Demonstrate pitch control.
Uses a planar multibody model
------------------------------------------------------------------------
See also PDDesign, Plot2D, TimeLabl, RK4, LoadAirfoilFile,
RHSVAWTPitchControl
------------------------------------------------------------------------
Contents
pitch =[ -34 59 5
-35 59 4
-35 58 4
-36 58 3
-36 58 3
-37 58 2
-37 57 2
-38 57 1
-38 57 1
-38 56 0
-39 56 0
-39 56 -1
-40 55 -1
-40 55 -2
-41 54 -2
-41 54 -3
-42 54 -3
-42 53 -4
-43 53 -4
-43 53 -4
-44 52 -5
-44 52 -5
-44 51 -6
-45 51 -6
-45 50 -7
-46 50 -7
-46 50 -8
-47 49 -8
-47 49 -6
-47 48 -3
-48 48 0
-48 47 3
-49 47 6
-49 47 8
-50 46 8
-50 46 7
-50 45 7
-51 45 6
-51 44 6
-52 44 5
-52 44 5
-53 43 4
-53 43 4
-53 42 4
-54 42 3
-54 41 3
-54 41 2
-55 40 2
-55 40 1
-56 39 1
-56 39 0
-56 38 0
-57 38 -1
-57 38 -1
-57 37 -2
-58 37 -2
-58 36 -3
-58 36 -3
-58 35 -4
-59 35 -4
-59 34 -5
-59 34 -5
-59 33 -6
-60 33 -6
-60 32 -7
-60 32 -7
-60 31 -7
-60 31 -8
-60 31 -8
-60 30 -9
-60 30 -9
-60 29 -10
-60 29 -10
-60 28 -11
-59 28 -11
-59 27 -12
-58 27 -12
-58 26 -13
-57 26 -13
-56 25 -14
-54 25 -14
-53 24 -15
-51 24 -16
-48 23 -16
-44 23 -17
-40 22 -17
-34 22 -18
-26 21 -18
-17 21 -19
-5 20 -19
0 20 -20
5 19 -20
17 19 -21
26 18 -21
34 18 -22
40 17 -22
44 17 -23
48 16 -23
51 16 -24
53 15 -24
54 14 -25
56 14 -25
57 13 -26
58 13 -26
58 12 -27
59 12 -27
59 11 -28
60 11 -28
60 10 -29
60 10 -29
60 9 -30
60 9 -30
60 8 -31
60 8 -31
60 7 -31
60 7 -32
60 7 -32
60 6 -33
59 6 -33
59 5 -34
59 5 -34
59 4 -35
58 4 -35
58 3 -36
58 3 -36
58 2 -37
57 2 -37
57 1 -38
57 1 -38
56 0 -38
56 0 -39
56 -1 -39
55 -1 -40
55 -2 -40
54 -2 -41
54 -3 -41
54 -3 -42
53 -4 -42
53 -4 -43
53 -4 -43
52 -5 -44
52 -5 -44
51 -6 -44
51 -6 -45
50 -7 -45
50 -7 -46
50 -8 -46
49 -8 -47
49 -6 -47
48 -3 -47
48 0 -48
47 3 -48
47 6 -49
47 8 -49
46 8 -50
46 7 -50
45 7 -50
45 6 -51
44 6 -51
44 5 -52
44 5 -52
43 4 -53
43 4 -53
42 4 -53
42 3 -54
41 3 -54
41 2 -54
40 2 -55
40 1 -55
39 1 -56
39 0 -56
38 0 -56
38 -1 -57
38 -1 -57
37 -2 -57
37 -2 -58
36 -3 -58
36 -3 -58
35 -4 -58
35 -4 -59
34 -5 -59
34 -5 -59
33 -6 -59
33 -6 -60
32 -7 -60
32 -7 -60
31 -7 -60
31 -8 -60
31 -8 -60
30 -9 -60
30 -9 -60
29 -10 -60
29 -10 -60
28 -11 -60
28 -11 -59
27 -12 -59
27 -12 -58
26 -13 -58
26 -13 -57
25 -14 -56
25 -14 -54
24 -15 -53
24 -16 -51
23 -16 -48
23 -17 -44
22 -17 -40
22 -18 -34
21 -18 -26
21 -19 -17
20 -19 -5
20 -20 0
19 -20 5
19 -21 17
18 -21 26
18 -22 34
17 -22 40
17 -23 44
16 -23 48
16 -24 51
15 -24 53
14 -25 54
14 -25 56
13 -26 57
13 -26 58
12 -27 58
12 -27 59
11 -28 59
11 -28 60
10 -29 60
10 -29 60
9 -30 60
9 -30 60
8 -31 60
8 -31 60
7 -31 60
7 -32 60
7 -32 60
6 -33 60
6 -33 59
5 -34 59
5 -34 59
4 -35 59
4 -35 58
3 -36 58
3 -36 58
2 -37 58
2 -37 57
1 -38 57
1 -38 57
0 -38 56
0 -39 56
-1 -39 56
-1 -40 55
-2 -40 55
-2 -41 54
-3 -41 54
-3 -42 54
-4 -42 53
-4 -43 53
-4 -43 53
-5 -44 52
-5 -44 52
-6 -44 51
-6 -45 51
-7 -45 50
-7 -46 50
-8 -46 50
-8 -47 49
-6 -47 49
-3 -47 48
0 -48 48
3 -48 47
6 -49 47
8 -49 47
8 -50 46
7 -50 46
7 -50 45
6 -51 45
6 -51 44
5 -52 44
5 -52 44
4 -53 43
4 -53 43
4 -53 42
3 -54 42
3 -54 41
2 -54 41
2 -55 40
1 -55 40
1 -56 39
0 -56 39
0 -56 38
-1 -57 38
-1 -57 38
-2 -57 37
-2 -58 37
-3 -58 36
-3 -58 36
-4 -58 35
-4 -59 35
-5 -59 34
-5 -59 34
-6 -59 33
-6 -60 33
-7 -60 32
-7 -60 32
-7 -60 31
-8 -60 31
-8 -60 31
-9 -60 30
-9 -60 30
-10 -60 29
-10 -60 29
-11 -60 28
-11 -59 28
-12 -59 27
-12 -58 27
-13 -58 26
-13 -57 26
-14 -56 25
-14 -54 25
-15 -53 24
-16 -51 24
-16 -48 23
-17 -44 23
-17 -40 22
-18 -34 22
-18 -26 21
-19 -17 21
-19 -5 20
-20 0 20
-20 5 19
-21 17 19
-21 26 18
-22 34 18
-22 40 17
-23 44 17
-23 48 16
-24 51 16
-24 53 15
-25 54 14
-25 56 14
-26 57 13
-26 58 13
-27 58 12
-27 59 12
-28 59 11
-28 60 11
-29 60 10
-29 60 10
-30 60 9
-30 60 9
-31 60 8
-31 60 8
-31 60 7
-32 60 7
-32 60 7
-33 60 6
-33 59 6
-34 59 5
-34 59 5];
Number of simulation steps
nSim = 2;
dT = 5e-3;
xPlot = zeros(15,nSim);
Initial state
x = zeros(8,1);
Blade angles
d = struct;
d.angle = [0 120 240]*pi/180;
d.rho = [0 .2 .2 .2;0 0 0 0;0 0 0 0];
d.d = [0 cos(d.angle);0 sin(d.angle);0 0 0 0];
d.wind.rho = 1.225;
d.wind.nu = 1.48e-5;
d.wind.v = [0;.0000001;0];
chord, half-length, span, thickness
d.blade = struct('c',0.0254*8, 'H', (0.3048/2)*5, 'Span',(0.3048/2)*10,...
'thknss',0.12*0.0254*8);
inr = 0.0230;
d.m = [0 1 1 1];
d.inr = [10 inr inr inr];
d.airfoil = LoadAirfoilFile( 'NACA0012.af' );
d.airfoil.alpha0 = 0;
Initialize time
t = 0;
omega = 6*pi;
wN = 55;
wD = 5*wN;
zeta = 0.7071;
xC1 = 0;
xC2 = 0;
xC3 = 0;
[aC, bC, cC, dC] = PDDesign( zeta, wN, wD, inr, dT );
thetaC1=zeros(1,nSim);
thetaC2=zeros(1,nSim);
thetaC3=zeros(1,nSim);
Run the simulation
for k = 1:nSim
phi = mod(omega*t,2*pi)*180/pi;
thetaC1(k) = interp1(0:360,pitch(:,1),phi)*pi/180;
thetaC2(k) = interp1(0:360,pitch(:,2),phi)*pi/180;
thetaC3(k) = interp1(0:360,pitch(:,3),phi)*pi/180;
yC = (x(2) - thetaC1(k));
uC1 = cC*xC1 + dC*yC;
xC1 = aC*xC1 + bC*yC;
yC = (x(3) - thetaC2(k));
uC2 = cC*xC2 + dC*yC;
xC2 = aC*xC2 + bC*yC;
yC = (x(4) - thetaC3(k));
uC3 = cC*xC3 + dC*yC;
xC3 = aC*xC3 + bC*yC;
d.t = [0;-uC1*inr;-uC2*inr;-uC3*inr];
xPlot(:,k) = [x*180/pi;d.t;thetaC1(k)*180/pi;thetaC2(k)*180/pi;thetaC3(k)*180/pi];
x = RK4(@RHSVAWTPitchControl,x,dT,t,d);
t = t + dT;
end
total_Torque=sum(xPlot(9,:));
ans =
0
v =
0
1e-07
0
ans =
2.0944
v =
8.6603e-08
-5e-08
0
ans =
4.1888
v =
-8.6603e-08
-5e-08
0
ans =
0
v =
-0.0095917
1e-07
0
ans =
2.0944
v =
-0.0095916
-5e-08
0
ans =
4.1888
v =
-0.0095917
-5e-08
0
ans =
0
v =
-0.0096036
1e-07
0
ans =
2.0944
v =
-0.0096035
-5.0002e-08
0
ans =
4.1888
v =
-0.0096037
-4.9998e-08
0
ans =
0
v =
-0.019231
1e-07
0
ans =
2.0944
v =
-0.019231
-5.0004e-08
0
ans =
4.1888
v =
-0.019231
-4.9996e-08
0
ans =
0
v =
-0.019231
1e-07
0
ans =
2.0944
v =
-0.019231
-5.0004e-08
0
ans =
4.1888
v =
-0.019231
-4.9996e-08
0
ans =
0
v =
-0.0222
1e-07
0
ans =
2.0944
v =
-0.0222
-5.0008e-08
0
ans =
4.1888
v =
-0.0222
-4.9992e-08
0
ans =
0
v =
-0.022217
1e-07
0
ans =
2.0944
v =
-0.022217
-5.0009e-08
0
ans =
4.1888
v =
-0.022217
-4.9991e-08
0
ans =
0
v =
-0.025204
1e-07
0
ans =
2.0944
v =
-0.025204
-5.0014e-08
0
ans =
4.1888
v =
-0.025204
-4.9986e-08
0
Labels
[t,tL] = TimeLabl( (0:(nSim-1))*dT );
xL = { '\theta_0' '\theta_1' '\theta_2' '\theta_3' ...
'\omega_0' '\omega_1' '\omega_2' '\omega_3' ...
'T_0' 'T_1' 'T_2' 'T_3'};
Plot2D(t, xPlot(1: 4,:), tL, xL(1: 4),'Pitch Control: Theta')
Plot2D(t, xPlot(13:15,:), tL)
Plot2D(t, xPlot(5: 8,:), tL, xL(5: 8),'Pitch Control: Omega')
Plot2D(t, xPlot(9:12,:), tL, xL(9:12),'Pitch Control: Torque')
