Annotation of OpenXM_contrib/gnuplot/demo/airfoil.dem, Revision 1.1.1.3
1.1 maekawa 1: #
1.1.1.3 ! ohara 2: # $Id: airfoil.dem,v 1.1.1.1.2.2 2001/12/04 14:12:53 lhecking Exp $
1.1 maekawa 3: #
4: # This demo shows how to use bezier splines to define NACA four
5: # series airfoils and complex variables to define Joukowski
6: # Airfoils. It will be expanded after overplotting in implemented
7: # to plot Coefficient of Pressure as well.
8: # Alex Woo, Dec. 1992
9: #
10: # The definitions below follows: "Bezier presentation of airfoils",
11: # by Wolfgang Boehm, Computer Aided Geometric Design 4 (1987) pp 17-22.
12: #
13: # Gershon Elber, Nov. 1992
14: #
15: # m = percent camber
16: # p = percent chord with maximum camber
17: pause 0 "NACA four series airfoils by bezier splines"
18: pause 0 "Will add pressure distribution later with Overplotting"
19: mm = 0.6
20: # NACA6xxx
21: thick = 0.09
22: # nine percent NACAxx09
23: pp = 0.4
24: # NACAx4xx
25: # Combined this implies NACA6409 airfoil
26: #
27: # Airfoil thickness function.
28: #
29: set xlabel "NACA6409 -- 9% thick, 40% max camber, 6% camber"
30: x0 = 0.0
31: y0 = 0.0
32: x1 = 0.0
33: y1 = 0.18556
34: x2 = 0.03571
35: y2 = 0.34863
36: x3 = 0.10714
37: y3 = 0.48919
38: x4 = 0.21429
39: y4 = 0.58214
40: x5 = 0.35714
41: y5 = 0.55724
42: x6 = 0.53571
43: y6 = 0.44992
44: x7 = 0.75000
45: y7 = 0.30281
46: x8 = 1.00000
47: y8 = 0.01050
48: #
49: # Directly defining the order 8 Bezier basis function for a faster evaluation.
50: #
51: bez_d4_i0(x) = (1 - x)**4
52: bez_d4_i1(x) = 4 * (1 - x)**3 * x
53: bez_d4_i2(x) = 6 * (1 - x)**2 * x**2
54: bez_d4_i3(x) = 4 * (1 - x)**1 * x**3
55: bez_d4_i4(x) = x**4
56:
57: bez_d8_i0(x) = (1 - x)**8
58: bez_d8_i1(x) = 8 * (1 - x)**7 * x
59: bez_d8_i2(x) = 28 * (1 - x)**6 * x**2
60: bez_d8_i3(x) = 56 * (1 - x)**5 * x**3
61: bez_d8_i4(x) = 70 * (1 - x)**4 * x**4
62: bez_d8_i5(x) = 56 * (1 - x)**3 * x**5
63: bez_d8_i6(x) = 28 * (1 - x)**2 * x**6
64: bez_d8_i7(x) = 8 * (1 - x) * x**7
65: bez_d8_i8(x) = x**8
66:
67:
68: m0 = 0.0
69: m1 = 0.1
70: m2 = 0.1
71: m3 = 0.1
72: m4 = 0.0
73: mean_y(t) = m0 * mm * bez_d4_i0(t) + \
74: m1 * mm * bez_d4_i1(t) + \
75: m2 * mm * bez_d4_i2(t) + \
76: m3 * mm * bez_d4_i3(t) + \
77: m4 * mm * bez_d4_i4(t)
78:
79: p0 = 0.0
80: p1 = pp / 2
81: p2 = pp
82: p3 = (pp + 1) / 2
83: p4 = 1.0
84: mean_x(t) = p0 * bez_d4_i0(t) + \
85: p1 * bez_d4_i1(t) + \
86: p2 * bez_d4_i2(t) + \
87: p3 * bez_d4_i3(t) + \
88: p4 * bez_d4_i4(t)
89:
90: z_x(x) = x0 * bez_d8_i0(x) + x1 * bez_d8_i1(x) + x2 * bez_d8_i2(x) + \
91: x3 * bez_d8_i3(x) + x4 * bez_d8_i4(x) + x5 * bez_d8_i5(x) + \
92: x6 * bez_d8_i6(x) + x7 * bez_d8_i7(x) + x8 * bez_d8_i8(x)
93:
94: z_y(x, tk) = \
95: y0 * tk * bez_d8_i0(x) + y1 * tk * bez_d8_i1(x) + y2 * tk * bez_d8_i2(x) + \
96: y3 * tk * bez_d8_i3(x) + y4 * tk * bez_d8_i4(x) + y5 * tk * bez_d8_i5(x) + \
97: y6 * tk * bez_d8_i6(x) + y7 * tk * bez_d8_i7(x) + y8 * tk * bez_d8_i8(x)
98:
99: #
100: # Given t value between zero and one, the airfoild curve is defined as
101: #
102: # c(t) = mean(t1(t)) +/- z(t2(t)) n(t1(t)),
103: #
104: # where n is the unit normal to the mean line. See the above paper for more.
105: #
106: # Unfortunately, the parametrization of c(t) is not the same for mean(t1)
107: # and z(t2). The mean line (and its normal) can assume linear function t1 = t,
108: # -1
109: # but the thickness z_y is, in fact, a function of z_x (t). Since it is
110: # not obvious how to compute this inverse function analytically, we instead
111: # replace t in c(t) equation above by z_x(t) to get:
112: #
113: # c(z_x(t)) = mean(z_x(t)) +/- z(t) n(z_x(t)),
114: #
115: # and compute and display this instead. Note we also ignore n(t) and assumes
116: # n(t) is constant in the y direction,
117: #
118:
119: airfoil_y1(t, thick) = mean_y(z_x(t)) + z_y(t, thick)
120: airfoil_y2(t, thick) = mean_y(z_x(t)) - z_y(t, thick)
121: airfoil_y(t) = mean_y(z_x(t))
122: airfoil_x(t) = mean_x(z_x(t))
123: set nogrid
1.1.1.3 ! ohara 124: set nozeroaxis
1.1 maekawa 125: set parametric
126: set xrange [-0.1:1.1]
127: set yrange [-0.1:.7]
128: set trange [ 0.0:1.0]
129: set title "NACA6409 Airfoil"
130: plot airfoil_x(t), airfoil_y(t) title "mean line" w l 2, \
131: airfoil_x(t), airfoil_y1(t, thick) title "upper surface" w l 1, \
132: airfoil_x(t), airfoil_y2(t, thick) title "lower surface" w l 1
133: pause -1 "Press Return"
134: mm = 0.0
135: pp = .5
136: thick = .12
137: set title "NACA0012 Airfoil"
138: set xlabel "12% thick, no camber -- classical test case"
139: plot airfoil_x(t), airfoil_y(t) title "mean line" w l 2, \
140: airfoil_x(t), airfoil_y1(t, thick) title "upper surface" w l 1, \
141: airfoil_x(t), airfoil_y2(t, thick) title "lower surface" w l 1
142: pause -1 "Press Return"
143: set title ""
144: set xlab ""
1.1.1.2 maekawa 145: set key box
1.1 maekawa 146: set parametric
147: set samples 100
148: set isosamples 10
149: set data style lines
150: set function style lines
151: pause 0 "Joukowski Airfoil using Complex Variables"
152: set title "Joukowski Airfoil using Complex Variables" 0,0
153: set time
154: set yrange [-.2 : 1.8]
155: set trange [0: 2*pi]
156: set xrange [-.6:.6]
157: zeta(t) = -eps + (a+eps)*exp(t*{0,1})
158: eta(t) = zeta(t) + a*a/zeta(t)
159: eps = 0.06
160: a =.250
161: set xlabel "eps = 0.06 real"
162: plot real(eta(t)),imag(eta(t))
163: pause -1 "Press Return"
164: eps = 0.06*{1,-1}
165: set xlabel "eps = 0.06 + i0.06"
166: plot real(eta(t)),imag(eta(t))
167: pause -1 "Press Return"
1.1.1.2 maekawa 168: reset
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to airfoil.dem CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gnuplot / demo