Viterna-Corrigan model: Code, sample run, & output
Matlab code file “makefinitedata.m”
%Function to convert airfoil data from infinite to finite aspect ratio
%>for symmetry airfoil from 0 to 90 degree angle of attack
%>Lanchester-Prandtl model was used before stalling angle
%>while Viterna-Corrigan model predict post stall coefficients
%
%INPUT
%infinitefile = filename of infinite 2D airfoil data
%>data in form of [degree Cl Cd] starting from 0 degree
%>to at least a point after maximum lift (stall) are needed
%A = airfoil aspect ratio
%
%OUTPUT
%finitedata = table of modified airfoil data with same format as input
%>data after stall are 'generated' at 1 degree increment till 90degree
function finitedata = makefinitedata(infinitefile, A);
infdata = load(infinitefile);
%Find effective angle of attack and Cd while Cl remain constant
%until stall occurs using Lanchester-Prandtl Eq.
i = 1;
dropping = 0;
%%Process first data then check if lift start dropping
while not(dropping) % need do while :(
%%Cl remain constant
finitedata(i, 2) = infdata(i, 2);
%%Cd increase or the lift slope decrease
finitedata(i, 3) = infdata(i, 3)...
+ (infdata(i, 2)) ^ 2 / pi / A;
%%Angle of attack increase or lift slope decrease
finitedata(i, 1) = (infdata(i, 1) / 180 * pi...
+ infdata(i, 2) / pi / A) / pi * 180;
dropping = (infdata(i + 1, 2) < infdata(i, 2));
i = i + 1;
end
%% i = first point after peak point
%Find finite lift and drag coefficient in post stall region
%using Viterna-Corrigon equation until 90 degree at 1 degree interval
alpha_S = finitedata(i - 1, 1); %all eq. refer to MODIFIED finite data!
C_LS = finitedata(i - 1, 2);
C_DS = finitedata(i - 1, 3);
if A > 50
C_Dmax = 2.01;
else
C_Dmax = 1.11 + 0.018 * A;
end
K_L = (C_LS - C_Dmax * sin(alpha_S/180*pi) * cos(alpha_S/180*pi))...
* sin(alpha_S/180*pi) / (cos(alpha_S/180*pi)) ^ 2;
K_D = (C_DS - C_Dmax...
* (sin(alpha_S/180*pi)) ^ 2) / cos(alpha_S/180*pi);
for alpha = (floor(alpha_S) + 1) : 1 : 90
finitedata(i, 1) = alpha;
finitedata(i, 2) = C_Dmax / 2 * sin(2 * alpha/180*pi)...
+ K_L * (cos(alpha/180*pi)) ^ 2 / sin(alpha/180*pi);
finitedata(i, 3) = C_Dmax * (sin(alpha/180*pi)) ^ 2 ...
+ K_D * cos(alpha/180*pi);
i = i + 1;
end
A sample run of “makefinitedata.m” in Matlab
» fin = makefinitedata('NACA12Re0300.data', 10)
fin =
0 0 0.0085
2.3644 0.1998 0.0107
5.9115 0.4998 0.0208
8.9029 0.7692 0.0369
12.4617 0.8563 0.0716
12.5648 0.8580 0.1163
13.0000 0.8457 0.1205
14.0000 0.8219 0.1305
15.0000 0.8033 0.1411
16.0000 0.7889 0.1525
17.0000 0.7779 0.1644
18.0000 0.7695 0.1771
19.0000 0.7634 0.1903
20.0000 0.7589 0.2041
21.0000 0.7560 0.2185
22.0000 0.7541 0.2335
23.0000 0.7532 0.2491
24.0000 0.7530 0.2652
25.0000 0.7533 0.2817
26.0000 0.7541 0.2988
27.0000 0.7551 0.3163
28.0000 0.7562 0.3343
29.0000 0.7574 0.3527
30.0000 0.7587 0.3716
31.0000 0.7598 0.3907
32.0000 0.7607 0.4103
33.0000 0.7615 0.4302
34.0000 0.7620 0.4503
35.0000 0.7621 0.4708
36.0000 0.7619 0.4915
37.0000 0.7614 0.5124
38.0000 0.7604 0.5336
39.0000 0.7589 0.5549
40.0000 0.7570 0.5764
41.0000 0.7545 0.5980
42.0000 0.7515 0.6197
43.0000 0.7480 0.6414
44.0000 0.7440 0.6632
45.0000 0.7393 0.6851
46.0000 0.7341 0.7069
47.0000 0.7283 0.7286
48.0000 0.7218 0.7503
49.0000 0.7148 0.7719
50.0000 0.7071 0.7934
51.0000 0.6989 0.8147
52.0000 0.6900 0.8359
53.0000 0.6805 0.8569
54.0000 0.6704 0.8776
55.0000 0.6597 0.8981
56.0000 0.6483 0.9183
57.0000 0.6364 0.9382
58.0000 0.6239 0.9578
59.0000 0.6108 0.9770
60.0000 0.5971 0.9958
61.0000 0.5828 1.0143
62.0000 0.5680 1.0323
63.0000 0.5527 1.0498
64.0000 0.5368 1.0669
65.0000 0.5204 1.0835
66.0000 0.5035 1.0996
67.0000 0.4861 1.1152
68.0000 0.4682 1.1302
69.0000 0.4499 1.1446
70.0000 0.4312 1.1585
71.0000 0.4121 1.1717
72.0000 0.3925 1.1843
73.0000 0.3726 1.1963
74.0000 0.3523 1.2076
75.0000 0.3317 1.2182
76.0000 0.3109 1.2282
77.0000 0.2897 1.2375
78.0000 0.2682 1.2460
79.0000 0.2466 1.2538
80.0000 0.2247 1.2609
81.0000 0.2026 1.2673
82.0000 0.1804 1.2729
83.0000 0.1580 1.2777
84.0000 0.1356 1.2818
85.0000 0.1130 1.2851
86.0000 0.0904 1.2877
87.0000 0.0678 1.2894
88.0000 0.0452 1.2904
89.0000 0.0226 1.2906
90.0000 0.0000 1.2900
» save('A10.NACA12Re0300.data', 'fin', '-ascii')
» plot(fin(:,1), fin(:,2), fin(:,1), fin(:,3))
Plotted graph of the finite airfoil data output
|
