1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
|
require 'narray'
require 'pgplot'
include NMath
include Pgplot
#-----------------------------------------------------------------------
# Demonstration program for PGPLOT vector field plot.
#
# Program to demonstrate the use of PGVECT along with
# PGCONB by illustrating the flow around a cylinder with circulation.
#-----------------------------------------------------------------------
twopi = 2*Math::PI
blank = -1.0e10
# number of points in the x and y directions
nx = 31
ny = 31
# cylinder radius
a = 1.0
# circulation strength
gamma = 2.0
# freestream velocity
vinf = 1.0
# max and min x and y
xmax = 3.0*a
xmin = -3.0*a
ymax = 3.0*a
ymin = -3.0*a
# point spacing
dx = (xmax-xmin)/(nx-1)
dy = (ymax-ymin)/(ny-1)
# compute the stream function, Cp, and u and v velocities
a2 = a**2
x = NArray.sfloat(nx,1).indgen!*dx + xmin
y = NArray.sfloat(1,ny).indgen!*dy + ymin
r2 = x**2 + y**2
r2[(r2.eq 0).where] = 1e-10
psi = vinf * y * (1-a2/r2) + gamma/twopi*0.5*log(r2/a)
u = vinf * (1 + a2/r2 - 2*a2*(x/r2)**2) + gamma/twopi * y/r2
v = vinf * x * (-2*a2*y/r2**2) + gamma/twopi * x/r2
cp = 1 - (u**2+v**2)/vinf**2
idx = (r2 < a2).where
u[idx] =
v[idx] = blank
#
# start drawing
#
pgbeg
pgenv( x[0], x[-1], y[0], y[-1], 1 )
pgiden
pglab('X','Y','Flow About a Cylinder with Circulation')
#
# contour plot of the stream function (streamlines)
#
tr = [ x[0]-dx, dx, 0.0, y[0]-dy, 0.0, dy ]
level = [ 1.0, 0.5, 0.0, -0.5, -1.0 ]
pgcont( psi, level, tr )
#
# draw cylinder
#
pgbbuf
pgsci(0)
pgsfs(1)
pgcirc(0.0, 0.0, a*1.1)
pgsfs(2)
pgsci(14)
pgcirc(0.0, 0.0, a)
pgsci(1)
pgebuf
#
# vector plot
#
pgsah(2, 45.0, 0.7)
pgsch(0.3)
tr = [ x[0], dx, 0.0, y[0], 0.0, dy ]
pgvect( u[1..-2,1..-2], v[1..-2,1..-2], 0.0, 0, tr, -1.0e10 )
pgsch(1.0)
#
# finish
#
|
