Korn System
Korn System
in Economics and
Mathematical Systems
Managing Editors: M. Beckmann and H. P. Kunzi
Control Theorv
(NASA-CR-142229) SUPERCRITICAL WING N75-1816'/
SECTIONS 2, VOLUME 108 (New York Univ.)
301 p HC $9.25 CSCL 01A
Unclas
00/01 12410
108
Springer-Verlag
Berlin • Heidelberg • New York
Lecture Notes in Economics and Mathematical Systems
(Vol. 1-15: Lecture Notes in Operations Research and Mathematical Economics, Vol. 16-59: Lecture
Notes in Operations Research and Mathematical Sy stems)
Vol. 1: H. Buhlmann, H. Loeffel, E. Nievergelt, Einfuhrung in die Vol. 30: H. Noltemeier, Sensitmtatsanalyse bei diskreten Imearen
Theorie und Praxis der Entscheidung bei Unsicherheit 2 Auflage, Optimierungsproblemen VI, 102 Seiten 1970 DM 16,-
IV. 125 Seiten 1969. DM 16,- Vol. 31: M. Kuhlmeyer, Die nichtzentrale t-Verteilung. II, 106 Sei-
Vol. 2: U. N. Bhat, A Study of the Queuemg Systems M/G/1 and ten. 1970 DM 16.-
GI/M1 VIII, 78 pages 1968. DM 16.- Vol. 32: F. Bartholomes und G. Hotz, Homomorphismen und Re-
Vol. 3: A. Strauss. An Introduction to Optimal Control Theory. VI. duktionen Imearer Sprachen. XII, 143 Seiten. 1970 DM 16,-
153 pages. 1968. DM 16,- Vol. 33: K. Hinderer, Foundations of Non-stationary Dynamic Pro-
Vol. 4: Branch and Bound: Eine Einfuhrung. 2 . geanderte Auflage. gramming with Discrete Time Parameter. VI, 160 pages. 1970.
Herausgegeben von F. Wemberg. VII. 174 Seiten. 1972 DM 18,- DM 16,-
Vol. 5: Hyvarinen, Information Theory for Systems Engineers VIII, Vol 34: H. Stbrmer. Semi-Markoff-Prozesse mil endlich vielen
205 pages 1968 DM 16.- Zustanden. Theorie und Anwendungen. VII, 128 Seiten. 1970.
DM 16,-
Vol. 6: H. P. Kiinzi, O. Miiller, E Nievergelt, Einfuhrungskursus in
die dynamische Programmierung. IV, 103 Serten. 1968. DM 16,- Vol 35: F Ferschl, Markovketten VI, 168 Seiten 1970 DM 16.-
Vol. 7: W. Popp, Einfuhrung in die Theone der Lagerhaltung. VI, Vol 36: M. P. J Magill. On a General Economic Theory of Motion.
173 Seiten 1968 DM 16,- VI, 95 pages 1970 DM 16,-
Vol. 8: J Teghem, J. Loris-Teghem. J. P. Lambotte. Modeles Vol. 37: H Muller-Merbach, On Round-Off Errors in Linear Pro-
d'Attente M/G/1 et GI/M/1 a Arnvees et Services en Groupes IV. gramming. VI. 48 pages. 1970. DM 16,-
53 pages 1969 DM 16.- Vol. 38: Statistische Methoden I. Herausgegeben von E. Walter.
Vol. 9: E. Schultze, Einfuhrung indie mathematischen Grundlagen VIII, 338 Seiten 1970 DM 22,-
der Informationstheorie VI, 116 Seiten 1969. DM 16,- Vol. 39: Statistische Methoden II. Herausgegeben von E. Walter.
Vol. 10: D. Hochstadter, Stochastische Lagerhaltungsmodelle. VI. IV. 155 Seiten. 1970 DM 16.-
269 Seiten. 1969 DM 18,- Vol. 40: H. Drygas, The Coordinate-Free Approach to Gauss-
Vol. 11/12: Mathematical Systems Theory and Economics Edited Markov Estimation VIM, 113 pages. 1970. DM 16,-
by H. W. Kuhn and G P. Szego VIII, IV, 486 pages. 1969 DM 34,- Vol. 41: U- Ueing, Zwei Losungsmethoden fur nichtkonvexe Pro-
Vol. 13: Heunstische Planungsmethoden Herausgegeben von grammierungsprobleme. VI, 92 Seiten 1971 DM 16,-
F Wemberg und C. A. Zehnder. II. 93 Seiten. 1969 DM 16,- Vol. 42: A. V. Balakrishnan, Introduction to Optimization Theory in
Vol. 14: Computing Methods in Optimization Problems Edited a Hilbert Space IV. 153 pages 1971 DM 16,-
by A. V. Balakrishnan. V, 191 pages 1969 DM 16,- Vol. 43: J. A. Morales, Bayesian Full Information Structural Analy-
Vol. 15: Economic Models, Estimation and Risk Programming: ses VI, 154 pages. 1971. DM 16.-
^Essays in Honor of Gerhard Ttntner. Edited by K. A. Fox, G. V. L. Vol. 44: G. Feichtinger, Stochastische Modelle demographischer
Narasimham and J. K. Sengupta VIM, 461 pages. 1969 DM 24,- Prozesse. XIII, 404 Seiten 1971 DM 28,-
Vol 16: H. P Kunzi und W Oetlli. Nichtlmeare Optimierung: Vol. 45: K. Wendler, Hauptaustauschschritte (Principal Pivoting).
Neuere Verfahren. Bibliographic IV. 180 Seiten 1969 DM 16.- II. 64 Seiten 1971 DM 16,-
Vol. 17: H. Bauer und K. Neumann, Berechnung optimaler Steue- Vol. 46: C- Boucher, Lecons sur la theorie des automates ma-
rungen, Maximumprinzip und dynamische Optimierung. VIM, 188 thematiques VIM, 193 pages. 1971. DM 18,-
Seiten 1969 DM 16,-
Vol. 47: H. A. Nour Eldin, Optimierung linearer Regelsysteme
Vol. 18: M. Wolff, Optimale Instandhaltungspolitiken in einfachen mit quadratischer Zielfunktion. VIII, 163 Seiten. 1971. DM 16,-
Systemen. V. 143 Seiten 1970 DM 16,-
Vol. 48: M. Constam. FORTRAN fur Anfanger. 2. Auflage VI,
Vol. 19: L. Hyvarinen Mathematical Modeling for Industrial Pro- 148 Seiten. 1973. DM 16,-
cesses VI, 122 pages 1970. DM 16,-
Vol. 49: Ch. SchneeweiB, Regelungstechmsche Stochastische
Vol. 20: G. Uebe. Optimale Fahrplane IX, 161 Seiten. 1970. Optimierungsverfahren. XI. 254 Seiten. 1971. DM 22,-
DM 16,-
Vol 50: Unternehmensforschung Heute - Ubersichtsvortrage der
Vol. 21: Th. Liebling, Graphentheorie in Planungs- und Touren- Zuricher Tagung von SVOR und DGU, September 1970. Heraus-
problemen am Beispiel des stadtischen StraBendienstes. IX, gegeben von M. Beckmann. VI. 133 Seiten. 1971. DM 16,-
118 Seiten. 1970 DM 16,-
Vol. 51: Digitale Simulation Herausgegeben von K. Bauknecht
Vol. 22: W. Eichhorn, Theorie der homogenen Produktionsfunk- und W Nef. IV, 207 Seiten 1971. DM 18,-
tion. VIM, 119 Seiten. 1970 DM 16,-
Vol 52: Invariant Imbedding. Proceedings of the Summer Work-
Vol. 23: A. Ghosal, Some Aspects of Queuemg and Storage shop on Invariant Imbedding Held at the University of Southern
Systems. IV, 93 pages. 1970 DM 16.- California. June-August 1970. Edited by R. E. Bellman and E. D.
Vol. 24: Feichtinger Lernprozesse in stochastischen Automaten. Denman IV, 148 pages 1971. DM 16,-
V. 66 Seiten 1970. DM 16.- Vol. 53: J. Rosenmiiller, Kooperative Spiele und Markte. IV, 152
Seiten. 1971. DM 16,-
Vol. 25: R. Henn und O. Opitz, Konsum- und Produktionstheorie.
I. II. 124 Seiten. 1970 DM 16,- Vol. 54: C C von Weizsacker, Steady State Capital Theory. Ill,
102 pages 1971 DM 16,-
Vol. 26: D Hochstadter und G Uebe, Okonometrische Methoden.
XII, 250 Seiten. 1970. DM 18,- Vol. 55: P. A. V. B. Swamy, Statistical Inference in Random Coef-
ficient Regression Models. VIM, 209 pages. 1971. DM 20,-
Vol. 27: I. H. Mufti, Computational Methods in Optimal Control Vol. 56: Mohamed A. EI-Hodiri, Constrained Extrema. Introduction
Problems. IV, 45 pages 1970. DM 16,- to the Differentiate Case with Economic Applications. Ill, 130
Vol. 28: Theoretical Approaches to Non-Numerical Problem Sol- pages 1971. DM 16,-
ving. Edited by R. B. Banerji and M. D. Mesarovic. VI, 466 pages. Vol 57: E. Freund, Zeitvariable MehrgroBensysteme. VII, 160 Sei-
1970. DM 24,- ten. 1971. DM 18,-
Vol. 29: S. E Elmaghraby, Some Network Models in Management Vol. 58: P. B. Hagelschuer, Theorie der linearen Dekomposition
Science. Ill, 177 pages. 1970. DM 16,- VII, 191 Seiten. 1971 DM 18.-
Control Theory
108
Springer-Verlag
Berlin • Heidelberg • NewYork1975
Editorial Board
H. Albach • A. V. Balakrishnan • M. Beckmann (Managing Editor) • P. Dhrymes
J. Green • W. Hildenbrand • W. Krelle • H. P. Kiinzi (Managing Editor) • K. Ritter
R. Sato • H. Schelbert • P. Schonfeld
Managing Editors
Prof. Dr. M. Beckmann Prof. Dr. H. P. Kunzi
Brown University Universitat Zurich
Providence, Rl 02912/USA 8090 Ziirich/Schweiz
Authors
Dr. Frances Bauer • Prof. Paul Garabedian
Dr. David Korn • Prof. Antony Jameson
New York University
Courant Institute of Mathematical Sciences
251 Mercer Street
New York, N.Y. 10012/USA
•t
Library of Congress Cataloging in Publication Data
This work is subject to copyright. All rights are reserved, whether the whole
or part of the material is concerned, specifically those of translation,
reprinting, re-use of illustrations, broadcasting, reproduction by photo-
copying machine or similar means, and storage in data banks.
Under § 54 of the German Copyright Law where copies are made for other
than private use, a fee is payable to the publisher, the amount of the fee to
be determined by agreement with the publisher.
© by Springer-Verlag Berlin • Heidelberg 1975. Printed in Germany
New York, N. Y.
November 1974
Chapter I. Theory 1
1. Introduction 1
2. Models of Shock Structure 2
3. Iterative Schemes for Three Dimensional Analysis. 11
4. Choice of Coordinates and Conformal Mapping 17
5. Two Dimensional Analysis with a Turbulent Boun-
dary Layer Correction 22
6. Design in the Hodograph Plane: A New Model of
the Trailing Edge 25
7. Design in the Hodograph Plane: Choice of Para-
meters 28
8. Bibliography 33 .
1. Introduction
designed have been tested with some success, and satisfactory agree-
ment of the results of our analysis with experimental data has been
and boundary layer effects. We hope that the data we have compiled
\
displacement thickness computed by a semi-empirical turbulent
shock waves are not defined sharply. Our contention is that the
the assumption that the reader has some familiarity with Volume I.
mathematics.
<*x>x= °'
|B - A | <_ C(b - a) .
The answer consists of two straight lines with the slopes C and -C
_ a+b B-A
X
0 ~ 2 2C~~ '
(See Figure la.) The problem has an analogy with transonic aero-
discuss the method of Murman and Cole [13] . Let equally spaced
mesh points be laid down on the interval [a,b] and denote by <j> .
the values of the potential <|> at these points. We call the jth
< an< >
point subsonic when 4>- + i ^--i 3 supersonic when $•+•> ^-.-i-
According to one version of the scheme of Murman and Cole our
form
tion
(a) Exact solution.
> >
*j - *j-2 *j ~ *j+2 °
hold at that mesh point. (See Figure Ib.) Moreover, there are
valid solutions containing a segment of shock points on which c)> .
remains constant. These smeared shock waves terminate with one
higher value and then a downturn leading to a supersonic point
and a shock point followed by subsonic points. (See Figure Ic.)
They need not fulfill any shock relations whatever, and they seem
to occur in the applications.
One way to. remedy the situation we have just described would
be to replace the scheme of Murman and Cole by a finite difference
analogue of the ordinary differential equation
(
*x}x = h
*xxx '
which is in conservation form and has been provided with an arti-
ficial viscosity term on the right. The small positive factor h
should be of the same order of magnitude as the mesh size. The
general solution of this equation is
X-X-
<j> = - h log cosh(—— ) +
=
d>
v $ h £ d>
x *xx ^x
tiated.
h
<*x>x= tx" ( E fe4>x> '
in which e is now differentiated. The appearance of e in the last
vative <f> should approach zero from the left at XQ , but may be
that
2 1 - e exp(x-xQ)/h
= C
^x 1 - exp(a-xQ)/h
nonlinear relation
1/2
B-A t fl-£ exp(x-x0)/hj dx
C j [l-exp(a-x0)/h J 2e-l '
In the limit as h -*• 0 this reduces to our earlier formula for x..
(*j+r*j)2 -(Wi)2 =P
J "pj-i"
where
p.. = max {0, (<t> j + 1 -<t>j_ 1 >} <<f>j+i~ 2 < f > j + l *'j-i> •
condition.
Now consider the problem of calculating the transonic flow
treating the natural boundary condition on <(> and the free surface
vanishes when the flow is subsonic but is positive when the flow is
by the decision process of Murman and Cole, while the highest order
derivatives appearing in the artificial viscosity are equivalent to
sional flow past an airfoil with a single valued stream function i|),
global considerations show that the total mass flux diji is actual-
ly conserved across the shocks even when they do not satisfy the
exact shock condition. Thus the method of Murman and Cole provides
a good approximation to the flow at nearly sonic speeds.
[(<f>x-cjd4; + p dy] .
The jump of the integrand across a shock wave is of the third order
of Murman and Cole that tends to yield shock waves behind which the
speed drops barely below the speed of sound through a jump roughly
10
tions of motion like those that have been described above. For the
similar to the ones obtained the old way. Some examples appear
occurs (cf. Chapter II, Section 2). Where the flow is expanding in
grids suggest that the truncation error remains quite small, on the
for Mach numbers that are about 0.02 smaller than those observed
specific problem. The method has been applied both in two dimen-
access to his experimental data for comparison with the theory [6].
speed [Link] across any shock wave, and we use the standard
(c2-q2)(j>ss + c2(A<|>-<j>ss) = 0 ,
where A<)> denotes the Laplacian of <J>. Since the direction cosines
of the stream direction are u/q, v/q, and w/q, the streamwise
second derivative can be expressed in the form
1 2 2 2
6
ss = -^=-*r
T
2 (u *xx
d> + v T<j>
yy + w yzz
4 + .
2uv<j>
y
xy + 2vw<t>
y
yz + 2uwd>
y
xz ) .
avoided.
(Ax)2
_ At, + .
XX AX Xt t.
neglect lower order terms, its principal part will have the form
on the split between new and old values in the difference scheme.
/U2 ,. . . . , ±_
(M -1; <p_^ ~ <l)™_, ~ 0 ~ i—o— ~ ap ~
+ 2 At
XX ZiX X"C
max Uu|,|v|,|w|J *t zt
+ lK
$ 1
i ,- H
i_v
,K"_^ ! 1
-i,
V-
_]~ ,"f"
_ JC_ -i-l
_ 1n ~k±r J3fT-'\
_ -\J-k -J-, 3 >K
Ax At
£ a. <j> . + r<j>
2 2 3 c2-q2 3
c V * - I <DXi$x.*XiX. = I h. ^- e g 3Xi -
h Ax
i = l*x.l i '
such as
relaxation factor.
The proposition has the advantage that it breaks up into
artificial time that specify the iterative scheme we use. The more
tion option for the two dimensional program with boundary layer
the extra terms in the equations that would result from the use of
of the flow over a yawed wing we have therefore used a square root
the wing about a singular line just inside the leading edge, which
x + iy = (X + iY)2 .
18
have a fast and accurate method of doing the conformal mapping. The
practice.
circle using polar coordinates r and to. The modulus h of the map-
explicit inversion.
s
19
tion of the map function around any circle exterior to the unit
z_ — z. = mf -=— , = 2iric
dz da _ . ~.
2 1 j da
Thus the mapping represents the wake as a gap with a constant thick-
n=0 a11
If a and s are the tangent angle and arc length of the contour in
N
ds
log T— = y (a cos noo + b sin noo) ,
doo rt n n
n=0
N
a - u) = £ (bn cos nu - an sin noo) ,
n=0
where
c
n = an + ib
n
-
formulas is of the order (1/K)
function
U
k = a2k - W2k + i
(a2k+l - U2k+l)
defined for 0 ~~
< k ~"
< K-l. Let U,
JC
be the complex Fourier trans-
VQ = 0 ,
-ito
Vk = Uk e , k > 0 .
Vk yield log (ds/du) at the shifted mesh points 2k+l and 2k+2,
ds
v
It 2k+l k '
2s s
U = cos -
S
0
and with splines to represent the contour has been found in prac-
Correction
« + (H + 2 - M2) i ^ = T
ds q ds
and the shape factor H and the skin friction T are determined
from semi -empirical formulas of Nash and Macdonald [14] . We
ness 6 . After that we update the map function in the unit circle
formulas for the shape factor H near the point where the boundary
layer separates.
According to the turbulent boundary layer method of Nash and
SEP
= -q-dl> -004 '
accurate and we have felt free to modify them. Thus over most of
over some prescribed interval near the trailing edge. Our idea is
the profile.
unit circle, even when the flow is computed at a mesh twice as fine,
layer we have drawn from the paper of Nash and Macdonald [14] or a
heavily aft loaded airfoil. We note that Bavitz [2] has also
the physical plane. New insight has been gained by experience and
cases with a boundary layer correction, for which the observed lift
was fifteen or twenty percent less than its predicted value. The
remain bounded on the upper surface. Heavy aft loading can still
This has the advantage of making the values of the pressure coeffici-
ip = 0 that proceed from the tail out to infinity and in effect delin-
eate the boundary layer wake. The new model of the trailing edge
thus obtained agrees with the one we have been using all along in
finite on the upper surface of the airfoil near the tail means that
tangent to the level curve of the speed q through the tail. There
are two different ways this can happen. First, we can impose a
profile ijj = 0 and with the angle of the flow monotonically increas-
ing as we pass from the upper surface to the lower surface. Both
gradient on the lower surface, and with the flow angles above and
aft loaded airfoil, and its success depends on the pressure coeffi-
cient C being nearly zero at the tail. Thus the speed at the
tail is almost the same as that at infinity and the flow angles
28
have been made, and they are listed in Section 8 of Chapter III.
We believe that the better model of the trailing edge which has
examples (cf. Section 1 of Chapter II), both before and after the
problem, and they furnish perhaps the best guide available to those
interested in the design method, which has turned out to be harder
able speed and slope, for the tail and lay down automation paths
through which the profile .ought to pass in the subsonic part of the
the location of the tail, and the more significant parameters defin-
ing the analytic function g(n) so as to arrive at a compatible
that had heavy aft loading failed to come up to their design speci-
fications in wind tunnel tests because we did not shape the profile
belief is that this fit should be carried far enough to ensure that
the inequality
SEP
6 dq
= - |dJi -004 '
conforming to the new criterion have more camber near the tail than
30
the computer program to design, with most of the runs using about
both our first example, the heavily aft loaded Airfoil 70-10-13
point with NCR = 5, required about 100 runs to perfect. The diffi-
tion conditions near- the tail. The problem of achieving smooth nose
the automation paths well short of the nose, and by choosing the
stagnation point so that they are more compatible with the auto-
mation paths.
sonic speeds attained within five percent of chord from the lead-
32
expected to reduce the drag creep that tends to occur just below
from the maximum possible design Mach number. This also tends to
suppress drag creep.
33
8. Bibliography
foils which we have been able to design. These are labelled with
the automation paths from Tape 6 have been included. This should
enable the reader to run the examples through Programs B and D and
to use them as starting points for new designs. For our newer and
geometry.
The newer airfoils are given first. The best are 79-03-12,
Airfoil 79-03-12 uses NCR = 5 (see pages 27-28) and has a low lift
(cf. Section 5). Airfoil 65-14-08 resulted from applying the same
criterion. Airfoils 70-11-12 and 65-15-10 are included largely for
series, and are included, not because they represent the best that
can currently be achieved, but because they have been tested (cf. [7,
at the NASA Langley Research Center. There are also plans for a two
70-10-13 at the NASA Ames Research Center, and for a two dimensional
space Corporation.
Our final example is a compressor blade which was designed in
c,
+ •* + V * VV+
t
+• •t +
+ + -t
^• + -t
o-o
•V
+ + *•,
.M •+
4-
.8
1.2
38
-J
CD
T
O
rv>
03
PO
a
-<
H
CD
CO
C~.I
PO
a)
39
07X23/7*
= -109
CIRCULATORY FLOW ABOUT A TRANjSONlC AIRFOIL
2 o
-.800 0.000 2
-1.000 0.000 2
2 0
.300 .050 2
,3ito -.062 2
TAPE 7
-6-109 4 -.12 .15 .08 1.40 .790 -.009 -.052 .120 1.50
22 1 2 5 6 10 13 It 17 18 33 34 37 38 42 49
50 53 54 57 58 &1 &2
AUTOMATION PftTHS
5 0
-.070 -.130 1
.190 -.335 -1
,2&0 -.350 3
.335 -.160 2
.350 -.090 2
3 0
-.910 -.300 -1
-.630 -.400 2
-.700 -.500 2
3 0
-.700 -.500 -1
-.600 -.500 1
-.510 -."t70 t
3 0
-.HID .390 -1
-.580 .410 2
-.840 .250 2
4 0
.120 .295 -1
.240 .320 2
.310 .260 1
.390 .160 1
3 0
.390 .160 -1
.390 .060 1
.355 -,0<t5 2
41
-1.21
-.81
-•Ml
"
[Link]
-41
•81
1.21
C\J
O
CO
.—i
O
II
>-
a
en
LO
ii
_i
LJ
O
00
46
08/21/73
= -67
CIRCULATORY FLDW A30UT A TRANSONIC ftlRFOlL
2 0
..800 0.000 2
.1.000 0.000 2
2 0
.300 .050 2
.*75 -.270 2
TAPE 7
-8 -87 <4 -.12 .15 .06 l.<*0 .780 0.000 -.116 .0=>5 1.50
25 1 2 5 6 9 10 13 It 17 18 22 33 3<* 37 38
tl H? 49 5(5 53 5i* 57 56 &1 &2
A U T O M A T I O N PATHS
4 0
..950 -.290 -2
..B80 -.400 1
..600 -.450 1
-.700 -.495 1
2 0
-.700 -.495 -1
-.560 -.495 5
4 0
..100 -.200 1
.245 -.220 1
.500 -.300 .1
.245 -.318 15
4 0
..100 -.200 1
.245 -.2,20 1
.245 -.316 -1
.170 -.310 15
5 .0
..075 .240 -1
..ISO .270 1
-.270 .320 1
..390 .390 1
-.530 .390 1
4 0
-.530 .390 -I
..S85 .340 2
-.790 .290 1
..930 .220 1
4 0
..075 .240 .1
..080 ,2itO 1
.200 .280 1
.270 .260 1
4 0
.270 .260 .1
.410 ,190 1
.490 .040 1
.470 -.100 1
48
-1.21
-.Ml
[Link]
.8.4.
4
4
1.21
O
CD
CD
O
II
>-
a
en
o
CD
CJ
O
CM
r-
50
07/23/71*
RUN= -20
CIRCULATORY FLOW ABOUT A TRANSONIC AIRFOIL
TAPC 6« PATH o
2 0
-.600 0.000 2
-1.000 0.000 2
2 0
.300 0.000 2
.t30 -.350 2
TAPE 7
-8-020 4 -.12 .15 ,08 [Link] .720 .001 -.109 .060 1.50
23 1 2 5 (, 9 10 13 14 33 34 38 tl 42 45 46
49 5o 53 5i» 57 58 6l 62
A U T O M A T I O N PftTHS
5 0
-1 .000 -.290 -1
- .910 -.450 1
> .73n -.540 2
- .710 -.550 5
- .650 -.550 5
I 0
- .940 .240 -1
- .730 .390 1
5 0
- .126 .394 -1
- .210 .416 2
- .410 .460 1
- .620 .450 1
- .740 .420 1
5 0
- .096 .390 -1
0 .000 .3B2 2
.103 .392 2
.200 .410 1
.280 .400 1
3 0
.230 .400 -1
."mo .160 3
.480 -.030 1
3 0
.480 -.030 -1
.500 -.170 4
.450 -.310 4
4 0
- .100 -.170 1
.170 -.370 -1
.220 -.310 2
.280 -.270 2
4 0
- .100 -.170 1
.280 -.270 -1
.340 -.280 2
.400 -.320 2
52
-1.21
-.81
-.41
[Link]
•41
1.21
-J
CD
.CO
-<
II
o
rv>
o
57
07/23/71*
SUN= -12
CIRCULATORY FLOW ABOUT A TRANSONIC AIRFOIU
TAPE 6* PATH
2 0
-.800 0,000 2
-1.000 0.000 2
2 0
.300 -.050 2
.455 -.360 2
TAPC 7
A U T O M A T I O N PfiTHS
5 0
-1.005 -.335 -1
-.920 -,405 2
-.815 -.492 2
-.715 -,555 3
-.644 -.563 4
r a
-.138 ,339 -1
-.210 .355 2
-.410 .385 1
-.620 .373 1
-.740 ,341 1
-.620 .293 1
-.680 .230 1
5 0
-.098 .335 -1
0.000 ,324 2
.102 ,322 2
.200 ,325 1
.230 ,320 1
5 0
.230 ,320 -1
.400 .224 3
.490 -.035 1
.515 -.185 3
.480 -.350 3
!> 0
-.100 -,225 1
.190 -.390 -1
.212. -.335 2
.295 -.300 2
.360 -.315 2
.420 -.360 2
59
-.81
-.41
4. •(• + 4- +
+• 4-1
[Link] 4
+
+ "I-
.4
1.21
CM
CJ
O
C\J
O
CO
CD
CD
CJ
CD
CD
64
07/23/71*
3UN= -138
TAPE 6* PATH
2 0
..600 0.000 2
-1.000 0.000 2
2 0
.300 -.100 2
.180 -."HO 2
TAPE 7
-7-138 <* -.12 .20 ,08 1.40 .TOO .012 -.202 .040 1.50
25 1 2 5 6 1Q 13 14 18 2l 22 33 3<t 37 38 42
45 46 49 50 53 54 57 58 61 63
A U T O M A T I O N PftTHS
6 0
- .100 -.200 1
.210 -.410 -1
.265 -.360 1
.310 -.330 2
.330 -.350 2
.440 -.400 2
3 0
..850 -.500 -1
- .720 -.560 1
- .570 -.550 4
4 0
..850 -.500 -1
..930 -.410 3
..990 -.310 3
.1.010 -.260 1
3 0
..430 .310 -1
..650 .290 1
- .620 .220 1
Z 0
.200 .240 -1
'.400 .310 1
5 0
.200 .240 -1
.300 .240 1
.400 .190 1
.500 -.040 1
.530 -.200 1
3 0
.300 -.100 1
.530 -.200 -1
.510 -.350 2
66
a
a
o
o
o
o
a
o
ii
71
07/23/74
RUN= -85
TAPE 6t PATH o
Z 0
..600 0.000 2
-1.000 0.000 2
a o
.300 .050 2
.520 -.320 2
TAPE 7
-6-065 H -.12 .20 .06 1.40 .700 ,0i7 -.210 0.000 1.50 7
23 1 2 5 6 lo 18 21 22 33 34 37 36 42 45 46
49 50 53 54 57 58 61 &2
AUTOMATION PATHS
3 0 '
- .840 -.510 -1
..700 -.550 1
- .550 -.545 4
4_ 0
.640 -.510 -1
..910 -.400 3
..990 -.320 3
-1 .010 -.260 1
5 0
.100 .520 -1
..100 .520 1
- .350 ,550 1
5 0
..430 .310 -1
..650 .290 1
- .830 .220 1
& 0
.200 .260 -1
.300 .230 1
.390 .160 1
."*50 .060 1
.430 -.040 1
.480 -.200 1
5 0
..100 -.300 1
.300 -.300 1
.220 -.420 -1
.320 -.390 1
.too -.360 1
73
-1.2 1
-.81 4
4-
4
4
[Link]
^ 4 t .
.8
1.2
CD
cn
o
n
r~
I!
-C
CD
CD
^c
11
o
cn
o
CD
co.
75
07/24/74
RUN= -41
TAPE 6t PATH 5
2 0
-.600 0.000 a
-1.000 0.000 a
2 0
.400 -.200 2
.515 -.450 2
TAPE 7
-8-041 4 -.12 .25 .08 1.40 .650 ,0i8 -.371 .023 1.50
25 1 2 5 6 9 10 13 It 17 18 33 3t 37 48 12
45 ft i*9 5o 63 Si* 57 58 61 &2
AUTOMATION PATHS
2 0
-.790 -.610 -1
-.640 -.600 4
2 0
-.790 -.610 -1
-.940 -.460 4
2 0
-.940 -.460 -1
-1.100 -.260 3
•4 0
-.350 .560 -1
-.460 .630 3
-.570 .610 2
-.700 .570 2
•3 0
.200 .290 -1
0.000 .335 2
-.280 .490 1
-3 0
.200 .290 -1
.350 .220 1
.440 .060 1
-3 0
.440 .060 -1
.460 -.030 5
.440 . -.120 5
-5 0
-.100 -.120 1
.150 -.220 -1
.250 . -.120 1
,300 -.120 2
.350 -.140 2
77
-1-21
-.81
[Link]
-41
-81
1-21
3-
O
II
LJ
O
O
O
>-
D
CM
r^
zr
ii
_j
LJ
O
LD
CO
82
08/2i*/73
C I R C U L A T O R Y FLOW A 3 0 U T A TRANSONIC
2 0
..800 0.000 2
-1.000 0,000 2
2 o
.300 .050 2
.563 -.320 2
TAPE 7
-7-111* f -.12 .25 .08 [Link] .650 .006 -.320 0.000 1.50
23 i 2 5 6 9 10 17 18 33 51 37 36 <*2 t5 46
i*9 5o 53 5<» 57 58 61 c,z
ftUT0«HAfIQN PATHS
2 0
.,790 ».&00 -1
..&10 -.590 H
3 0
..790 -,&oO «1
.,880 »,5feO 3
.,950 -,»70 3
3 0
.,950 -,i»7Q -1
• I, 020 -.390 3
.1.080 -,2&0 I
3 0
,200 ,%30 -1
0,000 ,<H»0 1
,,280
3 0
.,350 ,5QP .1
,.5?0 ,5&0 2
.,710 ,550 2
(4 0
,200 »»go .1
,370 ,320 1
,*SO ,120 1
.460 .020 2
i* 0
..100 -.100 1
.350 -.160 -1
.290 -,820 1
,150 -.260 1
84
-1 -6_
-.81
-.41
[Link]
•81
-*•
4
1-21
CD
CD
O
O
r~
ii
no
CD
CD
-<
II
O
CD
89
07/16/74
RUN= -274
CIRCULATORY FLOW ABOUT A TRANSONIC AIRFOIL
2 0
-.800 0.000 2
.1.000 0.000 2
2 0
.330 .050 2
.570 -.325 2
TAPE 7
-8-274 2 -.12 .30 .08 1.40 .600 ,0n7 -.330 .025 .50
21 1 2 5 6 17 18 33 34 37 36 4i 42 f5 49 50
5s Sit 57 56 61 62
AUTONATION PATHS
4 0
-.800 0.000 1
•1.100 -.250 -1
-1.080 -.370 2
-1.000 -.500 2
5 0
-.600 0.000 1
-1.000 -.500 -1
-.930 -.575 2
-.614 -.620 2
-.691 -.625 «•
6 0
-.691 -.625 -1
-.580 -.630 2
-.490 -.620 2
-.380 -.600 2
-.280 -.575 1
-.200 -.560 1
5 0
-.200 -.560 -1
-.100 -.545 1
-.054 -.540 1
.046 -.540 1
.110 -.545 1
3 0
-.250 .400 -1
-.650 • .500 2
-.600 .450 2
2 0
-.250 .400 -1
.100 .»15 2
5 0
.340 .038 -1
.480 .005 1
.530 -.125 1
.530 -.230 1
.475 -.255 1
6 0
-.100 -.365 1
.300 -.365 1
.495 -.350 -1
.385 -.440 1
.210 -.540 1
.110 -.545 1
91
++
^4
+ 4 4 4 +4 4 4 4 *
4
4-
o.o 1
-81
I-.21
GO
ro
o
n
f—
ii
en
CO
o
a
-<
ii
a
en
o
ii
o
U3
ro
93
07/18/71*
*UN= -255
TAPE 6» PATH 0
2 0
-.600 0.000 2
-1.000 0.000 2
2 0
.300 .050 2
.470 -.260 2
TAPE ?
-8-255 it -.12 .15 .08 1.40 .820 .On5 -.105 .055 ,50
23 1 2 5 6 10 13 It 17 IS 33 34 37 38 <*i 42
"*9 5Q 53 5t> 57 58 61 62
A U T O M A T I O N PftTHS
5 0
-.961* -.290 -1
-.930 -.350 1
".874 -.405 1
-.794 -.450 1
-.675 -.485 1
3 0
-.675 -.485 -1
-.550 -.470 2
-.460 -.450 5
5 0
-.800 0.000 2
-.875 .255 -1
-.832 .285 1
-.685 .350 1
-.560 .365 1
5 0
-.075 .250 -1
-.180 .260 1
-.270 .330 1
-.390 .360 1
-.520 .370 1
4 0
-.020 .250 -1
.060 .267 1
.200 .300 1
.270 .300 1
4 0
.270 .300 .1
,4io .190 1
.490 .040 1
.470 -.100 1
4 0
-.100 -.200 1
.245 -.220 1
.320 -.355 -1
.268 -.330 15
4 0
-.100 -.200 1
.245 -.220 1
.266 -.330 -1
.215 -.305 15
95
-1-21
C,
--81
4-
4-
-.41
o-ol
4- +
j.4
.41
-8.1
1.21
CD
CD
en
CM
UD
n
_i
LJ
—fr
O
LO
98
06/20/73
= -131
TAPt 6« PATH 0
2 0
..300 [Link] 1
-1.000 [Link] 1
2 0
.100 .200 1
.390 -.142 1
TAPE: 7
99,131 4 -.12 .40 .06 1.1*0 .750 .010 -.120 0.000 0.00
6 1 2 5 6 17 18
LISTING OF C O O R D I N A T E S F QR AIRFOIL 7 5 - 0 6 - 1 2
A\IG KftPPA
KAPPA MACH CP
_
54 .21022 ~ •05052 -2.50 .70 .9039 _ .3799
55 .19512 ™ • 04976 -3. 13 .76 .9032 .3783
56 .18044 ™ • 04839 -3. 30 .82 .90^5 . .3749
57 .16622 ™ • 04736 -4. 50 .90 .9036 ,367i»
56 .15247 ~ •04669 -5. 25 .99 .6994 - .357-4
59 .13920 ™ • 04536 -6. 04 1 .09 .6935 - .3439
60 .12644 ~ •04394 -6. 89 1 .20 .6852 -
. .3263
61 .11420 •• •04237 -7. 79 1.33 .8772 . .3049
62 .10250 —
• 04067 -e. 73 1 .49 .8555 _ .2797
63 _-
.09134 03636
"• -9. 75 1 .66 .8544 .2504
64 .08075 —
• 03693 -10. 34 1 .83 .8405 _ .2173
65 .07074 ~ •03491 -12. 01 2 .It .8251 _ .1804
66 .06131 ™ • 03230 -13. 29 2.49 .8032 .1397
67 .05250 —
• 03061 -14. 69 2.93 .7396 . .0950
69 .04431 02834
"• -1&. 25 3 .53 .7692 «• .046r
69 .03676 —
• 02602 -13. 03 4.37 .7468 .007e,
70 .02987 ™ • 02364 -20. 10 5.61 .7222 .0665
71 .02366 ™ • 02122 -22. 58 7.54 .6945 .1325
72 .01814 ™ • 01877 -25. 58 10 .71 .6623 .208?
73 .01336 —
t 01626 -29. 75 16 .17 .6232 .299o
74 .00932 " » 01370 -35. 37 25 .93 .5711 .4169
75 .00605 ™• 01105 -43. 14 33 .82 .4944 .581?
76 .00351 ™ • 00827 -52.45 44 .69 .3354 .789s
77 .00165 ~ •00536 -62. 46 62 .39 .2531 .9817
78 .00049 *" • 00240 -75. 49 70 .59 .1098 1 .1177
79 0.00000 •
00074 -86. 59 49 .26 .0423 1 .I44n
80 .00003 • 00399 -94. 60 43 .t3 .1373 1 .0590
91 .00056 •
00717 -104. 40 64 .t6 .3552 .8406
82 .00171 •
01015 -117. 71 74 .24 .5297 .5071
93 .00365 •01296 -130. 57 51 .62
84 .00640 01573 -139. 54 25 .91
.6490
.7141 _ .2393
•
.0859
85 .00993 t 01852 -144. 14 17 .39 .7651 .0386
Q6 .01415 *
02133 -146. 32 12 .04 .8118 . .1464
87 .01906 •
02415 -151. 72 9.11 .6532 . .2476
38 .02461 •02696 -154. 58 7.11 .8918 . .3394
89. .03060 •02973 -157. 07 5 .81 .9294 .425P
90 .03761 * 03245 -159. 30 4.87 .9635 _-
. .5075
91 .04503 *
03510 -161. 32 4.13 .9979 .5865
92 .05304 •03766 -163. 20 3 .74 1 .0331 . .6659
93 .06163 04009 -165. 07 3 .54 1 .0709 _ .749(t
•
94 .07080 U4237 -165. 98 3 .63 1 .1143 _-.842R
•
95 .08007 •
04444 -169. 05 3.37 1 .1659 .9501
96 .09101 •04629 -170. 74 2 .15 1 .2015 .1 .021*
97 .10215 •
04800 -171. 73 1 .27 1 .2122 -1 .0427
98 .11398 •
04963 -172. 56 1.01 1 .2139 -1 .046n
99 .12642 •05118 -173. 21 .62 1 .2120 -1 .0424
100 .13946 •
05256 -173. 79 .71 1 .2096 -1 .0357
101 .15306 •
05408 -174. 31 .63 1 .2044 -1 .0274
102 .16720 •
05543 -174. 79 .56 1 .1999 -1 .0185
103 .16183 •
05670 -175. 24 .51 1 .1952 -1 .0092
104 .19694 05791 -175. 67 ,t7 1 .1906 _ .9999
.
•
105 .21251 •05903 -176. 07 .43 1 .1360 .9907
106 .22850 •
06007 -176. 45 .40 1 .1315 .9817
107 .24486 •06104 -176. 82 .38 1 .1772 -
. .973(1
106 .26164 •06192 -177. 17 .36 1 .1730 .9646
-
101
-1.21
4
4
-•Ml
[Link]
4
4
•Ml
1-21
IT-
CD
II
>-
a
OD
UD
—»
o
LO
104
08/2H/73
RJ,N = -2"+2
CIRCULATORY FLOW A30UT A TRANSONIC AIRFOIL
TAPE 6. PATH 0
2 0
..300 0.000 2
-1.000 0.000 1
2 0
.300 0.000 1
.550 -.31+0 1
TAPE 7
-6-212 1 -.12 .25 .08 [Link] .750 .007 -.110 .060 .50
17 1 2 5 6 9 10 m 33 3<* 37 38 49 50 53 5>*
57 58
[Link] PATHS
4 0
..300 0.000 2
-.995 -.110 .1
..985 -.160 1
-.970 -.200 1
6 0
..918 -.260 .1
..8t8 -.390 1
-.763 -.445 1
..662 -.480 2
..562 -.495 2
-.1*78 -.495 2
4 0
-.093 -.210 1
.242 -.400 .1
.292 -.390 1
.353 -.365 1
15 0
-.098 .325 -1
..003 .315 1
.102 .305 1
.182 .300 1
.261 .270 1
.390 .183 1
.427 .140 1
.450 .060 1
.1*17 -.030 1
.443 -.095 1
.iff 2 -.150 1
.420 -.200 1
.382 -.230 1
.337 -.220 1
.297 -.208 1
7 0
..128 .330 -1
-.213 .355 1
..323 .365 1
..418 .380 1
-.523 .375 1
..623 .360 1
-.713 .335 1
6 0
..SOO 0.000 2
..980 .110 .1
..935 .200 1
-.900 .250 1
..640 .300 1
..760 .330 1
106
0-0
-Ml
-81
1.21
-O
o
o
n
r~
n
-o
U>
CO
O
ro
UD
ro
o
-C
109
OB/2<*/73
= -79
C I R C U L A T O R Y FLOW A B O U T A T R A N S O N I C ftlRFOlu
6. PATH 0
2 0
..BOO 0.000 2
-1.000 0.000 1
2 0
.300 0.000 1
.580 •.160 1
UPE 7
A U T O M A T I O N PATHS
3 0
-1.000 -.290 -1
-.910 -.170 1
..780 -.575 2
2 0
-.780 -.575 -1
..&20 -.560 5
1 0
..BOO 0.000 2
..980 .150 .1
..930 .260 1
..780 .120 1
5 0
-.128 .115 .1
..320 .170 2
..110 .180 1
.-.620 .180 1
.110 1
3 0
-.098 .115 -1
0.000 .105 1
.120 .115 1
3 0
.320 .350 -1
.160 .150 5
.4&0 -.030 1
5 0
..093 -.210 1
.220 -.120 .1
.330 -.380 2
.130 -.300 1
.150 -.200 1
111
/xxx X X X X X
1.01
-61
4 4
•21.
0.0
CO
O
O
f\3
II
cn
en
o
CO
O
ii
cn
113
compare the analysis of the flow past an NACA 0012 airfoil by these
can be seen to give essentially the same shock jump, but the fully
tive scheme does not give the full shock jump but agrees better
good agreement provides evidence that our new model of the tail
should eliminate the loss of lift which was experienced with the air-
the fine grid, but the flow is almost shock free on the crude grid.
wing. The program uses the full three dimensional difference scheme
ding Mach numbers and yaw angles. First we compare crude and fine
the airfoil is operating below its design point and two shocks
coarse grid in the unyawed condition. Away from the shock waves,
section was designed for Mach .79, drag rise is only just beginning
effects.
115
-1.21
-.8
-.4
•01
-Ml
1-21
-1.2..
-.8..
-.4 „.
.0..
.4 ..
.8..
1.2..
-1.21
-.81
-ol
-Ml
•81
1.21
-1.21
-.81
-01
•81
1.21
-1.21
-.81
-.41
.0
.41
•8
1.21
.020,.
,. DELS
.010
.005..
.000
.020^
UPPER SURFRCE DELS
015.. DELS
.010..
.005..
.000
.2 ,M .6
-1.21
--Hi
-01
I-21
•>••*•<
a
cu
a
oa
a
3-
•
I
0°
0°
i
o
CD
f
I
flIRFOIL 65-15-10
M = .650 . flLP = 0.000
CL = 1.4866 CO = -.0001 CM = -.2474
123
o
10
a
C\J
" >*
4-
X
4-
o
CO
fllRFQIL 65-15-10
M = 1.000 flLP = 0.000
CL = -50M8 CD = .1154 CM = -.3198
124
o
10
CO
I
a
3-
!
|.+
(_>°l
-,* x i
K4-
4
»
a
ID 4
flIRFQIL 65-15-10
M = 1.200 flLP = 0.000
CL = -M07M CD = .1048 CM = -.2738
125
- 2 . 4 -,
-1 .8..
-1.2-
-.6..
.0-
1.2-
-2.M _
-1.8..
-L.2..
-.6..
-0..
L.2..
v
UPPER SURFHCE PRESSURE LOWER SURFflCE PRESSURE
Correction
were found with transition set at PCS = .07, and with LSEP = 161.
PCH = .07, but with LSEP = 153. The experimental data on pages 147
and 148 are from two different series of tests of Airfoil 75-06-12
-1.21
-84
--Ml
•ol
-Hi
•81
-1
-1.2..
-.8-=^
-.4 ..
.0
1.2..
-1.2..
-.M ..
.0..
1.2..
1.2..
-1 .2..
-.8--
-8
1.2-
-1.21
-.81
--Ml
.01
-Ml
-81
1.21
-1.2-
-.4 __
.0..
.4 __
1.2..
-1.2-
-.8..
_ (q
I.2..
-1.2-.
-.8 ._
-.4 _.
1.2..
-1.2..
-.8..
-.4 ..
.8..
1.2..
-1.2-
-.8
1.2..
-1-2-.
-.8..
-.4..
.0..
.4 ...
1 .2..
-1.2-
-.8..
-0'..
.8..
1.2..
-.8..
.0-
1.2-
-1.21
-.Hl
•Ml
.81
1.21
-1.2
-.8..
-.4
.0..
.4..
.8..
1.2..
-.81
--Ml
.81
1.21
-.8..
-.4
.8..
-1.2..
-.8..
-0..
1.2..
-1.21
--Ml
1.21
-1.2
C,
-.81
-.4 J
.01
-41
.81
1-2
-1 .6
-1.21
-.81
•01
-41
.81
1.21
wave well forward on the wing section where the boundary layer is
-1-21
-.41
.81
-1.21
-.81
-.H
.8
1-2
-L.2..
-.8..
-.4 ..
.0..
1 -2
-1.21
-.81
.Ml
.81
1.21
-1.2..
-.8..
-.4 ..
.0..
.4 ..
1.2..
-1 .2..
-.8..
.8.
1.2..
5. Drag Polars
without any boundary layer correction, and they involve just the
The final drag polars are for three dimensional flows past
found that our evaluation of L/D compares favorably with the test
angles with the envelope obtained from calculations in which the yaw
optimal lift drag ratio does not change much with the design Mach
.035 i— = -800
-030 \—
= -790
= -780
."760
-750
FIGURE L
160
• 035 i—
= -85
• 030 \—
- .-75
.65
-55
.45
.35
FIGURE 2
161
• OMO i—
-035
-005 CL
-OMO i—
CL - -5
CL - .M
.005
FIGURE M
163
.030 —
THEORY
.025 A EXPERIMENT
.020
.015
.010
.005
CL
-.20 .00 .20 .40 .60 .80 1.00
.030 THEORY
A EXPERIMENT
.025
.020
.015
CJ
.010
.005
CL
-.20 0.00 .20 .MO' .60 .80
FIGURE 6
165
.030
THEORY
A EXPERIMENT
-025
.020
• 015
a
CJ
.010
.005
.030 THEORY
A EXPERIMENT
• 025
• 020
-015
Q
CJ
.010
.005
MO r-
THEORY
EXPERIMENT flT
3 YRW RNGLES
FIGURE 9
168
RIRFOIL 90-10-13
flLRFOlL 78-06-10
II
O
X
d
10 —
FIGURE 10
169
6. Schlieren Photographs
the flow at a fairly high Mach number below the design lift, with
two shocks on the upper surface and a shock on the lower surface.
Figure 2 shows the flow slightly below the design point with two
shocks quite far back, and Figure 3 shows the nearest approach to
the design flow for this airfoil with one fairly weak shock.
Airfoil 82-06-09; they are British Crown Copyright, and are repro-
Figure 5 shows the flow above the design point with a single shock
point. The design pressure gradient was too severe near the tail
on the upper surface and the flow was strongly separated in the
test.
170
The program which was written for the analysis of the flow
explained in Volume I.
of viscosity.
4. To obtain a redistribution of airfoil coordinates.
the FSYM value. The deck structure and data structure correspond-
data or any other data we have provided the user with the option of
will generally not need to be changed from case to case have been
necessary; default values will be used for this entire program. The
If these default values are [Link] required for the case under
mapping onto the circle performed and the boundary layer removed.
The printout from this program is similar to the printout from Pro-
gram F. Upon termination Tape 3 will contain the output data in the
YS, the surface slope and curvature, the pressure distribution CP,
NS1 cycles. This process is repeated every NS1 cycles until the
NCY = NS, where NCY is the running tally of flow cycles, or until
The slopes are needed to perform the mapping and are obtained by
After the surface slopes are obtained the mapping to the circle is
for the coordinates after the first mapping. Of these M+l mapped
coordinates, 108 points (NT) are saved. The points are obtained by
thinning out every other point of the upper and lower surfaces near
cruder grid, M*N = 80x15. After the flow, the ordinary differen-
input to this equation is the set of local Mach numbers at the NPTS
points of the circle plane which are obtained from the flow calcula-
176
value XSEP, SEP is set equal to its computed values even if 'SEP is
greater than SEPM. XSEP = .93 and SEPM = .004 are the default
mapped onto the circle and then the flow calculations are resumed.
ing the Mach number distribution after separation for input to the
some point along the upper surface to a base pressure BCP at the
x = XMON to the trailing edge. The straight line for the pressure
distribution after separation is determined each time by the value
making use of this option and have, therefore, set the LSEP default
value to M+l. This means that the pressure distribution derived
from the flow will not be modified on the upper surface. In diffi-
manner.
but the boundary layer computed at the NPTS points on the upper
surface and lower surface is plotted. In the plot the upper
at the top of the page. The quantities printed out after each KP
eight variables change during the flow calculation. The last five
correction subroutines and remain constant for the NS1 flow cycles.
The flow program has been modified so that the flow can be computed
correction are less than the tolerance ST. The ITYP parameter is
the coefficients of lift, wave drag, form drag, total drag and
if available, the airfoil and the sonic line are plotted. A plot
of the final displacement thickness for the upper and lower surfaces
appears on the last page. If ITYP < 4 no Calcomp plots are made.
3.2 operating system and data cards for the boundary layer program
which retrieve the program from the program library, store the air-
REWIND(T,TAPE3)
H. Execution of Program H
See Table 1 for the format of data for Tape 3. See Table 2
The first data card which may be used as input to the program
is:
[ $P NS=1, FSYM=4., EM=0.762, CL=0.362, IS=2, PCH=0.07,
However, for this case the input data can be simplified consider-
values set by the program. Thus, the first input card we use is:
After reading this card the program begins the mapping. The air-
circle. From these coordinates 108 (NT) are saved: every other
point in the first third of the original 161, every point in the
second third (around the nose) and every other point in the last
third, which includes the points near the upper surface trailing
edge. This was done to maintain the resolution at the nose and
108 points define the inner airfoil to which we add the boundary
[ $P NS = - 1, ITYP = 1$ I"
180
cruder mesh. The mesh size is then M*N = 80*15 for the flow calcu-
This card initiates the computation of the flow around the airfoil.
400 (NS) cycles of the flow on the crude grid size will be computed
and Mach numbers resulting from flow cycle NS1 are computed. Since
LSEP is not equal to its default (81 for the crude mesh, 161 for
surface from LSEP+1 to the trailing edge. Since the pressure BCP at
the trailing edge was not read in on a data card, the default
value BCP = 0.4 is used initially. At all boundary layer cycles
after the first, BCP is iterated and. underrelaxed using the mono-
by IS. This is the same parameter name used for the smoothings of
the original airfoil, but that smoothing is not. done for each new
outer airfoil. The amount of <$ to be added to the original airfoil
factor. After a spline fit at the NPTS points at which the equa-
defined at 108 (NT) points is then mapped onto the circle. The
airfoil and sonic line, and the plot of the last upper and lower 6.
Since some default values are used, Card 3 can be shortened to':
[ $P NS=400, LSEP=75, ITYP=4, KP=4$ ]
The mesh is restored to the finer grid, MXN = 160x30. The inner air-
redefined on the fine grid and the new airfoil is obtained at 161
ponding index for the fine mesh. All other required variables are
interpolated.
The fifth data card is:
[ $P NS=400, ITYP=1$ ]
400 cycles of flow are done on the fine mesh with a boundary layer
[ $P ITYP=0$ ]
The time required to compute 400 cycles of the flow and obtain
to obtain each new outer airfoil and map it onto the circle. In
total about 65 seconds are spent on the flow and 45 seconds on the
outer airfoil. 307 seconds are required for the 400 cycles calcu-
lation at a fine mesh size. 231 seconds are needed for the flow
Final Printout
^v COLS.
1-10 11-20 2'1 - 30 31 - 40
CARDS ^^
Title in Hollerith
1
(Columns 2-17 will be printed on plot)
2 FNU FNL EPSIL
3 Blank
4 Coordinates at nose
Points on upper surface
,*
FNU + 3 Coordinates at trailing edge
FNU + 4 Blank
Deck Structure
"V. COLS.
^V. 1 - 10 11 - 20 21 - 30 31 - 40
FSYM s^
3.0 u V X y
4.0 X y
5.0 X y 9°
Data Structure
^^^COLS.
CARDS^-\^^ 1-10 11-13
1 NP
^^^COLS .
1-6 7-12 13 -18 19 - 25 26 - 34
CARD!s-\^^
^^\COLS. 11 - 20
CARDS^-v^^ 1-10
3 XL CPX
• •
NP + 2 XL CPX
numbers of about 1.3 and yaw angles around 60°. At large Mach
the location of the singular line about which the square root
have any sharp bumps. The mapped coordinates are printed so that
sive span stations from the leading to the trailing tip of the
nates. If the wing sections are all similar only the profile for
the first span station is needed as input. The coordinates for the
flow direction
at infinity
vortex sheet
on the other hand, the sections are not similar, the program permits
flow each x,y plane is divided into three strips. Then horizontal
lines are relaxed in the middle strip, marching towards the body,
and vertical lines are relaxed in each outer strip, marching out-
plane.
Normally calculations are first performed on a coarse mesh,
and then on a fine mesh with twice as many intervals in each coordi-
the starting guess for the fine mesh. This procedure greatly
reduces the computer time required for a fine mesh solution. Using
195
the CDC 6600 it takes one second to sweep through about 4500 mesh
points. The time for one iteration cycle on a mesh with 72*12x16
Tape 6. Tapes 1, 2 and 3 are disc files used for internal storage
the code which does not use disc storage is also available. This
cards listing, the required data items. The complete set of title
cards provides a list of all the data which must be supplied, and
values for the parameters listed. The input parameters are given
data cards. All data items are read in as floating point numbers
are converted inside the program. The data deck for Airfoil
similar only the chord and twist angle are printed at the remaining
size, Mach number, angle of yaw and angle of attack are also print-
ed so that the case can easily be identified. Then for each itera-
tion the program prints the iteration number, the maximum correction
the flow equation together with the coordinates of the points where
center section, the relaxation factors Rel Fct 1, Rel Fct 2 and
Rel Fct 3 (see Glossary, Section 4), and the number of supersonic
points.
convergence criterion has been satisfied the section lift, drag and
drag, friction drag and total drag, the ratios of lift to form drag
and lift to total drag, and the pitching, rolling and yawing
wing and the pressure distributions over the upper and lower sur-
faces separately, with the leading tip at the bottom of the picture.
is refined.
REL FCT 2 The supersonic relaxation factor for the
velocity potential. It is not greater than 1.
and is normally set to 1.
REL FCT 3 The relaxation factor for the circulation.
It is usually set to 1., but can be increased.
BETA The damping parameter controlling the amount
of added <(> (see Chapter I, Section 3) .
It is normally set between 0. and 0.25.
STRIP Determines the split between horizontal and
vertical line relaxation and is the propor-
tion of the total mesh in which horizontal
line relaxation is used. Fastest convergence
is usually obtained by setting STRIP = 1.,
where horizontal line relaxation is used for
the entire mesh. If convergence difficulties
are encountered STRIP may be reduced to some
fraction between 0. and 1.
FHALF Determines whether the mesh will be refined.
FHALF = 0.: The computation terminates after
completing the prescribed number of iteration
cycles or after convergence for the input mesh
size.
FHALF 7* 0.: The mesh spacing will be halved
after NRELAX cycles have been run on the crude
mesh size. An additional data card must be
provided for the refined mesh giving the
numerical values requested by Title Card 2.
If FHALF < 0 the interpolated potential will
be smoothed |FHALF| times.
TITLE CARD 3 (Aerodynamic Parameters)
FMACH The free stream Mach number.
4 NC
4.
6 ISYM NU NL
0. 82. 80.
TE TE
7 XSING YSING
ANGLE SLOPE
0. -.047 .0085 .0164
NU
9 X Y (Lower Surface Coordinates Nose to Tail)
1.
NL
10 Z CHORD THICK ALPHA NEWSEC
IF ([Link].100.) 60 TO Z5
CL HAS SEEN INPUTTED, KEEP IT FIXED
NCr a 0
MODE = 0
YA = ,5*CL/CHD-OPHI
00 11» L = 1,M
00 Itif J = ItNN
PHIIL.J) = PHI<L,J>+YA*PHJR<D
DPMI = .5*CL/CHD
CLX = CL
2S CU = CLX
c CHANGE PARAMETERS WHICH DEPEND ON THE MACH NUMBER
EM = AMAXl(EM,.lE-«tO)
IF ([Link]) NCY = 0
Cl = C2+1./(£M *£M )
C6 = C2»EM *£«
CU = 1.+C6
C5 = l./(C6*C7)
aCRIT = (Cl+Cl)/(6Afl«tt+l.»
BET = SQRTd.-EH *EM )-l,
C CHECK FOR TERMINATE,RETRItVE, OR STORe INSTRUCTIONS
C IK HILL 3E -1 ONLY IF THEKE IS A NAMEuIST £RROK
IF (([Link].O).OR.([Link].-D) 60 TO I/O
CALL COSI
IF ([Link].O) 60 TO «fO
REWIND N3
IF (ITYP.6T.O) GO TO 30
WRITE(NS) COMC,[Link],3B,ARCOLa,ANSOLn«[Link],D£LOLD,R,[Link]
I ,OSUM,6AMMA,XI«ION,RBCP,RFLO t ROEL,8CPiMSl,KPiST
GO TO 140
30 ftEAO (M3) COnC«PHI«AA«3dtARCOL^«ANGOLO«XOL^fYOLO,0£LOLO,R,RS,RI
1 ,DSUH»SAMMA,XflOM,rlBCP,RFLO,ROeLi8CPi,MSl,KP»ST
CALL HAP
60 TO ItO
ifO C O N T I N U E .
IF ( N S . 6 T . O ) GO TO 70
NS = 0
C 60 TO CRUDE GRID IF [Link].O
IF [Link].O) CALL RE*IESH(-1)
60 TO ItO
70 IF ([Link].O) 60 TO 100
C (JO BACK TO FINER GRID
CALL 3EMESHU)
60 TO 1*0
100 XPHII = 0.
IF < [Link].O.) XPHII = 2t«CHD/RCL
XA = 1.-3./RFLO
ANSO = -RAO*BB(1)
TXT = 3H CL
IF («OD£.E3.0) TXT = 3HALP
C NO BOUNOARY LAYER CORRECTIONS ARE MADE FOR [Link].O.
IF ([Link].O.) MSI = 1000000
IXX = «+2
80 IXX = IXX-1
IF (XC(lXX-l).[Link]) GO TO 80
204
LC = 0
00 AT .10ST NS CYCLES
00 120 K = 1,NS
IF (HOO(LC«56) .NIE.O) 60 TO 105 •
WRITE (N2.210) TXT
LC a LC+1
105 CALL SW&EP
KEEP TRACK OF TOTAL NUMBER OF CYCLES
NCY = NCY + 1
ALPX = KAD*ALP
CLX= 2t*OPHI*CHD
YA = YA*XPHII
WRITE RESIDUALS ON N2 EVERY KP CYCLES
IF (MOOtK,KP).NE.O) 60 TO 110
LC = LC + 1
WRITE (N2,190> N C l f t Y R , Y A « 0 1 « 0 2 « I K , J K » N S P » C L A ( 2 - M O O £ ) t A N 6 0 » C P 1 t
1 3CP,SL
00 A 30UNORY LAYER CORRECTION tVERY Nsl CYCLES
ItO IF (MOO<K,NS1).NE.O) GO TO 125
IF ([Link]) GO TO 140
WRITE (N2.190)
LC = LC+1
FSYfl = b.
CALL GrUR8<Ol,[Link],fieCP)
AN60 = -RAO*BBU)
IF ([Link].o) OPHI = .5*tLX/CHO
CHECK TO SEE IF HE HAVE SATISFIED CONVERGENCE CRITERIA
185 IF (AflAXKASSl YK) «ABS( YA) [Link]) SO TO 140
120 CONTINUfc.
1<*0 ITYP = IA8SJITYP)
CL = CLX
LN a: RiM*l.E-6+.5
XPF = XPF*AMINO(l,IABS(Hi»-N'»))
XP = XP*XPF
CALL SECONOITI1E)
MTPE = HH
TXT = 3HALP
IF ([Link].O) TXT = 3H CL
150 WRITE (NTPE.200) E«ltTXT,CLA(MOOE-H) ,[Link],NS,TI«E,RFLO,RCL»ROELt
1 [Link]»ST,PCHiSEP«,XStP,NPTS,IS,LLiI2
IF ([Link]) GO TO 160
NTPE = N2
GO TO ISO
160 IF ([Link].2) CALL GTORB(01,02,CPlfBcP»SL«RDEL,RBCP)
EMX = EM
ITYP=1
tiO TO 10
170 ITYP = 4
IF ([Link].-l) WRITE (N^,220)
CALL 5TURB<Oi,02,CPl«BCPiSL,ROEL,RBCP)
TERMINATE PLOT
CALL PLOT(0.,0.t999)
CALL EXIT
180 FORMAT <7H READ P/>
190 FORMAT ( 5Xi It,<»E12.3, If, 13,16.2F10.
205
SUBROUTINE LSTERR
COM10.M /A/
IK = -i
CMO
SUdROUriNE RESTRT
CO*1*IO^ PHI(lb2«31)tFP(162*31)tA(31)tBr3I),C(31)«0(31)te(3l)
1 ,RP(31),KPP(31),R(31),RS(31),KI(31),aA(162),Bb(162),CO(162)
2 •SI(162),PHIR(162).XC(162),YC(162)«Fv|(lb2).ARCL(162)iOSUH(1
3 ,AMGOLO(lt2).XOLO(1&2),YOLD(162).ARCnLt) <162),UELOuU(162)
COMMOM /AX [Link],RAOt£«[Link],RN«[Link]
1 ,XA,[Link],OELTH«OEL.K,RA»[Link]»RA'»[Link],C2
2 ,Ci».C5«C6»C7«3ET,3£TAtFSY»l,XSfc:H,SEP|«1. rTLE(«f)
3 ,XK.J<«[Link].I 1 10DE,IStNFCtNCT t N«N,N
<* , [Link],.1H
SET UP CONSTANTS
TP = PI+PI
RAO = 180./PI
AUP = ALP/RAO
IF ( dM + 1) .[Link]. («•*•!).[Link]) NCT = 0
MM = .1 + 1
IF (LL.E8.0) LU = V2 + 1
NM = IM + 1
OR = -l./FLOATIN)
OT = TPXFLOATlM)
OCiM = COS(OT)
OS.M = SIM(OT)
OEUR = .5/OK
OEUTH = .5/OT
RA = OT/OR
RA4 = 01*OT
00 10 « = ItN
«(K) = l.+OR*FLOAT(K-l)
RS(K) = (RA*R(K))*(RA*3(K))
RI(K) = -,25*DT/R«)
10 CONTINUE
R(MN) = 0.
SET = SQRT(1.-EM*E"!> -1.
206
00 MAPPING
CALL AZRFOL
IF ([Link].1) CL = 8.*PI*CHO*SI(l)/(i.+BET)
OPHI = ,5*CL/CHO
SELECT NT OF THE Ml MAPPEU COOKDINATEs
MA = MM/3
MB = <W-2*((HA+l)/2)
IF(([Link]).OR. ([Link].O. ) ) JK = -1
J=l
DO 40 L = [Link]
OELOLO(L) = 0.
OSUfl(J) = 0.
ARCOLD(L)=ARCL(J)
IF< [Link]»l) GO TO 70
[Link] ) . O R . < J . G E . M 8 ) ) J=J+1
O S U M ( J ) = 0.
J=J*1
i»0 COMTINUE
70 NT = L
WRITE (Nif.100) NT
ioo FORMAT (iHO,m,i»5H POINTS WILL BL usep TO oEFiNt INNEK AIRFOIL )
CALL SPLIF(MMtflRCL,XC,PHI(1.3).PHKl.R)«PHI(1.7),3.0.«3.0.)
CALL [Link], XOLO,ftRCL,XC,PnI(1.3),PHI (1, b) , PHI ( 1, 7))
CALL SPLIF(MM,[Link],PriI(l,3)«PHI(li5)«PHI(1«7).3.0..3.U. )
CALL I^[Link], YOLU,ARCL,[Link]!(1,3).PHI(1,5)«PHI(1,7))
CALL SPLlF([Link](lf3)»PHI(l.s>«PHI(1.7).3.0.,3fO.)
CALL [Link]([Link],A.>JGOLO,ARCL,FM,PHl(l,3),PMI(l,5).PHI(l,7»
DO 60 L = l.M
DO 50 J = [Link]
50 PHKL.J) = H(j)*co(D*-oPHi»pHiK(L)
60 CONTINUE
FSfM = FSYM-12.
IS = 2
RETURN
CNO
COSI
SET THE SINES, COSINES, AND THE TERM AT INFINITY
COMMON PHI(162«31>,FP(162«31),A(31), 8(31), C(31>, 0(31), E(3t)
1 , RP(31),RPP(31) ,R(31),RS(3l),rtn31),aA(162).Bb(162),Co(168)
2 , SI (162) ,PHIR(162) .XC(162) ,YC(1&2) ,Fv|( 162 ) , ARCH 162) »OSUM(162)
3 , ANGOLO<162) ,XOLO<162) .YOLOUfa2) . ARCnLD(162) , DELOLL)( X62 )
COMMON /A/ Pl,TP,RAD,EMiALP,RNtPCH,[Link],CHD,OPHI,CL,KCL,YK
2
3 , [Link] IZ«ITYP, MOOE«ISiNFC,NCr,NRlM»NS»IOlM«N2»N3»N'»,NT,IXX
<* , NPTS,LL,I,LSEP,fl4
TPI = 1,/TP
ANS = ALP4-8B(1)
SN = SIN(ftMS)
CN = SQ«T (l.
00 10 L = liM
207
com = CN
SKL) = SM
PHIR(L) =<ANG+ATAN((BET*SN*CN)/(l.+BET*SN*SN)))*TPI
CN = CN*OCN-SN*OSN
SN = CO(L)*OSN+SN*OCN
ANG = ANG+OT
10 CONTINUE
CO(MH) = CM
CO(MH+l) = C0(2>
SI(PlM) = SN
SHrtM + lJ = SI(2)
RETURN
ENO
SUBROUTINE SWEEP
SWEEP THROUGH THE GRID ONt TIME
C01KON PHK162»31),FP(162»31),A(31),Bj31),C(31)tD<31),E<31)
1 , KP(31) , RPPOD ,R(31) ,RS131),rtl<31),aA<162) ,BB( 162),CO(162)
2 , SI (162), PHI R( 162) , XC(162) , YC(162) «F,v<(162) . ARCL (1&2) i OSUMl 162)
62) ,XOLO(162) ,YOLD(162) ,ARCol-0(l&2) ,OELOLL)(162)
/A/ Pi;[Link]»£1«ALP,RN»[Link]«CHDiDPHl f i;UtKCL»YR
1 , XA, [Link]«DT,OK,OE:LTH,OeLH,RfltUCM,[Link].t;pSXL,QCrtir,Cl,C;2
2 ,C<ttC5«C6«C7.B£T.3tTA«FSY."l,XSt.P,StP«.TTLE('H [Link]
3 ,IK,JK,IZ,ITVP»HOOE,[Link],NRN,Ns.IOIMtN2tN3,N4 t NT,IXX
4 t [Link],."!'*
YR = Oi
NSP = 0
00 10 J = 1«NN
PHI(M*I.J) = PHKli J)+OPHI
PHKMVI+ltJ) = PHI(2,J)*OPHI
E(J) = 0.
10 RPP(J) = 0.
SWEEP THROUGH THE GRID FRO* NOSE TO TfllL ON UPPEK SURFACE
TE = -2.
00 30 I = [Link]
CALL
00 30 J
30 PHUI-ltJ) = PHKI-ltJ)-RP(J)
UPDATE PHI AT THE TAIL FROM UPPER SURFACE
00 50 J = liN
PHI<H«I»*J) = PHKHMtJ)-E(J)
E(J) = 0.
RPP(J) = 0.
50 PHIU.J) = PHI(MH.J)-OPHI
SWEEP THROUGH THE GRID FROfl NOSE TO T a IL ON LOWER SURFACE
TE = 2.
1 = LL
80 I = 1-1
CALL flURciAN
00 60 J = 1<N
60 PHI(I+liJ) = PHI(I*lfJ)-RP(J)
IF ([Link].a) 60 TO 80
208
00 70 J = l.N
70 PHK2.J) = PHI(2,J)-E<J)
ADJUST CIRCULATION TO SATISFY THE KUTrA COMOITION
IF (RCL .Ed.O.) 50 TO 90
YA = RCL*<<HHI(«,l)-tPHl(2,l)+OPHI))*o£.LTH+Sl(l))
IF (MOJfc.E».l) GO TO 90
ALP = ALP-.b*YA
CALL COSI
60 TO 9b
90 YA = TP*YA/<1.+BET)
OPHI = UPHl+YA
9S 00 97 L = lt«
9? PHI([Link]) = OPHl*PHIR(L)
IF([Link].O) RETURN
00 100 J = l.N
00 100 L = 1,«
100 PHI(LtJ) = PHICL, J)+YA*PHIrUL)
RETURN
ENO
SUdROUTINE HURilAN
c SET UP COEFFICIENT ARRAYS FOR THE TRIQIAGONAL SYSTEM USED FOR LINE
C RELAXATION AND COMPUTE THt UPOATEO PHT ON THIS LINE
1 ,HP(31) ,RPP(31),R(31),RS(31),KI(31),aA<162) ,Bb(162>,(-0(162)
2 .SI(l62),PHIK(162) iXC(162) • YC ( 162 ) t F«| ( 162 ) . ARCH Ib2) tOSUM
i> tAMGOLU(!62) , XOLD < 162 > , YOLD(lb2 > , ARCnLO ( 1 &2 > , OELOLL) ( 162 >
CO^HOM /A/ Pl,[Link]«[Link],[Link],KCLtYR
1 ,XAtYAtTE<OT,OK,L)ELTHtOEL.-<iRA«lJCNtUS.M«RAl»«EPSlL<QCKir,CltC2
2 , C < f , c 5 t C 6 « C 7 « B E T « 3 E T A . F S Y , < 1 , x S E P , S E P n . T T L E ( < f ) .[Link] NSP
3 , I K , J r < « I Z » I T Y P . « O O E i I S » N F C . N C T t N R N , N B , IDIMt [Link] t N T « IXX
"* • [Link]
DO THE BOUNUARY
E(NN) = 0.
FAC = -,b*TE
IH = 1-1
IF ( F A C . L T . O . ) I* = 1+1
KK = 0
PHIO = P H I ( I , 2 ) - 2 . * O R * C O ( I )
PHlYP= P H I < I i 2 ) - P H K I t l )
PHIYY = PHIYP+PHIO-PHKI,!)
PHIXX = PHl(I + lil)+PHI(I-lil)-PHI(I«U-PH
PHIXM = PHKItltD-PHKI-lfl)
PHIXP = PHKI + lf 2 ) - P H I ( I - l t 2 )
CHECK FOR THE TAIL POIMT
IF (I..ME..1M) SO TO 10
= (Cl+Cl)*RS(l>
= -C(U+XA*C1-C1
0(1) = C l « ( P H I X X * R S t l ) * P H I Y Y - l - R A i » * C O ( I ) - E ( l ) )
SO TO <»0
10 U = PHIXVUDELTH-SIU)
BQ = J/FP(I«1)
209
QS = U*BO
CS = C1-C2*QS
BQ s BQ*QS*(FP(I-1»1)-FP(I+1,1»
X = RAt*(CS+QS)*CO(I)
C(l) = »CS+CS)«RSU)
0(1) = CS*RS(1)*PHIYY+RI(1)*BO+X
CMOS = CS-QS
PHIXT = 8ETA*ABS(U)+ABS(C«QS)
IF ([Link]) GO TO 30
FLOW IS SUPERSONIC, BACKWARD DIFFERENCES
KK = 1
PHIXT = PHIXT-CMQS
PHIXX* = RPP(l)
AU) = -<C(1)+PHIXT)
0(1) = Om+CMOS*PHIXX1-PHIXT*E<l)
GO TO 40
FLOW SJBCRITICAL, CENTRAL DIFFERENCES
30 Ad) = XA*CMQS -CU>-PHIXT
0(1) = D(1)+CM9S*PHIXX-PHIXT*E(1)
00 NON-bOuNOARr POINTS
*0 RPP(l) = PHIXX
DO 60 J = 2,N
PHIXX = PHHI+liJ)+PHI(I-l,J)-PHI(I,Jt-PHI(I,J)
OU = PHIXP
PHIXP = PHI(I+1,J+1)-PHI(I-1,J-H)
PHIXY = PHIXP"PHIX1*(E(J*-1)-E<J-1»*FAC
PHIXM = OU
OU = OU*OELTH
PHIYYI"! = PHIYY
PHIYM = PHIYP
PHIYP = PHI(I«J+1)-PHI(I,J)
PHIYY = PHIYP-PHIY.I
U = R(J)*DU-SI(I)
OV = R(J»*(PHI(I,J+1)-PHI(I, J-1))*OELR
V = OV*«(J)-CO(I)
RAV = R(J)*RA*V
BQ = l./FPd.J)
BQU = aa*u
US = 88U*u
UV = ([Link])«V
VS = 33»V*V
OS = US+VS
CS = C1-C2*QS
CMV/S = CS-VS
CMJS = CS-US
PHIXT = BETA*ABS(U)
PHIYT = BETA*A3S(RAV)
COMPUTE CONTRIBUTION OF RIGHT-rtANO SIQE FROM LOW OROEK TERMS
0( J) =RA<t*<
UV = .5*8QU*RAV
IF ( Q S . L E . Q C R I T ) GO TO 50
SUPERSONIC FLOW, USE B A C K W A R D OIFFEREistClNG
KK = KK+1
CMQS = CS-QS
210
FQ = l./OS
AUU = US*FQ
BUJ = rtS(J)*AUU
SW = VS*FQ
AW = RS(J)*BW
BUV = JV*FS
AUV = 39U*A8S(RAV)*FO*TE
PHINN = BVV*PHlXX-BUV*PHIXY4-BUU*PHIYr
B(J) = CS*BUU
PHIXT = PHIXT-CMQSKAUU+AUU-AUV) +CS*qVV
PHIYT = PHITT -CMQS*(At/V+AVV-AUtf)
C(J) = 8(J)+PHIYT
PHIXXd = «PP{J)
IF ([Link].o) GO TO US
PHIYYI"! = PHI(I.J + 2)-PHI(ItJ+l>-PHlYP
PHIXY1 = PHIYP+PHKIM.J)-PNl(If1iJ*l)
60 TO f6
i|5 PHIXYW = PHI(I«,J)-PHUI«tJ-l)-PHIYH
89 = 8(J)
B(J) = C(J)
C ( J ) = bQ
if6 PHISS = AUU*PHlXXM*AU»/*PHJXYMtAW*PHIrYM
A ( J ) = -O< J)+C{J)*PHIXT)
D ( J ) = D(J)+C«aS*PHISS4-CS*PHINN-E(J)*pHIXT
GO TO 60
SUBSONIC FLOW, USE CENTRAL DIFFERENCES
50 C(J) = RStJ)*C«VS
B(J) = C(J)*PHIYT
PHIXT = PHIXT+CMUS
A(J) = XA*C«US-8(J)-C(J)-PH!Xr
0(J) = U(J)+CMUS*PHlXX-Utf*PHIXY+C<J)*PHIYY-PHIXT*E(J)
IF ( V . L T . D . ) GO TO 60
8(J) = C(J)
C(Ji = CUJ+PHIYT
60 RPP(J) = PHIXX
MSP = IMSP-fKK
SOLVE THE TRIOIAGONAL SYSTEM
CALL T«10
RETURN
SUBROUTINE TRIO
SOLVE H OIMtNSIONAL TRIDIAGONAL SYSTEM OF EQUATIONS
COMMOM PHI(lfe2t31),FP(162«31)tA(31»,B(31),C(31),D(31)tE(31)
1 ,RP(31I,RPP(31),R(31),RS(3l),KI(31>,ftA(162),B«(162),CO(16z>
2 ,SI(162),PHIR(162)«XC(162),YC<162)«FM<162>«ARCH 162)»DSUM( 162)
i> ,AIMGOLU(l62) , XOLO (162 ) , TULD (1&2 ) ,ARCoLD(162) ,OELOLO(162)
COMMON /A/ Pl,TP,RftO,[Link],RN,PCH,XP,TC,CHO,DPHI,CL,[Link]
1 ,XA, YA«TE,OT,OR,OELTH,OELK,RAtOCN,US,M«RAt>,EPSlL,QCKIT,Cl,C2
2 ,C4,C5.C6«C7,8£.T,3tTA,[Link]»l«XStP,[Link](t) ,H,N.n,1,NN,NSP
3 ,IK,J*.IZ,ITYP«. w IOO£.IS,NFCtNCY,NRN,[Link],N2«N3»N l VtNT»IXX
1
•* , [Link],."! *
211
XX = l./AU)
RP(1) = E(l>
Ed) = XX«0(1)
00 ELIMINATION
00 10 J = 2«N
C(J-l) = C(J-1)«XX
XX s-. [Link](J>-B(J)*C(J-1))
RPtJ) = E(J)
10 E(J) = (0(J>-B(J)*E(J-1))*XX
00 BACK SUBSTITUTION
E«X = ABS(E(NM
00 20 J = 2«N
L = NM-J
E(L> = E(L»-C(L)*E(L+1)
20 EMX = AHAXl(EMX,ABS(E(L))>
FIMO THE; LOCATION OF THE MAXIMUM RESIDUAL
IF (E*[Link].A8S(YR)) KETURN
IK = I
DO 70 J = 1«N
IF (ABS(t(J)).EQ.E*IX) SO TO 74
70 CONTINUE
7* JK = J
YR = E(JK)
RETURM
END
SUBROUTINE REMESH(LSISM)
C 60 TO CKUDER GRID IF LSIGN IS -1
C SO TO FINER BRIO IF LSIGN IS +1
C01MOM PHI(162,31),FP(162.3l),A(31),B(31),C(31),0(31)<E(31)
1 ,RP(31).RPP(31),R(31),RS(31),RI(31),AA(162),BB(162),CO(162)
2 , SI (162) ,PHIR(162) • XC < 162 ) , YC < 162 ) . F"l < 162 ) , ARCLtlbi!) .DSUM(162)
4 , ANGOLO(l&2),XOLO(162),YOLU(16a),ARCoUO(lb2),UELOLL)(16i!)
COMMON /A/ [Link],RAO,E*[Link],[Link] t [Link]»CHO,OKHliCLtRCL»YR
1 , XA,YA,[Link],OK,OELTH,DELr<,[Link],DS^,KA'*,EPSlL,QCKll,Cl,C2
2 ,Cf ,C5»C6tC7»3ET!8t;TAtFSY?1,XSt;p,SEPM,TTLE(<») ,P1,NiMntNN,NSP
i , IKtJK, IZ«lTYP»«OJEiIStN(-CtNCY«NRN,N(;,IOIMtN2«N3«N<f,[Link]
1 , NPTS,[Link]£P,,1H
X = 2.**LSieN
NG = FLOAT(NG)/X+.2
H = FLOAT(«)*X +.2
M = FLOAT(N)*X*.2
LL = FtOAT(LL-l)*X*1.2
IF (LSI6N.6T.O) MM = H+l
IF ([Link].O) NN = Ntl
LSEP = FtOAT(LSEP-l)»X'H.2
PF = l./X
OELR = X*OELR
OELTH = X*0£LTH
OR = PF*OR
OT = PF*OT
OCN = COS(OT)
212
OSN = SlN(OT)
NCY = 0
I = LSI6N
HP = MM* i
CALL PEKMUT (R,NN,1)
CALL PEKMjT ([Link].l)
DO 5 J = itN
5 ftl(J) = -,25»OT/R(J)
CAUL PERMUT (OSUM,1P,1)
00 20 L = ItNN
20 CALL PEK«L|T ( PHI < 1 . L ) , IP, 1 )
00 30 L = liMP
30 CALL PERMUT (PHI (L, 1 ) ,MN, IOIM)
m = «n-x
IF ([Link]..5) 60 TO 80
00 HO L = 1.M.2
DSU«(L*D = . 5 * ( O S U H ( L ) * D S U M ( L * 2 ) )
00 40 J = ltNM,2
40 PHKL+lfJ) = .5*(PHI(L,U)*PHI(L-t-2.J) )
00 50 J =
DO 50 L =
50 PHKLtJ+1) = .5*(PHHL(J)+PHI(L,J-»-2))
60 CALL «AP
RETURM
END
20 L = L+l
60 TO SO
30 00 10 J = 1«NY.2
A(J) = AX(U>.
HO L = L + JX
00 50 J = [Link].2
A(J) = A X < L )
50 L = L+JX
60 L = 1
00 70 J = [Link]
AXIL) = ft(U)
70 L - L+JX
RETURN
CNO
SUBROUTINE GETCP(COF)
COMPUTE [Link] AND CM 3Y INTEGRATION AND OUTPUT MACrt DIAGRAM
COMMON PHI(162.31),FP(162»31),A<31),8(31),C(31).0(31)«E(31)
1 ,RP(51),RPP(31),R(31),RS(3l),KI(31),aA(162),B8(162),CO(16?)
2 ,SI(162),PHIR(162),XC(162).YC(162)iFM(162),ARCL(162).DSUM(162)
3 ,ANGOLD(l62),XOLD(162).yOLO(l&2),ARCr)UD(162),OELOLO(162)
COMMON /A/ PI,rP«RAO.£1«ALP,RN»PCH,[Link],[Link],Cl.,[Link]
1 ( XA,rA«TE<OT,OR t OELTHiOELR,RA«OCN|OSN«RA4«EPSIL«QCRIT,Cl,C2
2 ,C"*,C5,C6«C7.3ET,3ETA,FSTV|,XSEP,[Link]('*)««.NiM.«ltNN,NSP
3 , IK, jKiIZ,ITYP,[Link]«[Link],IMRN,(g S ,IOIM,N2»N3.N<*,NT,IXX
H , [Link].I.LSEPil1*
REAL «ACH,MtMACH
COMPLEX [Link]
DIMENSION «ACHM(l)-iCPX(l)tMN(l),IHACH»21)
EQJIVALtMcE (MACHN(l)«A(1>),(CPXI1),PnIR(l)),(f1N<1)«FP(1«31))
DATA IMACH/lHQ,lHR,lHs,lHTtlHu»[Link],lHX.lHYilHZflH0.1Hl,lH2>lH3
l,lrtit,lHS,lH6,lH7flH8,lH9,lH+/
DATA TX /i*HCDF=/
MACH(O) = SQRT(Q/'(C1-C2*Q))
inC(Q) = MINO([Link].«8)-H>
CLCO = 0.
CM = Ot
IF (([Link].O.).OR.([Link])) SO TO In
DY = YOLO(NT).YOLOd)
REMINO HI*
WRITE (flt,120) [Link] NRN,nn
10 00 20 L = ItMM
CP = CPX(L>
COMPUTE CP*DZ
TMP = CP*SQRT(FP(Ltl> )*CnPLX<COS<FM<|_» .SIN(FPKL) M
SUI UP [Link], AND CM
CLCO = CLCO+TMP
CKI = C«*(XC(L)-.25)*REAL(T«P).YC<L)*AT«AG<TMP)
WRITE PUNCH OUTPUT ON flf IF XP=0 AND [Link].80
IF «[Link]>0.).OR.([Link].80» GO TO 2n
Q = MACHN(L)*SQRT(Cl/(l.+C2*MACHN(L)niACHN(L))»
V = Q*SIN(FM(L»>
214
U = Q*COS(FM(L))
IF ([Link].O) GO TO 15
WRITE (nit,130) U,V,XC<L)«YOLO(L),CP
GO TO 20
15 WRITE (Pl<t,130> U i V t X C < U ) t Y C U > t C P
20 CONTINUE
CORRECT CLtCO FOR ANGLE OF ATTACK
CLCO = -(0 T*CHO)*CLCD*C«PLX(SIN(ALP).COS<ALP»
CM = OT*CHO«CH
WRITE CD,CL,CM ONTO NH
COd = REftt(CLCO)
CO = COW+COF
CL3 = AI«IAG(CLCD)
IF ([Link]) GO TO 85
IF ([Link].O.) GO TO 70
WRITE (Nit,90) E!*l,CL2»C*l,CDWtTX»CDF, Co
GO TO 30
70 WRITE (Ntt,90) EM. CL2i t«l, CUW
CONSTRUCT MACH NUM3CR DIAGRAM
WRITE (Nlt,l<to)
80 I =
I =
USE PRINT WIDTH OF 1Z FOR 1ACH NUMBER DIAGRAM
HB = .lil
HC = MAXO(1,MB/IZ)
HA = 1C+!1AXO(1,M9-IZ*MC)
WRITE OUT MACH NUMBERS AT INFINITY
WRITE (Nlf, 100) (Ii L = HA,M8,HC)
00 «ACH NUMBERS ONE LIME AT A TIME OOyN TO THE BODY
J = NN-«C
l»C KSJ = «(J)*R(J)
00.50 L = MA,MB,1C
U = (PHI(L+l,J)-PHI(L-liJ))*R<J)*0£LTH-SI(L)
V = (PHI<ttJ+l)-PHI<L.J-l»*DEt-R*RSJ -CO(L)
0 = (U*U+,/*V)/FP(Lt J)
1 = I1C««flCH(Q»
MN(L) = I«ACH(I)
50 CONTINUE
WRITE (Nit,100) (MN(L)tt = ."IA,MB,MO
J =• J-flC
IF ([Link].t) GO TO tO
DO THE LINE WHICH IS THE BODY
00 60 L = MA,M3,>1C
1 = InCfKACHNIL) )
60 MN(L) = IMACH(I)
WRITE (Nij,100) (MN(L)tL = [Link])
IF ([Link]) CALL GRAFIC(CD)
RETURN
65 RNX = .l*aINT(RN»l.E-5)
WRITE (Nit,ISO) [Link]
RETURM
90 FORMAT (lHl2XSHE!*l=[Link],ifX3HCL=[Link]«i»XsHCH=F6.i»i'*X«»HCOW=F7.5.ifXA4
1 fF7.5»'*x 3HCO=F7.5///)
100 FORMAT <3X«130A1)
180 FORMAT (3H M=t Fif.3.5Xt 3HCU=tF5.3t5Xt3HOY=fFi*.3«6X.»HT/C = «
215
1 Ft.3tl*X«2I5)
130 FORMAT U020)
1(»0 FORMAT (1HO//)
150 FORMAT (lHO//7X3HEM=»F4,3,<»X3HCL=,F&.<»,<tX*HT/C=»F»,3,»X3HCI"l=«
1 F6.
ENO
SUBROUTINE SRAFIC(CD)
COMPLEX ZP,ZG,SFAC,SIG
REAL MACHN
COMMON PHIU&2»31),FP<162«31»,A<31),B<31),C<31>,0<31),E<31)
1 ,RP(31),RPP(31),R<31>iRS(31),KI<31),4A<162),8B(162),Co(162)
2 ,SK162),PHIR(162)«XC(162),YC(162>,FM 162)»ARCLU62)»DSUM(16a)
3 , A N G O L D ( 1 6 2 ) ,X O L O U 6 2 ) , Y O U O ( l 6 2 ) , A R C O U O ( 1 6 2 ) , O E L O L O < 1 6 2 )
CO««OM /A/ PI«TP,RADt£1 1 Al.P,RN«PCH,[Link].«yR
1 , X A . y A t TE«OT,OR,OELTHtDEURtRA«[Link]»EPSIl..QC«IT,CltC2
2 tC4.C5iC6tC7»8£T«8tTA,FST!«l,XStP,SEP«.TTLE(«>).[Link],NN,NSP
6 , IK,jK,IZ,ITYP,MOOE,IS,NFCiNCr,NRNiNr,.IOIM,N2.N3,Ni+,NT,IXX
» , [Link].I*
DIMENSION CPX(1).MACHN(1).T(6)
EQUIVAUEIMCE (CPX(l) tPHIR(l) ) , (MACHNt 1 ) , A(1) )
DATA TOL/l.E-6/ , PF/-.4/ , SCF/5.0/,yOR/it.O/tSIZE/.l't/tSCO/ZOO./
MOVE THE ORIGIN TWO INCHES OVER AND TjO INCHES UP
CALL PLOT(2.0«2.5i-3)
YOR = AflAXl(3.5,.5*AINT(20.*EM-7.0»
PLOT CP CURVE AS A FUNCTIOM OF X
CPF = l t /PF
CCP = CPF*CPX(1)
CALL PLOT(SCF*XC(1),YOR+CCP,3>
00 10 L = 2,MM
CCP = A«IMl(8.5-YORtCPF*CHX(L»)
10 CALL PLOT(SCF*XC(L>iYOR+CCP,2)
OR A nl AMO LABEL THE CP-AXIS
CALL cPAXlS(-.5fYOR«l.-l./PF,7.5-YOR,pFi
COMPUTE AND PLOT CRITICAL SPEED
CALL SYMBOL ( -.5, YOR+CPF*CPX(NI"H-1) 12.*SIZE«15«0. ,-1)
PLOT BODY
CALL PLOT(SCF*XC<l).SCF*Yt(i),3)
00 20 L = [Link]
30 CALL PLOT(SCF*XC<L)»SCF*YC(L).2)
LABEL THE PLOT
ALPX = RAO*ALP
TXr=8HANftLYSIS
[Link].b.) TXT=6HTHEORY
XL=-,9
**«*NON-ANSI - SEE VOLU1E I, PAGE 209****
IFtFSY^.GE.b.) GO TO 30
ENCODEI60,19l,T) TTLE,[Link]
GO TO HO
30 LN=RN*I>E-6+.5
ENCODE(60,190. T ) T T L E t I t N t f J C Y i L N
40 CALL SYMBOL<-[Link],-1.0iSIZE,T.O.,56)
216
00 90 J = 1<N
ZP = SFAC
RJ = R(J)
OS = 3
IF ([Link].l) GO TO 82
u = (pHUL+i.J)-PHi<L-it J) »*RJ*OELTH-sx
v = <PHI<L«J+i)-PHi<LiJ-i>)*DELR*RJ*RJ-CX
Q = (U*U + V/*V)/FP(L»J)
Q = SaRT(Q/(Cl-C2*a) )
82 SIS = C<"IPLX(RJ*CO(L>,RJ*SI(L))
C COMPUTE «1-SIGMA)**U-£PSIL))SIGMA
SFAC = CEXP(EX*CLOS«[Link].)-SIGM/SIG
c SUM UP FOURIER SERIES TO OBTAIN CONJUGATE OF w
S = -36(1)
00 B<* < = liNFC
LT = MOO«L-l)*KtM)
S = S+RJ*( AA<K+1)*SHLT«-1)-BB<K*1)*CO<LT+1) )
RJ = *J*R(J)
IF ([Link]) GO TO 65
R<* CONTINUE
C COMPUTE THE ARGUHEiMT OF OZ/DR
Sb SFAC s -SFAC*CNIPLX(COS(S)iSIN(S))/CABS(SFAC)
C MULTIPLY THE ARGUMENT 3T THE MAGNITUDE TO OBTAIN OZ/OR
SFAC = SFAC*(CHD*S9RT(FP(L|J)))>(R(J)*R(J))
C PERFORM THE INTEGRATION
23 s Z9+FAC*SF4C
FAC = OK
IF ([Link].l.) GO TO 100
90 CONTINUE
100 ZQ = za-.5*OR*SFAC
ZP = za-.5*OR*(SFAC+ZP)
Rl = O-1.)/<Q-QS>
ZP = za+Ri*(ZP-za)
CALL PLOT <SCF*N£ALIZP),AHAXK-2.0iSCF*AIMAG<ZPM»2)
GO TO 120
110 IPEN s. 2 •
IF (MACHN(L-l).GE.l.) SO TO 70
120 CONTIMUt
C POSITION PEN AT BEGINNING OF NtXT PAGE
122 CALL PLOT UO.O,-2.b,-3)
IF {(FSY«[Link].7.).OR.([Link].6» RETURN
C PLOT THE BOUNDARY LAYE* DISPLACEMENT
MX = INOEXR ([Link])
CALL PLOT(2.,1.5,-3)
CALL SYM80L(1.36,-.[Link]»19HLOWER SuRFACE DELS <0.il9)
CALL CPAXIS ([Link].«[Link].,l./SCO)
c PLOT LOWER SURFACE
CALL PLOT (SCF*XC(1),SCO*USUH(1),3)
00 132 L s 2,MX
132 CALL PLOT (SCF*XC{L)iSCO*OSUM<L),2)
CALL PLOT(0.i4,5»-3)
CALL SrMBOL(1.36,-,65tSIZ£,19HUPPEH SijRFACE DELS tO.«19)
CALL CPAXIS (0.,O..[Link].,l./SCU)
C PLOI UPPER SURFACE;
CALL PLOT (SCF*XC(>IX>,SCO*DSUM(nX)«3)
218
00 I3f L = MX,M
CALL PLOT CSCF*XC(L+l)fSCU*OSUfl(L+l)o)
CALL PLOT(10.»-6.t-3)
RETURN
FORMAT <10X,I3)
150 FORMAT <3F6.3,F7.5,E9.1)
160 FORMAT «2F10.«H
17U FORMAT <Al2,«»H M=Ft.3,3XHHALP=F5.2t3x3HCL=»F5.3,3X3HCD=iF5.»)
190 FORMAT<tA!f,3XHH«*N=I3,lH*l2,3XtHNCy=IiMtX2HR = I2,BH MILLION)
191 FORMAT<tAt,3X«fHM*N=I3,lH*I2«3XtHNCY = Iln<*X12HNO VISCOSITY)
END
SUBROUTINE CPAXISJXOK,YOR«60T«TOP,SCFi
C DRAWS ANO LABELS THE CP AXIS
C [Link] IS THE LOCATION Oh THE ORIGIN OF THE AXIS
C BOT IS THE LENGTH OF THE AXIS BELOM TH£ ORIGIN
C SCF IS « SCALE FACTOR JSEO FOR LA3LLIMG
C SCF NEGATIVE FOR CP AXIS ANO POSITIVE FOR DELS AXIS
SIZE - .12-SIGN(.02,SCF)
C DRAW THt VERTICAL AXIS
CALL PLOT (XOR,YUR*TOP,31
CALL PLOT «XOR,YOR-BOT,2)
C DRAW HATCH MARKS AMU LABELS ONE INCH fiPART
N = 1 + INT(BOT)-HNT(TOP)
S = -AINT(BOT)*SCF •fl.E-12
XH = XOR-(3.*SIZE)/t7
YH = YOR-AINT(BOT)
00 10 I = ItN
CALL SYMBOL (XOR,YH,SIZE,15,0.t-1)
C ****NON-ANSI - SEE VOLUME I, PAGE 209****
IF ([Link].O.) ENCODE (10,25,A) S
IF (SCF.|_E.O. ) ENCODE (10,20, A) S
S = S+SCF
CALL SYMBOL <XH,[Link],A,0.,t)
10 YH = YH*1.
IF ([Link].O.) GO TO 30
CALL SYMBOL(XOR+.[Link]+2.5,. 1I*,1HC,0.,1)
CALL SYM80L(XOR+.25tYOR+2.3B..[Link],0.,l)
RETURN
C ORAW TH£ X-AXIS
30 CALL PLOT (XOR,[Link].3)
CALL PLOT (XOR*5.0,YOR-BOT,2)
CALL SYMBOL (XOR+5.5,TOR-.07,.It,1HX,o.«1)
YH = YOR-BOT-SIZE-SIZE
00 40 I = 1.5
S = ,2*PLOAT(I)
ENCODE <10«20«A) S
XH = YOR+FLOAT(I)-SIZE-SIZE
CALL SYMBOL (XH,[Link],A,0.,4)
i»0 CALL SYMBOL (XOR*FLOAT(I),YOR-dOT,SIZr,15.90. .-1)
CALL SYMBOL ( XOR+. 25, YOH*.3.0, . It, tHDELS, 0 . , t)
RETURN
219
25 FORMAT ( Ft.3)
80 FORMAT (F4.1)
END
I S H I F T ( X X X . Y Y Y ) = SHIFft X X X t Y Y Y )
N = MOOdABS(NRN)flOOO)
CALL REAOCP (I0,ai3il)
10(1) = I S H I F T ( l D ( 2 ) . A M O . M S f - 1 8 )
00 10 L = [Link].2
J = L/2+1
IF (LTAB(j)-IO(l)) IQtZOtlO
10 CONTIMJE
L = NU*-i
20 ENCODE (60.30.10) MAHE(L) i MAHECL+1) «N
IF (NRN.5T.1000) SO TO 50
CALL PLOTS (600,10)
KETURM
50 CALL PLOTSBL (600,10)
30 FORI«IAT(A10t5H — ,[Link].I3)
ENO
SUBROUTINE AIRFOL
C REAOS IN DATA FOR AIRFOIL AND DETERMINES TqE MAPPING
C FUNCTION 8Y COMPUTING FOURIER COEFFICIENTS
C IF ONLT X.Y COORDINATES A«£ PRESCSISEo SLOPES ARE COMPUTED
COMMON P H I ( 1 6 2 , 3 1 ) , F P ( 1 6 2 « 3 1 ) , A ( 3 1 ) , B ( 3 1 ) , C ( 3 1 ) , 0 ( 3 1 ) . E ( 3 1 )
1 ,RP(51).RPP(31),R(31),RS(31),RI(31),aA(162),BB(162),CO(162)
2 ,SI(162),PHIR(162),XC(162),YC(162),FM(162),ARCL(lfa2),DSUM(162)
3 .AN60LO(l62),XOLD(162),YOLO(l62),ARCoLO(162),OELOLO(162)
COMMON /A/ Pl,TP,RAO,E?[Link]»PCH,[Link],CHO,OPHI,CL,RCL«rR
1 ,XA,rA,TE,OT,Oft,DELTH,OELK,[Link]'RAl*.ePSIL,QCKIT,Cl.C2
2 ,C4,C3,[Link],FST»1,XSEP,SEP|«|,TTLE('f) ,n,N,P1«,NN,NSP
3 ,IK,JK»IZ,ITYP,HOOE,IS»NFC,NCY,NRN,NR,IDIM,N2,N3,N4,NT,IXX
•* . iMPTS«LL,Ifl-SEPt«lt
220
DIMENSION XX(l),YYU>.U(l).V(l)tH<l)tSP<l)»CIRCU)iTH<l)tTT(l)
1 , OSd).SS(l) tCX(l),SX(l),9SR(l) ,TITLt(15).2(l)
EQUlVALtNCEUXdl «FH(l,3 )),(YY(l),FP(1.5))%(U(l)iFP(l.D),
2 ( 1,13) ) i (TH( l)iFP< 1,15 )).<TTU)iFP< 1.17 ) ) • < DS( 1 > iFP( 1,19) ) .
3 (SS(1).FP(1.21) ) ,(CX(1),KP (1,23) )f(Sy(l).FP (1.25) It (USR(l),
H FPd.27) )»(Z<1),FP<1,29)>
53(82) = Q2*Q2
SMOOTH(«l,aa,Q3,Ot) = a2+Sa(SQ(SQ(Qt) i )*.25*<Ql-Q2-82+Q3)
OIS(Ql) = {81-ERR)*< <01-EKR)*(ai-ERRUCONST)
D A T A TOL, NT. ISYM, CONST, VAL/.HE-7.999«n«.2»'*HRUN /
DATA OXDsitOXOs2,OlTDsliOYDs2/it*0./ . KT/-1./
N*1P IS THE NUMBER OF POINTS IN CIRCLE PLANE FOR FOURI£R SERIES
LC = NFC
^MP = 2*LC
MC = NflP + 1
PILC = Pl/FLOAT(LC)
IF ([Link].6.) 60 TO 150
WRITE <Nij,it70)
REWINO N3
REAO (N3,i|10) TITLE
IF ([Link].3.) GO TO 100
READ in COORDINATES AS PROOUCEO BY PROGRAMS o AND F
EPSIL = 2.
XX(1) = 0.
NL = 2
REMIND N3
READ (PJ6.510) E M , C n O Y , T C » M R N
IMC = « 0 0 < I N T ( 1 0 0 . * £ l " H - . 5 ) U O O )
ICL1 = M O D ( I N T ( C L + . 0 5 ) f l O )
ICL2 = M O O ( I N T ( 1 0 . * C L + . 5 J i l O )
ITC1 = rtOOlINTUO.*TC + . 0 5 > » 1 0 >
ITC2 = M O O < I N T ( 1 0 0 . * T C * . 5 ) t l O )
ENCODE U 0 . 5 3 0 . T T L E ) I«IC, 1CL1.ICL2. I T c l » I T C 2
MODE = 0
IF ([Link].O) FSYM=2.
10 REAO (M3,500) U ( 2 ) , V ( 2 ) , XX ( 2 ) , YY ( 2 ) , FAC
IF (XX(ii).LT.l.) 60 TO 20
SAVE TAIL POINT ON LOWER SURFACE
= U(2)
XX(1) = XX(2)
YY(1) = YY(2)
GO TO 10
20 00 HO L = 3t999
READ (M3.500) U(L) i V(L) iXX(L) t YY (L) <F&C
****CHECK FOR END OF FILE****
IF (EOF(NS).NE.O) SO TO SO
IF (XX(L).EQ.l.) 63 TO 70
IF (XX(L).[Link](NL)) !XL = L
HO CONTINUE
AIRFOIL HAS BEEN EXTENDED IN PROGRAM 0
50 XT = It
70 NT = L
IF (XX(l).EQ.l.) GO TO 95
221
IF ([Link].O.) XT = l.+.6*OY
NRN = lABS(NRN)
C INTERPOLATE TO PUT THE TAIL AT X=XT
C LOWER SURFACE INTERPOLATION
1 = 1
L = 2
80 Rl = (XT-xX(|_+l))/<XX<L)-XX<L+lM
R2 = lt-Rl
rr(I) = R1*YY(L)+«2*YY(L+1>
U(I) = Rl*U(L)+R2*U(L+l)
V(I) = R1*V<L)+R2»V(L«-1)
XX(I) = XT
IF ([Link]) GO TO 150
C UPPER SURFACE INTERPOLATION
I = NT
L = NT-2
GO TO 80
C READ IN AIRFOIL DATA FROM CARDS
100 READ <N3,<»20) FNUtFNLtEPSlL
READ <N3,<»70)
NT = FNU+FNL-1.
NL = FNL
00 110 I = NLtNT
110 READ <N3,if20) U(I)tV(I)tXX(I),YY(X)
READ (iM3,!f70)
00 120 I = [Link]
J = NL+l-I
120 READ (N3,t»20) U(J)iV(J)tXX(J),YY(J)
00 130 J = l.H
ISO TTL£(J) = TITLE(J)
IF ([Link].) GO TO 150
00 ItO L = 1,NT
TH(L) = XX(L)/RAO
XX(L) = U(L)
mO YY(L) = tf(L)
GO TO 195
C HO PERIOD IN THE STREA1 FUNCTION
95 EPSIL = 0.
C DEFINE SLOPES SO THAT ARC LENGTHS CAN BE COMPUTED TO FIRST ORDER
150 IF (([Link].S.)) SO TO 170
00 160 I = [Link]
160 TH(I) - 0.
ISYH = 1
GO TO 20Q
C COMPUTE SLOPES FROfl VELOCITIES
170 TH(1) = A T A N C V U ) / U ( 1 M
QSR(l) = U ( 1 ) * U ( 1 > + V ( 1 ) * V < 1 )
DO 190 I = 2,NT
c CHOOSE NEAREST BRANCH FOR THE ARCTANGENT
OTrt = ATAN((U(I-1)*V(I)-U(I)*V(I"1))/(U{I-1)*U(I)+V{I-1)*V(I)))
TH(I) = TH(I-1)+OTH
190 QSR(I) = U(I)*U(I)+V(I)*V(I>
195 IF ([Link].1.) EPSIL = (TH(!)-(PI*Tn(NT)))/PI
IF ([Link].5.) EPSIL = (TH(1)*TH«2)-TH(NT)-TH(NT-1))/TP-l.
C COMPUTE ARC LENGTH TO FOUKTH ORDER ACCURACY
222
200 SP(1) =- 0.
00 210 I = 2 , N T
DUN = A M A X K . l E - 2 0 . . 5 * 4 B S C T H < I ) - T H ( I - i )))
OX = X X d ) - X X ( I - l )
or = rym-YYd-i)
210 S P ( I ) = SPCI-1)+SQRT<DX*OX+DY*OY)*OUM/SIN(DUM)
ARC = S P ( N T )
SN = 2 t / A R C
SCALE = ,25*ARC
CE = .5*(1.-EPSIL)
00 220 I. =• ItNT
220 S S ( L ) = A C O S U i - S N * S P < L > )
SS(NT) = PI
IF ([Link].O) GO TO 330
CALL SPUIF <NT,SS«THiU,V»Wi3t0.i3iO.)
IF ([Link].5.) 60 TO 232
WRITE <Ni»,<HO> TlTLEtVALtNRN
IF <[Link].N2) WRITE (N2,«0) TITLE. VAi_.NRN
PRINT OUT AIRFOIL OATA
WRITE (Ni»,430)
00 230 L = [Link]
VAL = TH(L)*RAO
SUM =-SN*u«L)/AMAXl(,1£-[Link](SS(L))>
IF (([Link].([Link])) SUM = V(L)*SI6N(SN,FLOAT<L-2))
230 WRITE (Ni».<»eo> XX( L)»YY(L) «SP(L). VALtsUMt V(L) ,W(L)
WRITE (Ntf,i»«»o)
MAKE INITIAL GUESS OF ARC LENGTH AS A FUNCTION OF CIRCLE ANGLE!
232 OX = (XX(NT)-XX(l))/TP
OY = (YY(NT).YY(l))/TP
00 240 1 = ItWC
ANGL = FLOAT(I-1)*PILC
CIRC(I) = ANGL
CX(I) = COS(ANGL)
SX(I) = SIN(ANGL)
YY(I) = 1.
IF <EEtNE.O.) YY(I) = (2.-2.*CX(I))**££.
FAC = SISN(1.+CX(I)«FLOAT(LC-I»
240 SP(I) - ACOS(,5*FAC)
SPJMC) = PI
CIRC(HCJ = TP
IF ([Link].6.) GO TO 2t<t
SCALE = ARC/ARCu(M.1)
00 2H2 L = iiMM
Z(L) = FLOAT(L-1)*OT
CALL SPLIF (MM,Z,A^CL,CO,SI,PHiR,3,0..3iO.)
CALL IiMIPL (N«P»CIf([Link] t Z»ARCL«CO f SI,PHIR)
00 2H5 L = 1,LC
BB(L) = CX(2*L-1)
AA(L) = -SX(2*L-1)
00 AT flOST 100 ITERATIONS TO FIND THE FOURIER COEFFICIENTS
00 320 K = 1,100
CALL lNTPi_(NMPtSPtTT,SS,TH,U,V«H)
00 250 I = It NIP
350 TT(I) = TT(I) + ,5*(CIRC(I)+EPSIL»(CIRC|I»-PD)
ENSURE CLOSURE
223
DUM = o.
SUM = o.
FAC = 0.
00 260 L = ItNMP
OW = OUM - T T ( L )
SUM = S U ! " I - T T ( L ) * C X ( L )
260 FAC = F A C + T T ( U * S X ( U
OU« = OU«/FLOAT(NMP>
DA s 1 . - E P S I L - < O X * S I N ( O U M ) + O Y * C O S ( O U M ) ) / S C A L E - F A C / F L O A T < L C )
08 = < O Y * S l N < D U M ) - O X * C O S ( O U M ) ) / S C A L E - s U M / F L O A T < L C )
00 270 L = [Link]
270 TT<L) = TT<L)+OA*SX(L)-OB*CX<L)
FIND THE CONJUGATE FUNCTION OS
CALL CONj(NMP,TT,DStXX f [Link])
00 290 I = 1.N.1P
SUM = OS(I)
290 OS(I) = YY(I)*EXP(SUN)
OS(MC) = OSU)
CALL S^LIF(MC,CIRC,DS,XX,XC,Z.-3,0.,3.0.)
SCALE = ARC/Z(MC)
EftK = 0.
00 310 I = 1,N.«IP
VAL = ACOS(1.-2.*Z(I)/Z{NO)
,ERR = AMAXUERrt,A8S(SP(I)-VAL)>
310 SP(I) = VAL
IF ([Link].5.) WRITE (N"H<*90) [Link]
IF ([Link]) GO TO 330
320 CONTINUE
WRITE (Nlt,<t50>
330 CALL FOUCF(NMP,[Link],B3,AA)
AA(1) = ARC
AA(2) = l.-EPSIL-(OX*SIN(BB(l)»+OY*CO«5(BB(l)))/SCAL£
B8(2) = (-OX«COS(B^(1))+OY*SIN(BB(1)))/SCALE
IF ([Link].5.) GO TO 3"t2
WRITE (Nlf,i»60) EPSIL, NflP
IF <(FSYM.N£.1.).AND.(FSYI"[Link].3.)) GO TO 341
00 34'4 L = [Link]
3U1* Z(L) = fiLOAT(L-l)*OT
CALL SPLIF(MC«CIRC«SP,U«V«rit3iO.«3«0. )
CALL l N T P L ( M M i Z t D S « C I R C » S P « U i V « W )
CALL SPLIF (NT,SS,9SR,J,VtW,l,0.,l,0.)
CALL IUTPL<MM,OS,[Link].J,v.*o
3i*l IF ([Link].120) GO TO 3<f2
WRITE ^N^t5^»0)
00 3«fO L = [Link]
3<*0 W1ITE (Nl»,490) A A ( L ) « B 3 ( L )
342 CALL HAP
KETJRN
350 IF ([Link].5.) GO TO 355
OxOSl = (XX(2)-XX(1»/SS(2)
OXOS2 = (XX(NT)-XX(NT-1))/(SS(NT)-SS(MT-1))
OYDS1 = (YY(2)-YY(1))/SS(2)
OYOS2 - (YY(NT>-YY(NT-l))/(SS(NT)-SS(MT-l))
355 CALL SPUIF(NT,SS,XX,U,SP,W,l,DXDSlil,oXDS2)
CALL SPLIF(NT,SS,YY,V,TT,OS,1,OYDS1,1.0YDS2)
224
IF [Link].O) 60 TO 397
OC = PI/FLOATINMP)
ERR = SS(NL)
DU.1 = OIS(0.)
FAC = PI/(OIS<PI)-OUM)
00 360 L = 1«MC
360 C I R C < L > = FAC*(DIS(FLOAT(L-l)*OC)-OUMj
CALL lNTpL([Link],SXtSS,XXtUiSP,M)
CALL INTPL<NHPtCIRC,CX,[Link],OS)
SXMC) = XX(NT)
CX(MC) = YY(NT)
SFAC = l./(XX(NT)-XX(ND)
XXNL = XX(NL)
00 370 L = [Link]
CX(L) = SFAOCX(L)
SX(L) = SFAC*(SX(L)-XXNL)
XX(L) = SX(L)
370 YY(L> = CX(L)
WRITE (Ni»,520> IS
IF ([Link]) WRITE (N2,520) IS
IF([Link].O) 60 TO 395
00 IS SMOOTHING ITERATIONS
00 390 K = 1,IS
00 380 L = £,N:1P
XX(L) = S*lOOTH(SX(L-l),SX(L),SX(L*l)«sX(L) )
380 YYCL) = SMOOTH(CX(L-l>,CX(L),CX(U*l)t«!X(L»
00 390 L r 2,NIP
SX(L) = XX(L)
390 CX(L) = YY(L)
395 NT = MC
CALL SPLIF(NT,CIRC,XX,J,SP.W,1»0.,1,0.)
CALL SPl-lF([Link],YYi [Link],DS«l«0.«l,n. )
397 ISY.1 = 0
IF ([Link].5.) GO TO 170
- Sp(l)
= TT(D
U(MT) = SP(NT)
V(NT) = TT(NT)
GO TO 170
1*10 FORMAT (lXl6At<I<t)
HZO FORMAT (5F10.7)
1130 FORMAT (35HOAIRFOIL COORDINATES AND CjRVATURES/lHO,6X.1HX,14X1HY
1 ,9X,10HARC LENGTHi7X3HAN6,6X5HKAPPA«iOXt2HKPillXi3HKPP//)
FORMAT <lHl«<*Xi3HERRil»X,2HOA.l<*X«2HDQ//)
FORMAT (32H FOURI£« SERIES 010 NOT COMtfERGC)
t»60 FORMAT (34HOMAPPINS TO THE INSIDE OF- A CIRCLE//3X11HDZ/DSIGMA =
1 50H «{l/SI6MA**2)*U-SIGrtA)**<l-EPSIL)*<EXP<W<SIGMA))//3X,
2<f2HW{SlBMA> = SUM((A(N)-I*8(N))*SISMA**(N-l>)//3X,7HEPSIL =
3 F5.3«20X«I4,25H POINTS AROUND THE CIRCLE )
H70 FORMAT (1H1)
1*60 FORMAT (Fl2.6,2Fli*.6,Fli*.3,Fl<f.'*,2Em.3)
(»90 FORMAT (3E15.6)
C ****CHANGE (1020) TO (20At) ON IBM 36n****
500 FORMAT U020)
510 FORMAT ( 5X»Fi*. 3, 8X,F5. 3t8X.F1*.3t 10X,F^.3, ItX. 15)
225
SUBROUTINE HAP
SUM UP FOURIER SERIES TO OSTAIN MAPPING FUNCTION
COMPLEX TTtTMP
COnnO.M PHI(162«31),FPU62«31)iA<31),B<31),C<31)«0(31)tE<3l)
1 iRP<31),RPP(31),R(3l),RS<31),RI<31),aA<162),8B<162),Co<162)
2 ,SI(162),PHIR(162).XC<162), YC (162)»Fvi < 162 )t ARCLU62) »DSUM< 162)
3 ,ANGOLD(162),XOLO<162).YOLO<1&2)tARCoLO<162),OELOLO(162)
C0110N /A/ P l , T P . R A O , E 1 « A L P , R N i P C H , X P . T C « C H D . D P H I « C L t R C L f Y R
1 , XA,YA,TE,OT,OR,OELTHiOELR,RAtUCNtOS,MiRA«MEPSILiQCRIT,Cl»C2
2 ,C<ttC5»C6«C7t3£TiBE:TA,FSY,1,[Link]£(lf)
3 ,IK,jKtIZ,ITYPtMOOE,ISnMFC,NCr«NRN,N(; t IOi
t , [Link]
****CrtftNSE TO l.E-6 FOR SINGLE PRECISION ISM 360****
DATA [Link]/-12..10.E-12/
NOTE THAT THE SQUARE OF THE HAPPING PlnOULUS is BEING COMPUTED
MX s fl/2
SET TH£ SINES AND COSINES
com = i.
S I ( l ) = 0.
00 5 L = 1«WX
CO(L+1) = CO(L)*OCN-SI(L)*OSN
COIIH-L) = CO(L+1)
SKL+1) = C O ( L ) * O S N + S I ( L ) * O C N
5 SI(flH-L) = -SKL+1) "
SET MAPPING HOOULUS FOR CUSP AT THE TftIL
00 10 J = liN
FP (1,J) = l.+R(J)*(R(J)-2.)
00 10 L = 1,MX
10 FP(L+1»J) = l.+R(J)*(R(J)-2.*CO(L+D)
IF ([Link].O.) GO TO 30
ADJUST IF THERE IS AN ANGLE AT THE TAIL
00 20 J - liN
FP(1,J) s FP{l,J)**(l.-EPSlL)
00 20 L = 1»MX
20 FP(L+1«J) = FP(L+liJ)**(lt-EPSIL)
NOW COMPUTE CONTRIBUTION FROM FOURIER SERIES
30 00 50 J = l.N
NFCX = l"IINO<NFCtl + INT(POH/ALOG10(R<J)-TOL)))
RJ = 2.*R(J)
K = NFCX
S = AA<K+1)
35 S = R(J)*S+AA(K)
K = K-l
IF (K.6T.1) SO TO 35
FP(1,J) = FP(liJ)*£XP(S*RJ)
00 50 L s [Link]
K = NFCX
226
LX = K*L
LT = M O U ( L X . « )
S = AA(K-fl)*CO(LT*l)
8 = B8<K*1>*SI(LT+1)
<*0 LX = L X - L
LT = . l O O ( L X i M )
S = R< J ) « S + A A < K ) * C O ( L T + 1 )
a = R(J»*a+68(rt)*SI(LTfl)
K = K-l
IF ([Link].l) GO TO 40
OUM = FP(L-H.J)
FPMM-LtJ) = E X P ( R J * ( S - 0 ) ) * D U M
50 FP<L+1»J) = E X P < R J * ( S + 3 n * O U M
00 65 L = 1,M
S = PI-83(D
00 60 K = 1»NFC
LT = 100«L-1)*K»M)
60 S = S +A A ( K + l ) * S I ( t . r + l ) - 8 B < « + l ) * C O ( L J > i )
ANG = F L Q A T ( L - 1 ) * O T
F P ( L i M M ) = 1.
65 FI*I(LJ = S-.5*(ANG+£PSIL*(ANS-PI) J
FW(.HH)- = F«(l)-(l.*[Link])*PI
00 70 J = 1»NN
FP(«M,J) = FP(1,J)
70 FP(H1H + 1.J) = F P ( 2 t J )
COMPUTt: ARC LENGTH AND BODY FROM THE CAPPING BY INTEGRATION
XMIN = 0.
YMIN = 0.
Y H A X = 0.
S = -S3*T<FP(1,1)>
TV.P = Ci»IPLX(S*COS(F«(l) ),S*SIN(FM(1) ) >
00 80 L = [Link]
Q = SO^T(FP(L«l) )
S = S+9
ARCL(L) = S
S = S+Q
TT = C'"l p LX(Q*cOS(F i «l(L) ) » Q * S I N ( F M ( L ) ) )
T«P =
XC(L) -
YC(L) =
IF ( K N - 2 ) 6 0 , 7 0 , 8 0
60 V = ( 6 « * V N - V ) / U
60 TO 90
70 V = \IH
60 TO 90
80 V = (DS*VN+FPP(I))/(1.«-FP<IM
90 B = V
0 = OS
100 OS = S ( J ) - S ( I )
U = FPP<I)-FP<I)*V
FPPP(I) = ( V - U ) / O S
FPP(I) = U
FP(I) = < F U ) - F < I ) ) / O S - O S * ( V + U + U ) / 6 .
V = U
J =I
1 = I-K
IF ([Link].H) 60 TO 100
FPPP(N) = FPPP(N-l) .
FPP(N) = B
FP(N) = OF+0*(FPP(M-l)4-B+B)/6.
IF ( K N . S T . O ) R£TLI«M
IF KM IS h|£6ATIVE COMPUTE. THE JNTE6RAL IN FPPP
FPPP(J) = 0.
V = FPP(J)
105 I = J
J = J»K
OS = S ( J ) - S d )
U = FPP»J>
FPPP(J) = FPPP(n + t 5 * O S * ( F ( I ) * F ( J ) - O S « O S * ( U * V ) / 1 2 . )
V = U
IF ([Link].N) 60 TO 105
RETURN
ENO
SUBROUTINE F O U C F ( N « G « X « A i B )
FOURIER COEFFICIENTS BY FAST FOURIER TRANSFORM
COMPLEX S.ElVi3P,X.6K
OI1ENSION G(1)«X(1)<
DATA PI/3.1«H59265358979/
L - = N/2
V = PI/L
EIV = CI"IPLX(COStV)(SIN(tf))
ENI = l./FLOAT(N)
CALL FFORM(LiG,X»A,B)
GK = 0.
I = 1
00 5 J = 1 < L « 2
X(J) = Cm|pLX(B(I)tA(D)
X(J+1) = X(J)*EIV
1 = Ul
K = L
00 10 J = l.L
230
QP = GK-CONJG<G<J>)
GK = S«+CONJG(G(J)>-QP*X(J)
A(J) = -REAL<GK)*ENI
B(J) = AIMAG(GK)*ENI
GK = 5(K)
10 K = K-l
A<L*1) = .8(1)
B(l) = 0.
8(L+1) = 0.
RETURN
END
SUBROUTINE FFORM(N»FiX,CN»SN)
C FAST FOURIER TRANSFORM
C INPUT ARRAY F UITH REAL ANO IMAGINARY PARTS IN ALTERNATE CELLS
C REPLACED BY ITS FOURIER TRANSFORM
COMPLEX F(l),XU)t*
DIMENSION CN(1),SN(1)
IF ([Link].2) RETURN
NS =• 1
NR = 2
NS = N
11 00 10 K = NR,N
IF <MOO(NQiK).EO.O) GO TO Zi
10 CONTINUE
21 NOsNQ/K
NS = NS*K
NR = <
IQ = 0
ID = 0
00 22 I = liNS
00 24 J = liND
L = MOO(IQ+J,N)
W s F(L)
» s Q
00 26 K =. 2«NR
L S' L*MO
n s MOO(M+IO,N)
26 W = W+F(L)*CMPLX(CM(M+1)
2t X(IO+J) = W
10 = IO+NO
22 IQ = IQ+NQ
N9 = NO
IF ([Link].l) GO TO 61
00 32 K =
32 F(K) = X < K )
61 DO 60 K r NR,N
IF (MOO<NQ«K).EQ.O) GO TO 71
60 CONTINUE
71 NDSNQM
NS * NS*K
231
NR = <
IQ = 0
10 = 0
00 72 I = 1,NS
00 7H J = 1,NO
L = MOO(IQ+J,N)
W - X(U
M = 0
00 76 K = 2«NR
L = L + ND
76 * z W+X(U>*CHPI.X(CN(n+l>tSN(H+lM
7f F(IO+J) = W
10 = IO+MO
72 IQ = I04-MQ
NO = NO
IF ([Link].l) GO TO 11
RETURM
FUNCTION INOEXR(X«ARRAf,N)
DIMENSION ARRAr(l)
S = A8SU-ARRAYIM))
00 10 L = l.M
IF <ABS(X-ARRAr(l».GT.S) SO TO 10
INOEXR = L
S = A3S(X-ftRRAY(L)>
10 [Link]
RETURM
ENO
REAL «IACH,MACHM,NE''«MACHS
OIMENSIOM HP(162),SEPP(162),CPP(162),THETAP(162),DELP(162)
1 «OELX(1),TO(1)
DIMENSION H(1),THETA(1),OELS(1),XX(1).YY(1),HACHN(1)
1 , SEPR ( 1 ) , CPX ( 1 > t OSOT ( 1 ) , S ( 1 ) , MACHS ( 1 ) , ANGNEw ( 1 )
EQUIVALENCE ( MACHN( 1 ) , A( 1 ) ) , ( H( 1 ) , FP ( 1, 6) ) , ( THETA( 1 ) ,FP(
1 t(XX(l)«FP(l,3 )),(YY(1),FP(1,5) ),(OELS(1),FP(1,10) )
2, ( ANGNEW (1 ), FP (l,2t»« (SEPR (D,FP(l,l4)),( CPX (l).PHIR(D)
3 , (S(l) , FP(1,16) ) , (MACHS(l) ,FPU,26) >. (DSOT(l) tFP(ltSO))
"» , (DELX(l),Fp(l,l2)),(TO(l),FPa,20))
232
CP(Q) a C5*«CH/(1.+C2*Q*0)
QSX(8| r <C<f-U.«-a/C5»**U./C7n/C6
MACHO) = SQRT « Q / < C l - C 2 * e n
DATA ISW/0/»COF/0./tXPUT/.5/fXFAC/100 . /
00 10 J = liNN
PHKMH.J) = PHI(1.J)+OPHI
10 PHI(MM+1«J> = PHI(2«J)+OPHI
IF (IS«[Link] CALL SOPLOT(NRN)
COMPUTE AND STORE CP CRITICAL
CPX(i"!imi = CPU.)
JSX SET TO 1 FOR FSYM=1. AND FSYM=3 IF FLOW HAS NOT BEEN COMPUTED
ISX = <NCY+l)*(irYP-3)*ABS<FSY«+10.)+.2
IF ([Link].D SO TO 30
Hif = N5
FSY« = 0.
XC(NM) = 1.
ALP = 0.
XSEP = AHAXKO. iXSEP-1.)
OS = A(«i*I)
00 20 L r 1,MM
XOLO(L) = XC(L)
rOLO(L) = YC(L)
MACHN(L) = MACH(A(D)
20 CPXIL) = CPt«ACHM(l.»
IF ((A8S(YC(MM)-rC(l».L£.l.e-5>«ANO.(IABS(NRN>.ST.999)> 60 TO 50
60 TO 110
30 00 40 L = 2tM
U = (PHI(L+l.l)-PHI(L-l.D)»DELTH-SI(u)
QS = < < J * U ) / F P < L i l l
MACHM(L) = NACH(QS)
HO C P X ( L ) = C P < M A C H N ( L M
MaCHN(.in) = ,5*(NACHN12)+flACHNlM) )
MACHN(l) = HACrlN(M«)
CPX(l) = cP(nACHNd))
CPX(MM) = CPX(l)
3S=aSX(CPX(MM))
IF t F S Y M . E Q . 6 . ) SO TO 60
If ( ( F S Y M . L E . 5 . ) . O K . U T Y P . L E . 8 ) ) 60 To 50
ADVANCE PLOTTER PAPER TO THE NEXT SLAV* PASE
[Link]..5) CALL PLOTtl2.0*FLOAT«I\|T((20.2+XPLT)/12.»,o.«»5)
XPLT = t5
50 CALL GETCP«CDF)
CALL 60PRIN (HP,THETAP,SE^[Link],jfTRANS)
IF ([Link].l) CALL EXIT
ISW = 1
RETURN
60 DO 70 L = [Link]
70 CPP(L) = CPX(L)
lF(([Link]«0).OR,([Link].6.)> 60 TO gO
FIND THE BASE PRESSURE
OELSP = 10.
CPO = CP(|viACHN{IXX-l))
00 80 L = IXXtN
CPN = CP(I"1ACHN(L) )
OELBP s AHIN1([Link],M-CI'0>
233
BO CPO = CPM
SCP = dCP+RBCP*0£LBP
90 ISW = 1
PCH = ftSS(PCH)
IF [Link]) SO TO 110
MODIFY THE HACH OISTRI3UT10N
CPO = CPl,1ACHN(LSEP»
SEPX = XC(LSEP)
SL = <BCP-CPO)/CXC<M«)-SEPX)
00 100 W = LSEPiNM
CPP(L) = CPO+SL*(XC(L)»S£PX)
100 HACHN(L> r HACH<QSX<CPP{L» )
110 KQMIN = 1
KQ*AX = 1
QHIN = MACHN(l)
QMAX = 8MIN
OARC = TP/FLOAT(NPTS-1)
00 115 L = 1,NPTS
115 H(L» = FLOAT (L-1)*OARC
H(MPTS) s TP
00 116 L = 1,M
116 YV(t) = FLOAT «U-1)*OT
yY(KH) = TP
CALL SPLIF (MH,YY,AKCL»OSOT,COtTDt3»0.»3,Ot )
CALL INTPL (NPTS,H,S,YY»AKCL,OSOT«CO,TO)
S(NPTS) = ARCL(HN)
CALL SPLIT (MM,ARCt»[Link],TO,^.0. ,3.0. )
CALL INTPHNPTS.S .IftCHS, ARCL,[Link],CO,TO)
CALL SPLIF <MM,ARCLiXC,DSOr,[Link]«3«0.«3f 0. >
CALL INTPL (MPTS,S«XXiARCL,XC,[Link])
00 120 L =. [Link]
IF ([Link]) KQNAX = L
IF <MACHS<L».[Link]«I:N) KQHIN = L
QMIN = AXUNi(«ACHStL>,3!"IIN)
QMAX = A*)AX1(MACHS(L).S«AX)
SEPR(L) - 0.
H(L) x 0.
OELS(L) = 0.
120 THETA(L) r 0.
' IF (PC*«LT.O.) GO TO l'+0
KQflAX = KaMIN+INOEXK([Link](KQ«IN+l) f MPTS»KQMIN)
IF ([Link]) CALL ABORT
CALL NASHPIC ( K O W f t X i N P T S )
XTRANS = PCH
IF ([Link].O) XTRANS = XX(KQNAX)
K830T = INO£xR(XTRAMStXX«KQl«IIN>
IF ([Link].1) CALL A90RT
CALL MASHMC (KQBOTtl)
THtTA(l)=FAC*THETA(2)+(l.-FAC»*TH£TA(n)
H(l)=FAC*H(2>+(l.-FAC)*H(t)
OELS<1)=H(1)*THETA(1)
COMPUTE THE SKIN FRICTION ORAS
Q = SQ^T(QS)
RT = (C1-C2*QS)/(C1"C2)
234
HBT = <H(NPTS>+1.)*(1.-C2*QS/C1)-1.
HB8 = (H.(i)+i.)*(l.-C2«aS/Cl)-l.
COF = 2.*THETA(NPTS)»Q**(.5*( HBT +5.))*RT**3
COF = COF*2.*THETA(l)*a**«t5*( HBB+5.))*RT**3
IF ([Link].1) GO TO 200
c MAKE DISPLACEMENT IONOTONE INCREASING ON THE UPPER SURFACE
00 170 L = [Link]
IF <D£LS(L+l).[Link]<D) OELS(L+1) = OELS(L)
170 CONTIMUE
C LOWER SURFACE - FIND WHERE DELS START"; DECREASING
C TREAT THE LOWER SURFACE LIKE THE UPPER SURFACE IF [Link].O
XPC = .fad
IF ([Link].O.) XPC = 2.
J = KQBOT
180 J = J-l
IF (DELS(j-i).[Link](J)) GO TO iBb
IF ([Link].2) GO TO 180
GO TO 200
185 IF (XX(J).[Link]) SO TO 190
OELStJ-D =• DELS(J)
GO TO 180
C OISPLACtflENT MUST STAY MONOTONE DECREASING
190 J = J-l
IF (DELS(J-1).[Link](J» OELS(J-l) = OELS(J)
IF ([Link].2) GO TO 190
C SMOOTH OELS IS TINES
200 IF ([Link].O) GO TO 220
00 210 I = liIS
OLO = OELS(l)
00 210 L = [Link]
NEW = OE-LS(L-l)
OELS(L-l) = .25*(OLO*NErt*NE«+o£LS(L))
210 OLO = NtM
220 XPLT = XPLT+.5
FAC=(S{NPTS-l)-S{NPTS))/(S(NPTS-l)-S(MPTS-2»
DELS(NPTS)=FAC*DELS(NPTS-2)<-(l.-FAC)«nELS(NPTS-l)
IF ([Link].D SO TO 260
YFAC = 10./S(NPTS)
OH = (H(KQMAX + 1)-H(KQBOT-H)/ FLOAT(2*KQMAX-KQMIN)
FAC = ARCOLO(NT)/S(NPTS)
[Link].1.2) CALL SYMBOL).33,8.7t,.14,55HOISPLACEMENT THICKNESS
1 AT EACH BOUNDARY LAYER 11£RATION,270.,55)
CALL PLOT (XPLT4-XFftC*0£LS(l) ,10.5,3)
00 230 L s [Link]
CALL ? L O T ( X P L T + X F A C * O E L S ( L ) , 1 0 . 5 - Y F A C « S < L ) , 2 )
IF ( ( L . G E . K O B O T ) . AMD. ( L . L E . K Q M A X ) ) H( L ) = HCL-D+DH
230 Y Y ( L ) = S ( L ) * F A C
YY(NPTS) = ARCOLO(.MT)
C OELX WILL 8£ ROUN04RY LAYER DISPLACEMENT AT NT POINTS
CALL SP(-IP(NPTS,YY,OELS,OSOT,CO,T0,3,n.»3,0.)
CALL INTPL«NT,ARCOLO,OELX,YY,DELS,DSDT,CO,TOj
C THE FOLLOWING ARE 3EINS COMPUTED FOR FUTURE PRINT OUT
CALL SPLIF(NPTS,S,L>ELS,OSUT,CO»TQ,3,0.,3,0.)
CALL !NTPL(Wl«l,ARCL,DELP,S,OELS,OSOT,Co»TD)
CALL SPLIF (NPTS,S,H,DSDT,CO,TQ,3.0.,^,0.)
235
CALL INTPL(NHtARCL,NP,S.H«OSDT»CO»TO)
CALL SPLIF(NPTS,StTHETA,OSOT,CO,TO,3,n.«3.0.)
CALL INTPL <MM,ARCL,THETAP,S t THETA,DSr>T,CO,TD)
CALL SPLJF(NPTS,S, S£P*,OSOT,CO,TO,3, n. ,3,0. )
CALL lNTPL<MM,[Link]«S«SEPR»DSOT,CotTD)
GET THE SLOPES FOR THE OUTER AIRFOIL flT CORRESPONDIMG POINTS
00 2«fO L = [Link]
OOEL = ROEL*(DELP(L)-OSU«(L))
OELP(L> s OOEL
OSJH(L) * DSUNID+ODEL
S(L) = FflC*ARCL(L)
SCIfl) = ARCOLO(NT)
CALL SPLJF(MM,S,Ff1,DSOT,CO,TDt3,0..3,0.)
CALL INTPL(NTtARCOLD,AJ4SNE«,S«FM,[Link],TO)
OELHAX = 0.
00 250 L = [Link]
OOEL - OELX(L)-OELOLO(L)
OELilAX = A M A X K O E L ^ A X t A B S C O O E L n
OY = OELOLD(L)fRDEL*ODEL
AMG = . S t l A N G O L O l D + A N S N E W I L ) )
XX(L)=XOLO(L)
250 OELOLO(L) = OY
ISS = IS
IS - -1
IF (ITrP.EQt99) CALL SOPRIN (HP, [Link])
CALL AIRFOL
IS = ISS
Fsrn = 7.
RETURN
260 00 270 L = ItM*
ARCOLO(L) = ARCL(L)
CPP«L) = CPX(L)
270 ANGOLO(L) = F«(L)
CALL SPLlF(NPTSiS«OELS,OSDT,CO.T0.3.0. .3,0. )
CALL INTPL(MM,ARCLtDSU1,StOELStOSDT,CO«TD)
CALL: SPLIF(NPTStS«SEPRtOSOTtCO«TO«3«0.«3iO.)
CALL I^TPL (MMfARCL»SEPP«SiSEPR,OSOT,(iO,TD)
CALL SPLIF (NPTS,S,TH£TA,[Link],TO,3.0.,3.0. )
CALL INTPL (MWtARCLtTriETAPiS,THETA,OSoT,COtTD)
CALL SOPRIM ([Link]^[Link])
NT = »1i"l
CALL GETCP<COF>
IF ([Link].-l > CALL PLOT (O.,0.,999)
CALL EXIT
ENO
SUBROUTINE GoPRIN<H,TH£TA«SEP»CPP,OEL,XTR)
REAL IftCHN
COMHON PHI(162,31),FP(162,il),A(31),B(31),C(31),D(31),E(31)
1 ,RP<31>,RPP<31),R<31),RS<31),KI<31),AA<162),8B(162),CO<162)
2 ,SI(l62),PHIR(162),XC(162),YC(162),F*l(162),ARCL(162),DSUH(162)
236
3 ,ANGOL0(l62),XOLO<162},YOLO(l62),ARCoLO(162),D£LOL0<162)
COMWOM /A/ PI»TPtRAD,E1tALP«[Link],DPHI,CL«RCL«YR
1 ,XA«YA«TE«OT«DR,OELTH,0£LRtRAiQCN,OSN»RA<HEPSlLiOCRIT,CltC2
2 , C < l , C 5 t C 6 » C 7 t 3 E T . S E T A , [Link] Hit NNiNSP
3 , IK,jK»[Link],MOO£«[Link],N S .IOII'1,N2tN3tN<».NTiIXX
IPEN = 2
IF «NOO<|_+3t55).EO.O) WRITE (N4,330> TON
IF (XOLU(l).[Link]) 60 TO 90
TRANS = 1H
IF (MACHN(U).[Link]) TRANS = lOHSTAGuATION
IF «XOUOU+1>.[Link]).OR.<XOLO<U-l).[Link])) GO TO 65
IF ((XOLOIL+2).ST.XTK1.0R.(XOLU(t-a).[Link]))TRANS= 10HTRANSITION
85 WRITE <N<t,340) XOUO(L)»YOUD(L)[Link],FPPP(L),CPPCU,TRANS
SO TO 100
90 WRITE (Nlf.350) XOLO(U)«YOUOtL)«YS»FPP<UtFPPP(C)«CPPtU)tTHETAjI)
1 .SEP(L)
100 CONTINUE
IF ([Link].O.) NRN = -IABS(NRN)
XP = -A8S(XP)
RETURN
120 WRITE (NH.310)
WRITE (Ni»,300) IOFF
I = 1
YSEP = ABS(XSEP)
IF ([Link].O.) YSEP = 2.
00 150 L = 1,M«
IF (noo<LtSS).eo.o> WRITE <Nt«3oo> IOM
IF (XC(L).[Link]) 60 TO 130
TRANS = 1H
IF (HACHN(L).[Link]) TRANS = 10HSTA&MATION
IF «XCIL+1>. [Link]). OR. (XC(U-U. ST. XrR» SO TO 125
1 = -1
YSEP = ABS(XSEP)
IF ((XC(L+2).[Link]).OR.<XC(L-2).ST,XTR» TRANS = 10HTRANSITION
185 WRITE <Ni»,290) U » X C < U ) « YC(U) tFpP(U) «FPPP(U) iflACHNd.).
1 CP(L),CPP(L),Z,[Link].U
GO TO 15Q
130 BL(1) = 1H
8L(2) = 1H
BL(3) = 1H
BL(<*)a 1H
IF ([Link]) BL(1) = 2 H L S
IF(<SEP<[Link]).AND.<SEP<L+I).[Link])) BL(2)= 2HCS
IF ([Link]) 3U(3)= 2HLM
IF((XC(L).[Link]).AND.(XC(L+I).[Link])) 3L(tl = 2HLP
WRITE (Nit,2eo> BL .UtXC(U).YC<t).FpP(L).FPPP(U)trtACHN(U)t
1 CP(L)fCPP(L),THETA(L),OSUH(L)tSEP(L).rt(L)iO£L(L).L
ISO CONTINUE
GO TO 10
260 FORMAT(Il<H2F9.5 f 2F8.2,2F9.<H
280 FORMAT ( 3 X , 4A2,I5.2F9.5.F9.2«F8.2,F8.<f ,2F9.«f ,F9.5,F9.5.F9.S,F7.2,
1E9.2.I3)
290 FORMAT (Il6t2F9.StF9.2«F8.2tF8t4i2F9.K«2F9.S,8XtA10*7X,IS)
300 FORnAT<Il,l<fXlHL5X2HXS«7Xi2HYSt7Xt3HAMG«i»Xt5HKAPPA|i»XtifHnACH6X2HCP
1 ,6X3HCPl-.tX5HTHETA,5X!tHOEl.S<6X3HSCP,6XlHH,6X2HOO,6XlHl./)
310 FORnAT(lHl«l5X«ifOHLOWER SURFACE TAlt TO UPPER SURFACE TAIL )
330 FORMATdHl/ 17X26HLISTIN6 OF COORDINATES FOR.2X.4A4)
330 FORKATdl /11X1HX, 3X« HY.7X,2HTS,6X, 3HANG. tXt5HKAPPA,6X,2HCP. 5X,
1 5HTHETA,5X,3HSEP/)
FORflAT (F14.5,2F9.5tF8.2tF8.2,F9."*,i*X.A10)
238
2 tC't.C5«CSiC
3 .IK. jK,I2t lTYP..lOOE,ISiNFCi NCT f I
^ , MPTS.H..[Link].1'*
DIMENSION MACHS(l)«H(l),TH£TA(X),SEPR(l),S(l)iOELS(l)tXX(l)
EQUIVAUENCE(|V|ACHS(l)«FP(li26»MH(l),FP(l«&»t(THETA(l)iFp(l(e))
EQUIVALENCE ( SEPR < 1 ) ,F?< 1 • It) ) • < DELSI i > ,FP( 1 , 10» . < S( 1) ,FP
EOUIVALEMCE (XX< 1 ) ,FP< 1,3) )
REAL y)H,HHSQ,NJ,MACHS
OAT A T R t R T H O i T E i , T E 2 < S E P M A X « P I l < 1 I N , P I M A X /.3"*24»320. , 5.E-3,5.E«5,
GAM1 = .5/C2
CSIINF = C<*
INC s ISISN(1,K2-K1)
YSEP = ABS(XSEP)
IF (([Link].O.).AND.([Link].O» YSEP =1.
SEPHAX = SEPH
GE = 6.5
L = Kl
OS = A3S(S(L)-S(L-INO)
10 LP = L*INC
,MH s .S*(H1ACHS(L)*.%IACHS(LP»
iiHsa = HH*MH
CSIH = l.+C2*MHSQ
OSOLO =• OS
OS = ABS(S(LP)-S(D)
OQOS = («ACHS(UP)-1ACHS(L»)/(OS*MH*CSTH)
T = CSIINF/CSIH
RHOH = T**SAM1
NU = T*U,+TR)/(RHOH*(T+TR) J
RTH = RN*MH/(£**NU)
IF ([Link]) SO TO 30
THETAH= RTHO/RTH
THT = THETAH
30 FC = 1.0+.066*lHSa-.008*HH*MHSQ
FR = l.
239
BLOCK OATA
COilMON P H I ( 1 6 2 « 3 1 ) « F P ( 1 6 2 « 3 1 ) » A ( 3 1 ) , B f 3 1 ) , C ( 3 1 ) 1 0 ( 3 1 1 « E { 3 1 )
1 ,«P(31),RPP(31),R(31),RS(3l),Hl(31),aA(162),BB(162),CO(162)
2 ,31(162),PHIR(162),XC(162).YC(162)tFv) f162)tARCL(162).DSUM(162)
i ,AMGOLO(162),XOLO(162),rOLO(162),ARClLO(162),UELOLO{162)
COVI^OM /A/ [Link].E«liALP,RNtPCH,[Link],OPHItt(..RCL»YR
1 ,XA,YA»TE«OT,OR,OELTHiDELRtRAtQCN»OS(g'[Link] l Cl.C2
2 tC4tC3tC6tC7«3eT«aEI'AtFsr.*1 t XS£PtS£[Link]( l »-)«n«NtnfltNN«NSP-
3 ,IK,J<iIZtITYP,!«!OOe,[Link],NCY,NRNiNG,IDIH,N2tN3,N<»,NT,IXX
«t t NPTS,[Link]*l<*
****IQI« MUST 3E SET TO THE FIKST DIMENSION OF PHI****
DATA PI/3.14153265358979/ , EH/.75/ , ALP/0./ t CL/lOO./ .
1 pCH/,07/ t FSYfl/1.0/ t RCL/1.0/ , BETA/0.O/i RN/20.E6/ ,
240
MX =: NX +1
MY = NY +2
MZ == NZ +1
MIT = FIT(NM)
KIT = MIT + 2 -MB
IF ([Link].3I KIT =• 10
JIT = 0
COV = COVO(Nfl)
STRIP = STRIPO(NM)
BETA = BETAO<NM>
WRITE UWRIT,112)
112 FORMATU9HOCHORDWISE CELL DISTKIBUTION IN SQUARE ROOT PLANE,
1 46H AND MAPPED SURFACE COORDINATES AT CENTER LINE/
2 15HO X ,15H SURFACE HEIGHT)
LZ = NZ/2 *1
DO lllf 1=2,NX
11* WRITE (IWRIT.610) AO (I) ,SO (I,L* J
WRITE (IWRIT.116)
116 FORMATdSHO TE LOCATION , 15H POWER LAW )
WRITE (IWRIT.610) [Link]
WRITE (1WRIT.600)
WRITE UWRIT.118)
118 FORMAT(46HONORMAL CELL DISTRIBUTION IM SQUARE ROOT PLANE/
1 15HO Y )
KY =• NY +1
00 120 J=2«KY
120 WRITE UWRIT,610) BO(J)
WRITE (IWRIT«122)
122 FORMATdSHO SCALE FACTOR, 15H POWER LAW )
WRITE UWRIT,610) SY,AY
WRITE (IWRIT.600)
WRITE (IWRIT.12<t)
12f FORMAT(27HOSPANWIS£ CELL; DISTRIBUTION/
1 15HO Z )
00 126 K=2,NZ
126 WRITE (IWRIT.610) CO(K)
WRITE (IWRIT,128)
126 FORMATdSHO TIP LOCATION, 15H POWER LAW )
WRITE (IWRIT.610) ZMAX,AZ
WRITE (IWRIT.600)
WRITE UWRIT,132)
132 FORMAT(19HOITERATIV£ SOLUTION/
1 <t3HOSTRlP WIDTH FOK HORIZONTAL LINE RELAXATION)
WRITE (IWRIT,6iO) STRIP
WRITE (IWRlT,13if)
13<f FORMATdSHO NX «15H MY ,15H NZ )
WRITE (IWRIT,6<fO) NX,NY,NZ
CALL SCCOND(T)
WRITE UWRIT,660) T
WRITE UWRIT,136)
136 FORMATdSHO MACH NO t!5H YftW ,15H ANS OF ATTACK)
WRITE UWRIT,610) FMACHtYAtAL
WRITE UWRIT, 138)
136 FORMATUOHOITERATION,15H COKRECTIOM ,4H I ,tH J ,tH K ,
1 13H RESIDUAL <4H I <4H J t<fH K ,
245
2 10H [Link] REL FCT 1,10H REL FCT 2.10H REL FCT S,
3 10H SETA ,10H SONIC PTS)
NIT = NIT +1
JI1 f JIT +1
PI = PlO(NM)
P2 = paO(NM)
P3 = P30<NM)
IF' ([Link]) PI = 1.
IF ([Link].10) P3 = 1.
c UPDATE POTENTIAL BY RELAXATION
C EACH ITERATION IS ONE STEP IN ARTIFICIAL TIME
C EQUIVALENT TIME DEPENDENT EQUATION IS
c (it -ii**a>*GSS +GMM +GNN +TER,«IS IN [Link] ANO GT
CALL MIXFLO
C 10 = 0 INDICATES DISC FAILURE,RETURN TO PREVIOUS ITERATION
IF ([Link].O) GO TO 151
JO = 0
REMIND Nl
REWIND N2
C UPOATEO VALUES ARE STORED IN DISC FILES 1,2,3 IN ROTATION
C SET FILE NUMBERS FOR NEXT ITERATION
N = Nl
Nl = N2
N2 = N3
N5 s N
C WRITE NUMBER OF ITERATIONS NIT,
C LARGEST CORRECTION DG AND ITS LOCATION IG,JG,KG,
c LARGEST RESIDUAL FR ANO ITS LOCATION IR,JR,KR,
C CIRCULATION EO,RELAXATION PARAMETERS P1,P2,P3 ANO BETA,
C ANO NUMBER OF SUPERSONIC POINTS MS
WRITE (IWRIT,650) [Link],IG,JG»KG,FR,IR,JR,KR,EO(LZ),
1 P1,P2,P3,BETA,NS
C EVERY KIT CYCLES SAVE CURRENT VALUES ON TAPE i*
C TO ALLOW RESTART IN CASE OF MACHINE FAILURE
IF ([Link]) GO TO 251
IF ([Link](OG).[Link]»[Link](OG).LT.10.) GO TO 11*1
C STOP ON ITERATION COUNT OR IF tRROR HEtTS TOLERANCE
c OR IF ITERATIONS DIVERGE.
GO TO 1«>1
c Jo = i INDICATES SUCCESSIVE DISC FAILURES,GIVE UP
151 IF ([Link].l) GO TO 1
REWIND Nl
REWIND "2
JO = 1
C RESET FILE NUMBERS FOR PREVIOUS ITERATION
N = N3
N3 = N2
N2 = Nl
Nl = N
GO TO 1"*1
C GENERATE AND WRITE AERODYNAMIC PARAMETERS FOR EACH SPAN STATION
C READ FROM THE OISC ANO PROCESS SLICES OF THE G ARRAY *=OR FIXED Z,
C REPRESENTING VALUES OF POTENTIAL ON X.Y PLANES
C CONTAINING SUCCESSIVE WING SECTIONS
161 LX = NX/2 *1
246
CALL SECOND(T)
WRITE (IWRIT.660) T
WRITE (IWRIT.600)
C READ FIRST THREE SLICES OF POTENTIAL ftRRAY FROM DISC- FILE
00 162 L = l»3
BUFFER IN CNltl) (6(1,1,L)»G(MX,MY,L))
c RETURN TO PREVIOUS ITERATION IN EVENT OF DISC FAILURE
IF' (UNIT(Nl).GT.O.) 60 TO 151
162 CONTINUE
K a: a
C INCREMENT Z
171 K = K +1
IF ([Link]) GO TO 191
c SHIFT SLICES OF POTENTIAL ARRAY
00 172 J=1,MY
00 172 IsltMX
G(I,J,1> = 6(ItJ«2)
172 6(1,J,2) - G(I«J,3>
C READ SLICE OF POTENTIAL ARRAY FROM OIsC FILE
BUFFER IN (Nl.l) (611,1,3)«G(MX,NY»3)|
c RETURM TO PREVIOUS ITERATION IN EVENT OF DISC FAILURE
IF (UNIT(Nl).GT.O.) GO: TO 151
IF ([Link].KTE2) GO TO 171
Z = SCALZ*CO(K)
11 = ITEKK)
12 = ITE2(K)
C CALCULATE SURFACE SPEED SV,MACH NUMBER SM,PRESSURE COEFFICIENT CP
C ANO COORDINATES X»Y OF WING SECTION
CALL VELO tK,2,SV,SM,CP»X,Y)
CHORD(K) = X(Il) -X(LX)
C CALCULATE SECTION LIFT,DRAG ANU MOMENT COEFFICIENTS
CALL FOKCF (11,12tX,Y,CP,[Link](R),n..SCL(K)tSCO«K),SCM(K))
IF ([Link]) <»0 TO Ig5
WRITE (IWRIT.600)
WRITE (IWRIT.182)
182 FORMAT(2i»HOSECTION CHARACTARISTICS/
1 15HO MACH NO ,15H YA« «15H ANG OF ATTACK)
WRITE (IWRIT.610) FMACH,YA,AL
WRITE ([Link])
IS* FORMAT(15HO SPAN STATION,15H Cl- «15H CO ,
1 15H CM )
185 WRITE (IWRIT.610) Z,SCL(K),SCO(K),SCM(K)
C IF KPLOT = 0 LIST ANO [Link] CP
IF ([Link].O) CALL CPLOT (II,I2,FMAcH,X,Y,CP)
IF ([Link].2) GO TO 171
C IF KPLOT s 2 GENERATE CALCOMP PLOT OF SECTION CP
CALL GRAPH (IPLOT,II,I2»X,Y,CP,TITLE,[Link],
1 Z,SCL«)«SCO(K),CHOROO)
IPLOt = 0
GO TO 171
c CALCULATE TOTAL LIFT .DRAG AND MOMENT COEFFICIENTS
191 CALL TOTFOR(KTE1,KTE2,CHORO,SCL,SCO,SCM,CO,SCALZ«.23,
1 CL,CDl»CMP,CflR,CMY)
C01 s CYAMtCOl
CO s COO +C01
247
VL01 = 0.
IF (ABS(C01).GT.l.£-6) VL01 = CL/C01
VLO = Oi
IF <ABS«CO)tGT.l.E-6) tfLO = CL/CD
WRITE (IWRIT.600)
REWIND Nl
CALL CHARTY
WRITE <IWRIT,600)
WRITE (IWRIT.192)
192 FORMAT121HOWING CHARACTARISTICS/
1 ' 15HO MACH NO ,15H YAW ,15H ANG OF ATTACK)
WRITE UWRIT.610) F M A C H , Y A * A L
WRITE <IWRIT,19<H
19t FORMATU5HO CL ,15H CO FORM ,15H CO FRICTION ,
1 15H CO ,15H L/D FORM ,1SH t/0 )
WRITE (IWRIT.610) C L » C 0 1 » C O O » C O » V L 0 1 , V L D
WRITE <IWRIT,196)
19& FORNATUSHO CM PITCH i!5H CM ROLL »15H CM YAW )
WRITE (IWRIT.610) CMP»CMRtCMY
REWIND Nl
IF ( K P L O T . L T . l ) GO TO 201
IF KPLOT GT 0 GENERATE THREE DIMENSIONAL CALCOMP PLOT
CALL THREEO<IPLOT*SV,S«*CP.X,YtTITL£,[Link]»CO,CHOROO>
IPLOT = 0
10 = 0 INDICATES DISC FAILURE,KETUKN TO PREVIOUS ITERATION
IF (IOtEQ.0) GO TO 151
STOP ON OPERATOR COMMAND
201 IF ([Link].l) GO TO 501
IF (FHALF(NM).EQ.O.) SO TO 1
REFINE GRID IF FHALF ME 0
NX = NX +NX
NY = NY +NY
NZ = NZ +NZ
RECALCULATE MESH LOCATIONS ON KEFINEO GRID
CALL COORO (NX,NY«NZ,xrEOtZTIP«XMAX ( ZqAX,SY<SCALtSCALZ,
1 AX t [Link],AO,Al,A2,A3,BO f Bi,B2,B3,CO,Cl,C2,C3)
INTERPOLATE UNWRAPPED SURFACE ON REFIiyED GRID
CALL SURF (NOfNE,MC,NX«NZtKSYfl«NP,KTElfKT£2,ITCl«IT£2,XV,
1 YAW,SCAL«SCALZiZS«XS,YS,SuOPT,TRAIL«
2 SO,ZOiAO,[Link]«Ul,[Link]»X,Y»lND)
IND = 0 INDICATES SPLIME FAILUKE DUE TO BAD DATA,SIl/E UP
IF ([Link].O) GO TO 291
INTERPOLATE POTENTIAL ON REFINtD GRID
CALL REFIN
10 = 0 INDICATES DISC FAILURE,KETURN TO PREVIOUS GRID
IF (10.£0.0) GO TO 221
REWIND Nl
REWIND N2
NSMOO = -FHALF(NM)
IF FHALF LT 0 SMOOTH INTERPOLA'ED POTENTIAL ABS(FHALF) TIMES
IF ([Link].l) GO TO 211
00 202 N=l,NSMOO
CALL SMOO
10 = 0 INDICATES DISC FAILURE,RETURN TO PREVIOUS SKID
IF ([Link].O) GO TO 221
248
REWIND Nl
202 REWIND N2
C RESET FjLE NUMBERS
211 N = Nl
Nl = N2
N2 = N3
N3 = N
C INCREMENT NUMBER OF MESH
NM = NM +1
NIT = 0
SO TO 111
C RESTORE PREVIOUS GRID
221 .MX = NX/2
NY = NY/2
NZ - NZ/2
C RECALCULATE MESH LOCATIONS UN PREVIOUS GRID
CALL COORD (NX,NY,NZ,XT£0,ZTIP,XNAX,Z.v|AX,SY,SCAL,SCALZ,
1 AX,[Link],[Link].A2«A3,BO,Bi,B2.B3,CO»Cl,C2,C3)
C INTERPOLATE UNWRAPPED SURFACE ON PREVIOUS GRID
CALL SURF (NO,N£,.MC,[Link]»[Link],KTE2«IT£l,ITE2,IV,
1 YAW,SCftL,SCALZ,ZS,XS,YS,SLOPT,TRAIL,
2 SO,ZO» AO,CO,[Link], 1)1,02,03»X,Y»IND)
c INO = o INDICATES SPLINE FAILUKE DUE jo BAD DATA,GIVE UP
IF ([Link].O) GO TO 291
SO TO 151
c WRITE THREE COPIES OF; INFORMATION NEEDED TO RESTART ON TAPE
251 Kl a. KTEl -1
K2 = KTE2 +ITE2(KT£2) -NX/2
00 252 M=l,3
WRITE (4) NX,NY,NZ,N«,KliK2,NIT
00 262 K=liMZ
BUFFER IN <NI,D <G<I,I,I>,G<MX,MY,IM
c RETURN ro PREVIOUS ITERATION IN EVENT OF DISC FAILURE
IF (UNrT(Nl).GT.O.) GO TO 281
262 WRITE (if) «G<I,J,1),I=1,MX)«J=1,MY)
REWIND Nl
WRITE (t) (EO(K).K=KltK2)
ENDFILEi if
252 CONTINUE
REWIND t
c ALLOW OPERATOR TO STOP CALCULATION
CALL SSWTCH([Link])
IF ([Link].l) GO TO 161
JIT a 0
IF ([Link](OG).[Link](OG).LT.10.) GO TO
SO TO 161
281 REWIND H
SO TO 151
291 WRITE (IWRIT,600)
WRITE (IWRIT,292)
292 FORMAT(2<fHOBAD DATA«SPLlNE FAILURE')
SO TO 1
C TERMINATE CALCOMP FILE
301 IF ([Link].O) CALL PLOT(0.,0,i999)
STOP
249
500 FORMAT(IX)
510 FORHAT<8E10.7>
530 FORMAT<20At)
600 FORMAT(lHl)
610 FORMATtF12.H,7F15.'H
620 FORMAT(8El5.5>
630 FORMAT(1HO»20A<»)
640 FORMATt16,7115)
650 [Link])
660 FORMATU5HOCOPIPUTING TIM£,F10.3»10H SECONDS)
END
NU = FNU
NL = FNL
N = NU +NL -1
C FSYH = 1 INDICATES SYMflETRIC PKOFILE
C FOR WHICH ONLY THE UPPEK SUKFACE IS R^AD
C NU AND NL ARE NUMBERS OF UPPER AND LOuER SURFACE POINTS
REAO (IREAO,500)
C REAO TRAILING EDGE INCLUDED ANGLE AND SLOPE.
C AND COORDINATES OF SINSULAR POINT
REAO (IREA0.510) [Link].XSIN6,YSING
REAO (IREA0.500)
c REAO UPPER SURFACE COORDINATES
00 12 I=NL»N
12 REAO (IREAD.510) XP(I).YP(I)
L = NL +1
IF ([Link].O.) GO TO 15
REAO (IREAD.500)
C REAO LOWER SURFACE COORDINATES
00 14 1=1,NL
REAO (IREAD.510) y/[Link]
J = L -I
XP(J) - VAb
!*• YP(J) = OUfl
GO TO 21
15 J = L
DO 16 I=NL»N
J = J -1
XP(J) = XP(I)
16 YP(J) = -YP(I)
21 WRITE (IWRIT.600)
WRITE (IWRIT.22) ZS(K)
22 FORMATU6HOPROFILE AT I = .F10.5/
1 15HO TE ANGLE ,15H TE SLOPE ,15H X SING
2 15H Y SING )
WRITE (IWRIT.610) [Link],YSING
WRITE UWRIT.2IM
1
2 * FORMATdSHO X ,15H Y
00 26 1=1,N
26 WRITE (IWRIT.610) XP(I)«YP(I)
C SCALE AND ROTATE PROFILE
31 SCALE = CHORO/(XP(l) -XP(NL))
XX = XP(NL) *(XSINS -XP(NL))*THICK
YY = YP(NL) *(YSING -TP(NL))*THICK
CA = COS(ALPHA)
SA = SIN(ALPHA)
00 32 1=1.N
XS(I.K) = SCALE*<(XP(I) -XX)*CA +THICK*<YP(I) -YT)*SA)
32 Y S < I , K > = SCAUE*(THICK*(YP(I) -YY)*CA -<XP(I) -XX)*SA)
SLOPT(K) = THICK*SLT -TAN(ALPHA)
TRAIL(K) = THICK*TRL/RAD
NP(K) = N
NP is NUMBER OF POINTS DEFINING PROFILE
CHORDO = [Link])
XTEO = AMAXl([Link](liK)l
CHORDO AND XTEO ARE MAXIMUM CHURD AND KEARMOST EXTENT OF WING
251
IF ([Link].O..[Link].O.) ISYM = 0
ISYM s 1 INDICATES SYMMETRIC WING
WRITE (IWRIT,52) ZS(K)
52 FORMAT(27HOSECTION DEFINITION AT Z a .F10.5/
1 15HO CHORD .ISHTHICKNESs RATIO,15H ALPHA
WRITE (IWRIT,610) CHORD,THICK,AL
K = K +1
IF ([Link]) GO TO 11
20 = .5*(ZS(1> +ZS(NC))
00 62 K=1,NC
62 ZS(K) = ZS(K) -ZO
ZTIP = ZS(NC)
ZTIP IS TIP LOCATION AFTER WIN<i HAS BEEN C£NTEBEO AT Z s. Q.
RETURN
500 FORflAT(lX)
510 FORHATCflElO.T)
600 FORNAT(lHl)
610
END
KY = NY +1
OX = 2. /NX
OY = l./NY
OZ = 2./NZ
SELECT POWER LAWS
AX = .5
AY = .5
AZ = .5
SY = .5
SY SCALES Y SPACING RELATIVE TO X SPACING
SCAL = XT£0/(.50001*XHAX*XMAX)
SCALZ = ZTIPm.000001*ZMAX)
V2 = (DX/DY)**2
Wl = SCAL/SCALZ
W2 = «W1*OX/OZ)**2
GENERATE X MESH
00 12 1=2, NX
252
00 = (I -1)*OX -1.
B = 1.
IF (ASS(OO).GT. X M A X ) GO TO 13
DO = DO
01 = 1.
02 = 0.
GO TO IH
IF ([Link].O.) B = -1.
A = 1. - ( ( 0 0 -B*X«AX>/(1. .XMAX))**2
C = A** X
D = (AX + AX -l.)*(l. -A)
00 = 8*X 1AX +(00 - B * X H A X ) / C
01 = A*C (1. 4-0)
02 = -(A +AX)*(DD -B*XMAx)
*(3. +0)/((1. -X«AX)**2)
AO(I) = 00
Aid) = ,5*D1/OX
A2U) = 01*01
12 A3(I) = ,5*DX*02
GENERATE Y MESH
00 22 J=2,KY
00 (KY -J)*OY
A 1. .00*00
C A**AY
0 (AY + AY -!.)*(!. A)
01 A * C / ( (1. +0)*SY)
B0( J) SY*00/C
Bl(J) ,5*Ol/OY
82(J)
22 B3(J) = -AY*00*OY*<3. +0)*A)
GENERATE Z ilESH
00 32 K=2,NZ
DO (K -1)*OZ -It
B 1.
IF (ASS(OO) >GT< ZHAX) GO TO 33
00 DO
01 1.
02 0.
GO TO i
IF (DO, LT.O.) B = -1.
A = 1. -((00 -B*Zf1AX)/(l. .ZMAX»**2
C = A**AZ
0 = (AZ +AZ -lt)*(lt -A)
DO = B*ZMAX +(00 -B*Z«AX)/C
01 = A*C/(1. +0)
02 = -(AZ 4-AZ)*(DD -B*ZMAX)
*(3, +0)*A*(1. -ZflAX)**2)
CO(K) = 00
CKK) = .5*01*W1/DZ
C2(K) = 01*01*^2
32 C3(K) = ,5*OZ*02
ENO
253
71 00 72 I=I1»I2
72 IV(I.K) = 2
C SEARCH FOR POINTS ON VORTEX SHEET AT I INDICES OFF «IN6 SURFACE
M =11-1
00 7» 1=2,»
c DETERMINE: STREAMWISE PROJECTION ON SINGULAR LINE
ZZ = Z -TYAW*AO(I)*AO(I)
C SET IV TO INDICATE VORTEX POINT
C IF PROJECTION IS 8EYONO PROJECTION OF UPSTREAM TIP
IF ([Link])) IV(I,K) = IVO
7H CONTINUE
« - 12 *1
00 76 I=M,NX
C DETERMINE STREAMWISE PROJECTION ON SINGULAR LINE
ZZ = Z -TYAW*AO(I)*AO(I)
C SET IV TO INDICATE VORTEX POINT
C IF PROJECTION IS BEYOND PROJECTION OF UPSTREAM TIP
IF ([Link](KTEl)) IV(I,K) = IVO
76 CONTINUE
KTE2 = K
GO TO 11
C Z IS BETONO LAST SPAN STATION
C SEARCH FOR POINTS ON VOKTEX SHEET
61 00 82 1=2,NX
C DETERMINE STREAHWISE PROJECTION ON SINGULAR LINE
ZZ = Z -TYAW*AO<I)*AO(I)
C SET IV TO INDICATE VORTEX POINT
C IF PROJECTION IS WITHIN PROJECTION OF DOWNSTREAM TIP
IF ([Link](NC>.[Link](KT£l)) ItMItK) = IVO
82 CONTINUE
60 TO 11
91 N = KTE2
IF [Link].) GO TO 93
C PROJECT DOWNSTREAM SIDE EOSE POINTS. ON SINGULAR LINE
10 = ITEKKTE2) +1
DO 92 I=IOiLX
N = N +1
92 ZO(N) = SCALZ*CO(KTE2) -TYA«*AO(\)*AO(I)
93 I = ITEl(KTEl)
ZO(KTEl-l) = SCALZ*CO(KTE1-1) - T Y A w » A O ( I ) * A O ( I )
ZOIN+D = SCALZ*CO(KTE2+1)
C LOCATE POINTS JUST BEYOND EDGE OF WINg OR VORTEX SHEET
DO 102 K=2iNZ
00 102 I=2iNX
IF (IV(I,K).GT.O) 60 TO 102
IF (IV(I+1<K+1),GT.O.O«.IV(I-1»K4.1).GT.O) IV(I»K) = IVl
IF (IV(I*1.K-1),[Link]<.IV(I-l,K-l).GT.O) IV(IiK) = IVl
102 CONTINUE
RETURN
END
257
SUBROUTINE ESTIM
C GENERATES INITIAL ESTIMATE OF POTENTIAL
C SUCCESSIVE SLICES OF THE' S ARRAY , REPRESENTING VALUES OF POTENTIAL
C ON X-Y PLANES AT SUCCESSIVE VALUES OF Z,ARE GENERATED
C ANO STORED ON TWO OISC FILES TO PROVIoE BACK UP
C IN EVENT OF SUBSEQUENT OISC FAILURE
COMMON G(193,26«<t) «SEPl,
1 A0(193) ,SEP2,A1<193) ,SEP3. A2( 193) ,SEP4» A3t 193) ,SEP5,
2 80<2&),SEP&»Bl<26)»SEp7,Ba<26),SEpe,B3(26>,SEP9»
3 CO(33),SEP10,C1(33),SEP11,C2(33),SEP12,C3<33),SEP13,
f SOU93,33),SEPl'*,EO<129),sEP15,ZO<129),S£Pl&,
3 IV< 193,33), SEP17, 1 TEK33).SEP18,ITE2(33),SEP19,
6 NX,[Link],KTEl,KTE2»KSYN,ScAL,SCALZ«
7 YAW, [Link], ALPHA, CA,[Link],N1«N2»N3, 10
MX = NX +1
KY = NY +1
MY = NY +2
MZ = NZ +1
C SET THE G AKRAY TO ZERO
00 12 1=1,193
00 12 J=l,26
00 12 K=l,i*
12 G(I,J,K) =• 0.
K =1
C SET VALUES OF POTENTIAL AT DUMMY POINTS BEHIND BOUNDARY
21 00 22 1=2, NX
G(I,KY*1,D s: 0.
IF (IV(I,K>.LT.2) GO TO 22
C IV = 2 INDICATES POINT ON WING SURFACE
C SET POTENTIAL BELOW SURFACE TO SATISFY BOUNDARY CONDITION
DSI = SOU + l.K) -SO(I-l.K)
OSK = SO(I.K-fl) -SO(I,K-1)
SX = A1(I)*OSI
SZ = ClfK)*OSK
U = CA*AO(I) 4-SA*SO(I«K)
W = SYAW
FH = AO(I)*AO(I) +SO<I,K)*SO(i,K)
V = Bl(KY)*(l. +SX*SX +FH*Sz*SZ)
1 -(CA*SO(I«K) -SA*AO(I) *U*SX +FH*W*SZ) /V
22 CONTINUE
WRITE SLICE OF POTENTIAL ARRAY ON TWO OISC FILES
BUFFER OUT(N3,1) ( G( 1, 1,1 ) ,G(MX,MY,1 ) )
GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N3).GT.O.) GO TO fl
BUFFER OUT(Nl.l) (G( 1, 1,1) ,G(MX,MY,l) )
GIVE UP IN EVENT OF OISC FAILUKE
IF (UNIT(NI).GT.O.) GO TO 41
INCREMENT Z
K = K +1
IF ([Link]) GO TO 21
SET TRAILING JUMP £0 IN POTENTIAL TO ZERO
Kl = KTE1 -1
K2 = KTE2 tITE2(KT£2) -NX/2
DO 32 K=K1»K2
258
32 EOfK)
SET 10 TO INDICATE SUCCESSFUL COMPLETION
10 = 1
RETURN
SET 10 TO INDICATE DISC' FAILURE
41 10
RETURN
END
SUBROUTINE MIXFLO
c UPDATES POTENTIAL at RELAXATION USING ROTATED DIFFERENCE SCHEME
C EQUIVALENT TIME DEPENDENT EQUATION IS
C (1. .M**2)*GSS +GPIM +GNN +TERMS IM [Link] AND GT
c SUCCESSIVE SLICES OF THE G ARRAY .REPRESENTING VALUES OF POTENTIAL
c ON x.y PLANES AT SUCCESSIVE VALUES: OF [Link] READ
C FROM ONE DISC FILE,UPOftTEU,AND WRITTEN ON A SECOND DISC FILE
C THREE SLICES ARE REQUIRED FOR COMPUTATION
C A FOURTH SLICE IS USED AS A BUFFER FOR DISC OPERATIONS
C INPUT AND OUTPUT BY SUFFER IN AND BUFFER OUT PROCEED IN PARALLEL
C WITH COMPUTATION
C IF THE BUFFER OPERATION IS NOT YET FINISHED,
C THE IF UNIT TEST DOES NOT RETURN CONTROL TO THE CENTRAL PROCESSOR
C UNTIL ITS COMPLETION,PREVENTING PREMATURE PROCESSING
COMMON) G(193.26,<»),SEPl«
1 AO(193)iSEP2,AK193),S£P3.A2(193),SEPt,A3(193),SEP5,
2 B0<26) .SEP6.BK26) .SEP7.B?(26).SEP8,B3<26),SEP9,
3 CO<33),SEP10,CK33),S£P11,C2(33),S£P12,C3(33),SEP13,
H SO<193,33),S£Pm,EUU29),sEPl5,ZoU29),S£Pi6,
5 IV(193.33),SEP17,IT£l(33).SEP18,ITE2(33),SEP19,
6 NX«[Link]«[Link]«ScAL«SCALZ.
7 YAW,[Link],CA,SA,FMACH,Nl<N2.N3,lO
COMWON'FLO/ GK1(193,26),BUF1,GK2U93,?6),BUF2,
1 SXX(19S),BUF3,SXZ(193)«BUF4«SZZ(193).BUF5.
2 SXI193).BUF6«SZ(193).[Link](193), t)UF8,Kl(193),8UF9,
3 C(193).BUF10.D(193)«3UF11.G1(193),BUF12«G2(193).
H STRIP,Pl,[Link],[Link],JR,KR,DG«[Link],[Link]
COMMON/SWP/ G10(26).SP41,G20(2&),SPA2.G30(26),SPA3tGtO(26).
1 IlfI2.K.L«MO«[Link],[Link]«Qlt82*[Link]
LX = NX/2 +1
MX = NX +1
KY = NY +1
MY = NY +2
TYAW = .5*SCAL*SYAW/CYAW
OX = 2./NX
Tl = DX*OX
AAO = l./FMACH**2 +.2
Ql = 2./PI
Q2 = 1./P2
FR = 0.
IR =0
JR =0
KR =0
259
06 = 0.
IS 30
JS s. 0
KG =0
NS = 0
C FR,IR,JR AND KR ARE VALUE AND LOCATION OF LARGEST RESIDUAL
C [Link] AND KG ARC VALUE AND LOCATION) OF LARGEST CORRECTION
C NS IS NUMBER OF' SUPERSONIC POINTS
C START AT THIRD ROW IF FLOW IS SUPERSONIC AT INFINITY,
C REQUIRING CAUCHY DATA
Kl = 2
IF ([Link].l.) Kl = 3
c DEFINE CENTRAL STRIP OF x-r PLANE FOR HORIZONTAL LINE RELAXATION
c EXTENDING FROM i = 11 TO i = 12 WITH UIOTH DEFINED SY STRIP
C STRIP = 0. ELIMINATES THE CENTRAL STRiP
C STRIP = 1. ELIMINATES THE OUTE« STRIPS
F = ABS(.5*STRIP*NX)
L a. F
IF ([Link]/2) L = L -1
11 = LX -L
12 = LX +L
IF ([Link].O) 12 - LX -1
C READ FIRST THREE SLICES OF POTENTIAL ARRAY FROM FIRST DISC FILE
DO 2 L=li3
BUFFER IN (Nl»l) <G(1,1,L)tG<MX,MY»L))
C GIVE UP: IN EVENT OF DISC FAILURE
IF (UNIT(Ml).GT.O.) GO TO lOl
2 CONTINUE
C SAVE OLD VALUES OF POTENTIAL AT UPSTREAM Z STATIONS
C TO GENERATE CORRECT MIXED SPACE-TIME (DERIVATIVES
DO 4 [Link]
DO f 1=1,MX
G(I,J,t) = G U « J t l >
GKKI.J) = 6(1,J,l)
<» GK2(I,J) = 6(1,J.I)
K =2
L = 2
NO = KTE1 -1
IF ([Link]) GO TO 11
C ADVANCE' AN EXTRA SLICE IF THE FLOW IS SUPERSONIC AT INFINITY
C WRITE FIRST SLICE1 OF UPOATEO POTENTIAL ARRAY ON SECOND DISC FILE
BUFFER OUT(N2,1) (6<lil«<M i6(MXtNYi<M )
C GIVE UP IN EVENT OF DISC FAILURE
IF (UNITIM2).GT.O.) GO TO 101
C READ FOURTH SLICE OF POTENTIAL ARRAY FROM FIRST DISC FILE
BUFFER IN (Nltl) (G( 1,1,<t) tG(HX,flYt«M )
GO TO 5l
C WRITE SLICE OF UPOATEO POTENTIAL ARRAY ON SECOND DISC FILE
11 BUFFER OUT(N2«1> <GU,1,t),G<MX,HY»4) )
Z = SCALZ*CO(K)
DO 12 J=1,MY
GIO(J) = G(I2,J,2)
G20(J) = 6(12-1,J,2)
G30(J) = 6(11,J,2)
12 GifO(J) =
260
SUBROUTINE YSWEEP
ROW RELAXATION
COMMON G(193,26.'M.S£P1,
1 AO(193),SEP2,A1<193),S£P3.A2(193).S£P'»,A3(193).SEPS,
2 BO<26>,SEP&»81(26).S£P7.B2(26),SEP8,B3(2&).SEP9.
3 CO{33),SEP10.C1<33).SEP11,C2<33),S£P12,C3(33),SEP13,
H SO(193.33),SEPm,EO<129),sEP15,ZO(129),SEP16,
5 I V( 193,33), SEPl7,IT£l<33).[Link]<33),SEP19,
6 NX,NY,!MZ,KT£l,KTE2»KSYI«l,ScAL,SCALZ.
7 YAW,CYAW,SYAW,ALPHA,CA,[Link],N1«N2.N3« IQ
COMMON/FLO/ GKl(193,26)tBUFl,GK2U93,3&).BUF2,
1 SXX(193),3UF3,SXZU93),BUP<t,SZZ(193).BUF5»
2 SX(193),BUF6,SZ(193),[Link](193),8UF8,RH193),BUF9,
3 C < 193 ), 8UF10, D ( 193), [Link] 193 ),8UF12, 62(193),
4 STRI P, PI, P2.P3, BETA, FR,IR,jR,KR, 06,15, JS, KG, NS
COMMON/SWP/ 610 ( 26 [Link], 620 (2&1.SPA2, 330(26), SPAS, 640 ( 26 ),
1 Il,I2,K,L,NO,LX,MX,Kr,MY,Tl,AAO,Ql,Q2,Z,TYAW
Jl = 2
IF ([Link].l. ) Jl = 5
C(I1-1) = 0.
0(11-1) = 0.
00 12 1 = 11, 12
RO(I) 1.
Rl(I) 1.
Gl(I) G(I,J1-1,L)
12 62(1) G (1.1. L)
J Jl
13 12
31 BC -Tl*81(J)*Cl(K)
00 32
AS -Ti*Ai(D*ai( j>
AC T1*A1(I)*C1(K)
YP SO(I,K) tSO(J)
A 1, -RO(I) +AO(I)*AO(I) +YP*YP
H RO(I)/A
FH RO(I)*A
OGI -GU-l.J.L)
06J = 6<I»J+1»L) -61(1)
D6K -GKl(I.J)
0611 6(1+1,J,L) -6(1,J,L) -G<I.J,L) +6(1-1,J,L)
+A3(I)*OGI
OGJJ 6(I,J+ltL) -6(1,J.L) -G(I.J.L) +6(1.J-l.L)
-B3(J)*06J
OSKK G(I,J,L+1) -6(1,J,L) -6(1.J.L) +6(I.J,L-1)
+C3(K)«D6K
DGIJ 6(1+1.J+l.L)
-6(I+1.J-1.L) +6(1-1»J-1,L)
262
AXZ = e.*s*u«*AC
BXX = (03 -UU)*A2dI
BYY = (F*Q9 -VV)*B2(J)
BZZ = (FH*OQ -W*)*C2(K)
BXY = .(QQ*SXd) +UV)*(AB +ABl
BYZ = -(Q9*HZ *VM)*(8C +BC)
BXZ 3 -UH*(AC +AC)
AQ = AA/BQ
OELTAG = BXX*OGII *BYY*OGJJ +BZZ*OGKK
1 +BXY«OGIJ +BYZ*OGJK 4-BX?*OGIK
OGII = Gd«J.L) -Gd«,J,L) -6(lM,J,L) *GdMM,J,L)
1 *A3(I)*06I
OGJJ a Gd.J.L) -6d,J-l,L.) -G(I«J-1,L) +62(1)
1 -B3( J)*DGJ
OGKK = G(I«JtU -6(I,J,t-l) -G(I.J,L-1) +GK2(ItkJ)
1 +C3(K)«DGK
DGIJ = G(I.J.L) -5(11, J,L)
1 -Gd»J»l.L) +G(IM«J»1.L)
OGIK = G(I.J.L) -G(1,J,L-1)
1 -G(IM,J,L) +G(IH,J,L-1)
OGJK = Gd.J.U) -G(I,J,L-1)
GSS = AXX*DGII *AYY*D6Jd +AZZ*QGKK
1 >AXY*OSIJ +AYZ*OGJK -t-AX7*OGIK
B = .5*(AQ -l.)*(AXX t-AXX *AXY +AXZ)
BP = AQ*BXX -d. -s)*»
BH 3 AQ*BXX -d. +s)*a
3 3 .AQ*(BXX *BXX +Q2*(BYY +BZZH
1 +(AQ -1.)*(2.*<AXX *AYY +AZZ) +AXY +AYZ +AXZ)
R = (AQ -l.)*SSS +Aa*DELTA6 +R
35 IF (ABS(R) •[Link](FR) ) Go TO 37
FR s R
IR 3 I
JR = J
KR = K
37 R 3 R -AYT*(Gld) -Gd,J-l,|))
1 -AZT*(5K1( I. J) -Gd«J.C-i))
8 3 B -AXT -AYT -AZT
BM 3 BM +AXT
B 3 l./(8 -BM*Cd-l)(
C(I) = B*BP
32 0(1) 3 B*(R -BM*Od-l) )
CG 30. ^
I = 13
DO HZ «=•!!iIS
CG 3 od) -cd)*cG
IF (ABS(CG).[Link](OG) ) GO TO <f3
OG = CG
IG = I
JG = J
KG = K
Hi 62 d) = Gl(I)
61 d) 3 Gd.J.L)
GK2(I,J) a GKld.J)
GKld.J) = Gd.J.L)
264
SUSROUTINE XSWEEP
COLUMN RELAXATION
COMMON G(193,26,i*).SEPl,
265
1 AO<193).SEP2,AK193)«SEP3.A2<193),SEPH,A3«193),SEP5»
Z BO<26),S£P6,Bl<26)tS£P7«Ba<26>*SEP8,B3(26>,SEP9,
3 CO(33),SEP10,CK33>,S£P11.C2<33),SEP12,C3(33),S£P13,
H SO(193,33),SEPl«f,EO(129),s£Pl5,ZO<129),SEP16,
3 IV<193,33),SEP17,ITEK33).SEPie,ITE2<33),SEPl9.
6 NX,NY,NZ,KT£l,KTE2»KSYM,ScAL,SCAI.Z*
7 YAW,CTAW,STAH,AUPHA,CA,[Link],Nl»N2,r43,10
COMMON/FLO/ GK1(193,26) .BuFl,GK2d93»?6) »BUF2,
1 SXXU93>,BUF3,SXZ(193),BUPt»SZZU93),BUF5«
2 S X ( 1 9 3 ) tBUF6,SZ(193) ,BUF7,RO(193),BUF8,«K193),BUF9»
3 C(193),BUF10,D(193),BUF11.G1(193),BUF12,S2(193),
"> STRIP,Pl,P2,P3,BETA,F«,[Link],KR,DG t IS,JS,K6,NS
,SPftl,G20(2b),SPA2.G30(26),SPA3»GiVO(26),
N = NO
Jl = 2
IF (FMACH .GE .1 ) Jl = 3
C(Jl-l) = 0.
D(J1-1) = 0.
S s 1.
II = 1
I = 12 +1
DO 12 J=2,KY
R0( J) = 1.
*1<J) = 1.
GKJ) = Gl
1Z G2tJ) = 62
=
21 IP I +11
-II
J2 KY
— LT .[Link]) 02 = NY
IF (IVd,K).
LV S IA! S ( 1 -IABS(IV(I«K»)
RO(KY) = AHINO(UViIABS(IV(I«K) ))
RHKY) = UV
AC = T1*A1(I)*CKK)
00 32 J=JltJ2
AB = -T1*A1(I)*B1(J)
BC = -T1*B1( J ) * C 1 ( K )
YP = SO(I,K) *30(J)
A = 1. -RO(J) +AO(I)*AO(I) +YP*VP
H = RO(JJ/A
FH = RO(J)*A
OGI = S*(6(IP«JtL)
OGJ = Gd.J+l.L) -G(I.J-l.L)
OGK = GU.J.t+l) -GKl(I.J)
OGII = G(I+ltJ,L) -G«a,J»L) -GlI.J.L) +GII-1.J.L)
+A3(I)*OGI
OGJJ = G(IiJ+l.L) -G(ItJtL) +G(I.J-1.L)
-B3(J)*OGJ
OGKK = G<I,JtL+l) -Gd.J.L) -G(I,J,D +G(I,J,L-1)
+C3(K)*DGK
OGIJ
OGJK
1 -G(ItJ+l<L-l) *G(I.J-l.L.l)
SY = -B1(J)*UGJ
U = A1(I)*OGI -SX(I)«GY +CA*AO(I) +SA*YP
V = 6Y +SA*AO(I) -CA*YP
W = RO(J)*(Cl(i<)*DGK -SZ(I)*sY +SYAW)
QXY = H*(U*U +V*V)
QQ = QXY + tf*W
AA ='DIM(AAO«.2»UQ)
HZ = FH*SZ(I)
F =1. +SX(I)*SX(I> +SZ(I)*HZ
AV = v -U*SX(I) -W*HZ
UU = H*U*U
VV = H*AV*AV
WW = FH*W*W
UV = H*U*AV
VW = AV*W
UW = U*W
AXX = R1(J)*(AA -UU)
AZZ = FH*AA -WW
R = .(AXX*SXX(I) +AZZ*SZZ(I) .(UW +UW)*SXZ(I))*SY
1 -T1*H*(CA»(U«U -V*V) 4>(sA +SA)*U*V
2 -QXr*(U*AO(D *V*yP))
AXT a A8S(U*A1(I»
AYT = ABS(AV*B1(J))
AZT x ABS(FH*W*C1(K>)
A = R0(J)*BETA*AA/AMAXl(AXT.A r r«AZT t (1. -RO(J)))
AXT = A*AXT
AYT = A*AYT
AZT = A*AZT
IF ([Link]) GO TO 33
AXX = AXX*A2(I)
AYY = (F*AA "VV)*B2(J)
AZZ = AZZ*C2(K)
AXY = -R1(J)*(AA*SX(I) +UV)*(Ag +AB)
AYZ = -(AA*HZ +(/W)*(BC +BC)
AXZ = -UW*(AC +AC)
8P a AYY
BM r AYY
B = -AYY -AYY -Q1*<AXX +AZz)
R = AXX*OGII tAYY«OGJJ +AZZ*OGKK
1 +AXY*OGIJ +AYZ*DGJK +AX^*DGIK +R
SO TO 35
33 NS = NS +1
AXX = UU*A2(I)
AYY a VV*B2(J)
AZZ = UW*C2(K)
AXY = 8.*S*UV»AB
AYZ = 8.*VW*3C
AXZ = 8.*S*U^*AC
BXX = (QQ -JU)*A2(I)
8YY = (F*OQ -VV)*B2(J)
6ZZ = (FH*QQ -MW)*C2(K)
BXY = .(QQ*SX<I> +UV)*(AB +AB)
8YZ = >(QQ*HZ +V«)*(BC +BC)
267
FH = AO(I)*AO(I) «-SO(I.K)*SO(i»K)
V = BKKY)*<1. «-SX(I)*SX(I) +FH*SZ(I)*SZ(I)I
[Link]+ltL> = e<I«KY-l,L>
-<CA*SU(I«K) -SA*AO(J) +U*SX(D
IF [Link])) 50 TO 61
M =. NX +2 -I
C = G(«»KY,L) -G(I,KY»L)
NO = NO *1
E0( NO) = EO(NO) +P3*(E -EO(NO))
N = NO
GO TO 61
51 IF ([Link]) GO TO 61
E = 0.
IF UVU,K).NE.l) 60 TO 57
zz = Z -TYAW*AO(I)*AOU)
53 IF <zz
.GE .ZOCN-1M GO TO 55
N = N -1
GO TO 53
55 R = (ZZ -ZO(N-1))/(ZO(N) -ZO(N-D)
£ = R*£0(N ) +(1. -R)*££MN-1)
57 n = NX +2 -I
G( I f Ky +1 ? L) = B|M,KY-l.L) -E
G(M,KY+lt L) = G(I,i<Y-l,k) +E
GK2(H,KY) = GKK [Link])
GKKM,KY) = <J<M«KY«L)
G(H f KYtL) = lid,KY,L) •«•£
61 IF ([Link]|GO TO 71
IF ll. £Q. 2) RETURN
I =• i +11
GO TO 21
71 S =. -i.
II = -i
I = ii -i
00 72 J=2 »KY
GKJ) = G30(J)
72 G2(J) = GtO(J)
GO TO 21
END
b IV<193.33)«SEP17,ir£i(33).SEPl8,ITE2(33)tSEP19.
6
7
OIHENSION SV <1> t SM(1),CP<1> t X <1),Y <i)
ITEl AND ITE2 ARE LOWER AND UPPER TRAILING EDGE POINTS
269
11 = ITEl(K)
12 = IT£2(K)
SURFACE INDEX IS NY +1
J = NY *1
01 = ,2*FMACH**2
Tl = l./(.7*FMACH**2)
00 12 I=I1»I2
DETERMINE MULTIPLIER H OF SQUARE ROOT TRANSFORMATION
YP = SO(I.K) 4-30(0)
H - AO(I)*AO(I) +YP*YP
DETERMINE SLOPES SX«SZ OF SHEAKIN6 TRANSFORMATION
DSI = SO(I+1<K> -SO(I-1«K)
DSK = SO(I,K+1) -SO(IiK-l)
SX = A1(I)*OSI
SZ = Cim*OSK
DETERMINE FIRST DIFFERENCES OF POTENTIAL
OGI = G(I+ltJtL) -GU-1«J,U
DGJ = G(I«J+1«L> -GdtJ-ltL)
OGK = GdtJ«L+l) -6(I,J«L-1)
DETERMINE VELOCITY COMPONENTS U,tf,H
U = A1(I)*OGI +SX*B1(J)*OGJ +CA*AQ(I> +SA*YP
v = -Bi{j)*uiiJ +SA*Aom -CA*YP
W = Cl(K)*OliK +SZ*B1( J)«OGJ -t-SYAri
DETERMINE SURFACE SPEED [Link] NUMBER [Link] COtFFICIENT CP
QQ =0.
IF ([Link].l.E-6) QQ = (U*U +V*V)/H +w*W
0 = SQRT(QQ)
IF ([Link].O.) Q = -Q
SV(I) = Q
QQ = 1. +Q1*<1. -QQ)
SM(I) = FMACH*Q/SQRT(QQ)
CP(I) = T1*(QQ**3.5 -1.)
DETERMINE SURFACE COORDINATES X,Y
X(I) = .5*SCAL*(AO(I)«*2 -SO<IiK)**2)
12 Y(I) = SCAL*AO(I)«SO(I.K)
RETURN
END
SET K PROPORTIONAL TO CP
K = 30.*(CPO -CP(D) +<t.5
SET ELEMENT K OF LINE TO + SYMBOL
LINE(K) = KODE<2)
WRITE UWRIT.610) X I I ) « T ( I ) t C P C I ) t L I N t
22 LINE(K) = KODE(l)
RETURN
610 FORHAT(3FlO.<ttlOOAl)
END
SUBROUTINE TOTFOR(KTEliKT£2«CHOROiSCi.>SCO»SCM»[Link],XM,
1 CLtCO»CrtP«Ci"IR«CMY)
C CALCULATES TOTAL LIFTtORAG AND MOMENT COEFFICIENTS
C IN DIRECTION NORMAL TO LEAOIN6 EOSE
C BY TRAPEZOlOftL INTEGRATION OF SECTION FORCE COtFFICIENTS
C SPANWISE FORcE IS «OT CALCULATED
C CMP IS PITCHING MOMENT COEFFICIENT REFERRED TO MEAN CHORD
C CMR IS ROLLING MOMENT COEFFICIENT REFERRED TO SEMI-SPAN
c CMY is YAWING MOMENT COEFFICIENT REFERRED TO SEMI-SPAN
DIMENSION CHORD<l)lSCL(l)fSCU(l).SC)*|{l)tCO(l)
SPAN = SCALZ*<CO(*TE2> -CO(KTEl))
271
CL = 0.
CO = 0.
CMP = 0.
CMR = 0.
CMY = 0.
S = 0.
N s KTE2 -1
00 12 KcKTEl .N
oz s .5*SCALZ*(CO<K+1)
.5*51 -CO(K))
z = .5*SCALZ*(CO<K+1)
.5*Si +CO(K»
CL r CU *SCL(K)*CHORO( K) )
CO = CO +LlZ*(SCO(K+l)*CHO«0(K*l) + SCO(K) »CHOKO( KJ )
CMP r CMP +OZ*<SCM(K+1)*CHOKO(K+1)**2 -fSCM(K)*CHORD(K)**2)
CMR = CMR +Z*DZ*(SCL(K+l)*CHORp<K-H) +SCU(K)«CHORD(K) )
CMY = CMY +Z«OZ*(SCO(K + l.)*CHORr)tK + l> +SCDIK) *CHORD(K) )
12 S = S *i *OZ*(CHOHO((<+1) +CHORO(K»)
CL s CL/S
CO S CO/S
CMP S CMP*
CMR =
(C1R *CMR)/(S*SPAN)
CMY = (COY
*CMYI/(S*SPAN)
RETURN
ENO
SUBROUTINE CHARTY
GENERATtS MACH NO CHARTS IN PLANE OF. uING PLANFORM
COMMON G(l93«26.<f).SEPl.
AOU93)«SEP2.A1U93),SEP3.A2(193) .SEPt. A3( 193) ,SEP5.
80(26)iSEP6i81(26)<SEp7.B2(26).SEpa«83(26),SEP9<
CO(33),SEP10,C1(33),SEP11.C2(33),SEP12»C3(33),SEP13.
SO(193.33),SEPm,EO(129),sEP15,ZO(129),SEPl6,
IV(l93i33),SEPl7,IT£i(a3).SEPlB.ITt2(33),SEPi9.
NXtNY.'NZiKTEl,KTE2tKSYM,[Link]
YAW,CYAW«SYAW,ALPHA,[Link]«Nl«N2<N3tIO
DIMENSION LV(33)
IWRIT = 6
LX = NX/2 +1
MX = NX +1
KY a NY +1
MY = NY +2
01 = ,2*FMACH**2
00 2 K=2,NZ
LV(K) = HY
IF (IV(LX.K).LE.O) LV«) = KY
2 CONTINUE
WRITE (ItlRlT.12)
12 FORMAT«t2HOUPP£R SUKFACE MACH NO CHART IN WING PLANE)
LI = 1
IM = NX
21 00 22 L=l,3
BUFFER IN <N1«1) (6(1,l.L)»6(MX,MY.L))
IF (UNIT(Nl).GT.O.) GO TO 101
272
22 CONTINUE
K a 1
31 K a K +1
L a LV(K)
N = 1
II a NX/2 +1
KI 0 a
JJ 2 a
KJ 1 a
33 I II =
J a JJ
IF ([Link].L) GO TO 35
J s NY
35 YP = SO(I.K) +30(J>
H s AO(I)*AO(I) +YP*YP
OSI a SO(I+ltK) -SO(I-l.K)
OSK s SO(I,K+1) -SO(I,K-1)
sx a A1(I)*OSI
sz a C1(K)*OSK
OGI s G(I+l»J«2t -G(I-1,J,2)
OGJ a GIItJ+1.2) -G(I,J-1,2)
OGK a 6(1, J, 3) .6(1, J,l»
IF ([Link].L) GO TO 37
H s NX +2 -I
OGI = .5*1061 *S(M-1,J«2) -G(vi+l,j,2) )
OGJ a ,5*(OGJ +3(M,J-1«2) -G(v|,j+i,2) )
OGK s ,5*(OGK +S(M,J,3) -G(M,j,l))
37 U a A1(I)*OGI *SX*Bl( J)*OGJ +CA«AO(I) +SA*YP
V = -B1(J)*OGJ +SA*AO(I) -CA*YP
W s C1(K)*OGK +SZ*BK«J)*OGJ +SYAW
QQ a (U*U +tf*V)/H *W»«
F = 1. +Q1*(1. -au)
Q a 0.
IF ([Link]. 0.) Q = SQRT(90/F)
IF (LI*(I -II) +J -JJ) 41,45,43
41 QO s Q
I a I +LI
J s KY
60 TO 35
43 Q a ,5*(Q +QQ)
L a HY
45 N a
N +1
IVIN.K) s
100.*F,1ACH*a
IF (II. EQ.!« ) GO TO 51
IF ([Link] .KY ) GO TO 47
KI a LI
KJ a 0
47 II a II +KI
JJ a JJ +KJ
GO TO 33
51 IF ([Link]) GO TO 61
DO 52 1=1 ,MX
00 52 J=l ,MY
52 6 ( 1 , J , 2 ) = 6(1,J,3)
273
SUBROUTINE REFIN
C INTERPOLATES POTENTIAL AT MESH POINTS OF REFINED GRID
c SUCCESSIVE SLICES OF THE G ARRAY .REPRESENTING VALUES OF POTENTIAL
c ON x-r PLANES AT SUCCESSIVE VALUES OF Z,ARE READ
C FROM ONE DISC FILE,UPDATED,AND WRITTEN ON A SECOND DISC FILE
COMMON 6(193,26,[Link],
1 AOU93) «SEp2,Al(193),SEp3,A2(193) ,SEpf ,A3(195) ,SEp5,
2 B 0 ( 2 6 ) ,SEP6,81(26),SEP7,82(26).SEPa,83(26),SEP9,
i CO(33),SEP10,C1(33),SEP11,C2(33),SEP12,C3(33),SEP13,
«* SO<193,33),S£Pl'»,EU<l29),sePl5.ZOU29),S£Pi6,
5 IVU93i33)tSEP17«ITEl(33).SEP18,ITE2<33),SEP19,
6 NX,NY»N2,KT£i,KTE2»KSYMiScAL,SCALZ»
7 YAWiCYAW,SYAWiALPHH*CA«[Link],Nl»N2,N3,10
MX = NX +1
KY = NY +1
MY = NY +2
HZ = NZ +1
MXO = NX/2 +1
HYO = NY/2 +2
HZO = NZ/2 +1
K =1
C INTERPOLftTt POTENTIAL ARRAY G
C READ SLICE OF POTENTIAL ARRAY FROM FIRST DISC FILE
11 BUFFER IN (NI.D (G(i,i«i),G(Mxo,mo,i))
c GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(Nl).GT.O.) GO TO fOl
c SHIFT I/ALUCS TO LOCATIONS IN NEW GRID
J = NY/2 +1
JJ = KY
21 I = MXO
II = MX
274
c INTERPOLATE IN Y
00 52 1=1, MX
00 54 J=2,NY,2
54 G(ItJfl) = ,5*(G( It J+l»l) +G(I» J-l«i»
52 G(I*MY«1) = 0.
c WRITE SLICE OF INTERPOLATED POTENTIAL ARRAY ON SECOND DISC FILE
BUFFER OUT<N2,1» (5(1, 1 .1 » ,6<MX,MY,i ) |
C GIVE UP IN EVENT OF DISC FAILU«E
IF (UNIT(N2).6T.O.) GO TO 401
C INCREMENT Z
K = K +1
IF ([Link]) SO TO 11
REMIND Nl
REWIND <M2
C REAO FIRST TWO SLICES OF POTENTIAL ARRAY FROM SECOND DISC FILE
BUFFER INI (N2,i) (Gti,i,i) ,G(MX,MY,D >
C GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N2).GT.O.) GO TO *01
BUFFER IN (N2«l) < G ( 1 , 1 , 3 ) , 6 ( MX , MY , 3 ) )
C SIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N2).GT.O.) GO TO 401
C WRITE FIRST SLICE OF POTENTIAL ARRAY ON FIRST DISC FILE
BUFFER OUT(N1,1) ( S( 1, 1,1) »G(HX,MY,1 ) )
C GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(Nl).GT.o.) 60 TO fOl
K =1
C INCREMENT 2
111 K = K +1
C INTERPOLATE IN Z
00 112 J=1,MY
00 112 1=1, MX
112 G(I,J,2) = .5*(G(IiJ«l> +G(I«J,3))
C WRITE TWO SLICES OF INTERPOLATED POTENTIAL ARRAY
C ON FIRST QlSc FILE
00 122 L=2,3
BUFFER OUT(Ni.l) (G( 1, 1»L) ,G(HX,HY,L) >
C GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(NI) .GT.O. ) GO TO 401
122 CONTINUE
IF ([Link]) GO TO 201
c SHIFT SLICES OF POTENTIAL ARRAT
oo 132 J=I,MY
00 132 I=1«MX
132 G(I,J,1) = G(I,J,3)
275
IF ([Link].1XO) 60 TO 221
I - II
231 I =1 -1
E =0.
IF (IV(I.K).NE.i) SO TO 237
C IV = 1 INDICATES VORTEX POINT
c INTERPOLATE JUMP E IN POTENTIAL
c TO SET POTENTIAL AT ou-iflY POINT BELOW VORTEX SHEET
ZZ = Z -TYAH*AO(I)**2
233 IF ([Link](N-l)) GO TO 235
N = N -1
SO TO 233
233 R = (22 -ZO(N-1))/<ZO<N) -ZO<N-1M
E s R*EO(N) *U. -R)*EO(N-l!
237 n = NX +2 -I
G<I,KY+1,2) = G(H,KY-1.2) -E
G(M«KY+1»2) = G(I,KY-1,2) +E
IF (IV(I.K).NEt-l) 60 TO 241
C IV = -1 INDICATES POINT JUST BEYOND EflGE OF WING OR, VORTEX SHEET
c KENORMALIZE POTENTIAL ON EITHEK SIDE OF CUT AT MEAN VALUE
G(I,KY«2> = .5*6(1,KY»1) + .25K 6( I ,KY,3) +G(PI,KT,3) )
IF <IVU,K+1).LT.1)
1S(I,KY»2) s ,5*G(ItKY«3) +.25*{6(I,KY,1) +<»<M«KY,1))
G{rt,KY»2) = 6(I,KY,2)
G(IiKY-l«2) = .5*(G(I,KY,2) +6(I,KY-2«2))
G(M,KY-1,2) = .5*<G(M,i<Y,2) +b(M,KY-a»2) )
2H1 IF ([Link].2) 60 TO 231
IF ([Link]) GO TO 261
c SHIFT SLICES OF POTENTIAL ARRA*
00 252 J=l«f"lY
00 252 i=l,HX
G(IiJ,l) a Gdi J.2)
252 6(1,J,2) = G(I,J,3)
C WRITE SLICE OF UPDATED POTENTIAL ARRAY ON SECOND DISC FILE
BUFFER ouT(N2,i> (GU,I,I),G<MX,MY,I) ,
c GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N2).GT.O.) GO TO 401
C READ SLICE OF POTENTIAL ARRAY FROM FIRST DISC FILE
BUFFER IN INI,D (S(i,i,3),ti([Link]«3»
c GIVE UP IN EVENT OF DISC FAILURE
IF (UNlTdMl).GT.O.) GO TO 401
C INCREMENT Z
K = K +1
GO TO 211
261 EOCNO + D = 0.
C WRITE LAST TWO SLICES OF POTENTIAL ARRAY ON SECOND DISC FILE
00 262 L=2,3
BUFFER OUT(N2,1) (6(1,1«L)«G<HX,HY,L)j
C GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N2).GT.O.) GO TO 401
262 CONTINUE
REWIND NI
REMIND N2
C COPY FINAL VALUES OF POTENTIAL ON FIRgT DISC FILE
00 502 K=[Link]
277
SUBROUTINE SflOO
C SMOOTHS POTENTIAL
C BY REPLACING THE VALUE AT EACH POINT gY A WEIGHTED AVERAGE
C OF THE VALUES AT NEIGHBOURING POINTS
C SUCCESSIVE SLICES OF THE G ARRAY .REPRESENTING VALUES OF POTENTIAL
c ON X-Y PLANES AT SUCCESSIVE VALUES OF [Link] READ
C FROM ONE DISC FILE,UPDATED*AND WRITTEM ON A SECOND DISC FILE
COMMON 6(193.26,t).S£Pl,
1 AO(193),SEP2,A1<193),SEP3,A2(193),SEP<»,A3U93),SEP5,
2 BO(2&).SEPS,Bl<26).SEP7,Ba(26),SEP8.B3(2&),SEP9»
i CO(33)tSEPlOiCl(33)iS£Pli.C2(33),St.Pl2. C3( 33),SEPl3.
4 80(193.33)«[Link](129),sEP15«ZO(129).SEK1&,
5 IV(193,33),SEpl7,irEl(33).SEplB,ITE2<33),SEP19.
6 [Link],Ni:,KTEl,KTE2tKSYM,ScAL,SCALZ«
7 YAW,CYAU,SYAW,ALPHA,CA,SA,FnACH,NltN2tN3,IO
MX = NX +1
KY = NY +1
MY = NY +2
«Z = NZ +1
c SET SMOOTHING PARAMETERS
PX = l,/6.
PY = 1./6.
PZ = 1./6.
C READ FIRST THREE SLICES OF POTENTIAL ARRAY FROM FIRST DISC FILE
00 2 L=l,3
BUFFER IN (Nl«l) (S(l.l.L).6(HX,MY.L))
c GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N1).6T.O.) GO TO 5l
2 CONTINUE
C WRITE FIRST SUlCE OF POTENTIAL ARRAY ON SECOND DISC FILE
BUFFER OUT(N2d> (G( 1« 111) tG(MX< MYtl) )
C GIVE UP IN EVENT OF DISC FAILURE
IF <UNir<M2).ST.O.) GO TO 51
K =1
C INCREMENT Z
UK = K +1
C GENERATt SMOOTHED VALUES OF POTENTIAL FOR MIDDLE SLICE
278
00 12 J=3,NY
00 11 1=2,NX
1* 6(IiJ,f) = (1. -PX -PY -PZ)»G(I,J.2)
1 +.5*PX*(G(I+i<Ji2> 4-6(1-1,J,2))
2 +.5*PY*(6(ItJ+l»2) +G(I,j-l»2M
3 +.5*PZ*<G(ItJ«3) +G(ItJ«t»
6(1,J,1) = G(1,J,2)
12 6(MX,J,1) = G<nx,J«2)
LEAVE BOUNDARY VALUES UNCHANGED
00 16 1=1,MX
6(1.1,<U = G(I,1,2)
6(1,2,11 = G(I,2,2)
G(IiKY,4) = G(I,KY,2)
lb G<[Link]«<n = [Link].2)
WRITE SLICE OF UPOATEO POTENTIAL ARRAY ON SECOND DISC FILE
BUFFER OUT(N2,1) < G( 1, !.<*) .G(nX,HY»t) }
GIVE UP IN EVENT OF DISC FAILURE
IF <UNIT(N2).GT.O.> GO TO 51
IF ([Link]) GO TO 31
SHIFT SLICES OF POTENTIAL ARRAY
00 22 J=1,MY
00 22 1=1,MX
G(ItJtl) = G(I,J,2)
22 6(1,J,2) = G(I,J,3)
READ SLICE OF POTENTIAL ARRAY fROH FIRST DISC FILE
SUFFER IN <NI,I> (S'(iii«3)«G(nxtnrt3>i
GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(Nl).GT.O.) GO TO 51
GO TO 11
WRITE LAST SLICE OF UPOATEO POTENTIAL ARRAY ON SECOND DISC FILE
31 BUFFER OUT(N2,1) (G{1,1«3),G<MX,MY,3)>
GIVE UP IN EVENT OF UlSC FAlLUKE
IF (UNIT(N2).GT.Q.) 60 TO 51
REMIND NI
REWIND N2
COPY FINAL VALUES OF POTENTIAL ON FIRsT DISC FILE
00 42 K=1,MZ
BUFFER IN (N2«l) (G(1,1,1),G(MX,NY,1))
GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(N2).GT.O.) GO TO 51
BUFFER OUT(N1,1) (G(1,1,1),G(WX,«Y,1))
GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(NI).GT.O. ) GO TO 51
i+2 CONTINUE
SET 10 TO INDICATE SUCCESSFUL COMPLETION
10 =1
RETURN
SET 10 TO INDICATE DISC FAILURE
51 10 - 0
RETURN
ENO
279
SUBROUTINE SPLIF(«»N»S,F,FP,[Link],[Link],VN»MOD£»FOM,INO)
C CUBIC SPLINE
C SPLINE IS FITTED TO DATA ARRAY F AT NoUES S
C FROM INDEX H TO INOEX ,M
C KM =• i»2 OR 3 INDICATES THAT FIRST,SECOND OR THIRD DERIVATIVE
C IS GIVEN VALUE VM AT POINT l"l
C KN = j.,2 OR 3 INDICATES THAT FIRST,SECOND OR THIRD OEKIVATIVE
C IS GIVEN VALUE VN AT POINT N
C IF MODE = 0 NOOAL VALUES OF FIRST,SECOND AND THIRD DERIVATIVES
C OF SPLINE ARE STORED IM FP,FPP AND FPPP ARRAYS
C SO THAT FITTED VALUE AT A DISTANCE H BEYOND A NODE IS
C F +FP*H +FPP*H**2/2. +FPPP*H**3/6
C IF HOOt GT 0 FPPP IS GIVEN THE NOOAL VALUES OF THE INTEGRAL OF F
C INSTEAD OF ITS THIRD DERIVATIVE, STARTING KITH THE VALUE FOR
C THEN THE THIRD DERIVATIVE CAN tfE RECOVERED AS
C (FPPII+l) -FP(I))/(S(I+1) -S(I>)
C INO IS SET EQUAL TO 0 IF S IS NOT A MONOTONE ARRAY
DIMENSION S(l),FU),FP(l),FPPm,FPpP(l)
INO = 0
K = IABS(N -M)
IF (K -1) 61,81,1
IK = (N -M)/K
I = M
J = M +K
OS = S(J) -S(I)
0 = OS
IF (OS) [Link]
11 OF - (F< J) -F<I) )/OS
IF (KM -2) 12,13,11
12 U = .5
V = 3.*(DF -Vfl)/DS
GO TO 25
13 U = 0.
V = VM
GO TO 25
1H U = -1.
V = -DS*VM
GO TO 25
21 I = J
J = J +K
OS = S(J) -S(I)
IF (D*OS) 81,81,23
23 OF = (F(J) -F(I»/OS
B s l./(OS +DS +U)
U = B*OS
V = B*(fa.*OF .V)
25 FPU) = U
FPP(I) = V
U = (2. -U)*OS
V = 6.*DF +DS*V
IF (J -N) 21,31,21
31 IF (KN -2) 32,33,34
32 V = (6.*VN -VI/U
GO TO 35
33 V = VN
280
GO TO 3S
3t V = (OS*VN +FPPJI))/(1. +FPfI)
35 8 = V
0 = OS
m os r S(J) -S(I)
u = FPP(I) -FP(I)*V
FPPP(J) = (V -U)/OS
FPP(I) = U
FP(I) s (F(J) -F(I))/OS »OS*(V +U
V = U
J = I
I = I -K
IF (J -H) <t 1,51, <fl
51 I = N -K
FPPP(Nl = FPPP(I)
FPP(N) = B
FP(N) = DF .»-D*(FPP«I) +B +8)/6.
INO - 1
IF (MODE) 81 ,81,61
61 FPPP(J) = FOM
V = FPP(J)
71 I s J
J = J +K
US S S(J) -S(I)
U S FPP<J)
FPPP( J) FPPP(I) +.5*OS*(F(I) *FfJ) -DS*OS*(U
-
V S U
IF (J -N) 71,61,71
ei RETURN
ENO
SUaROUTlNE INTPL(«l»NI,Sl,Fl»fl»N»S»ff»FP»FPP,FPPP,l'100E)
c INTERPOLATION USING PIECEWISE TAYLOR SERIES
C AS GENERATED BY CU3IC SPLINE OR ITS IMTEGRAL
c VALUES F,FP»FPP ANO FPPP OF FUNCTION ANO ITS FIRST,SECOND
C ANO THIKO DERIVATIVES A HE GIVEN AT NODES S
C FROM INDEX « TO INDEX ,M
C INTERPOLATED VALUES FI ARE GENERATED AT POINTS SI
C FROM INDEX «I TO INDEX NI
C IF MOOt 6T 0 A CORRECTION IS ADDED
C FOR A PIECEWISE CONSTANT FOURTH DERIVATIVE
C SO THAT INTEGRAL OF CUdIC SPLINE IS EVALUATED EXACTLY
DIMENSION SI(l),FI(l),S(l)»F<l),FP(i),FPP(l),FPPP(l)
K = IA3S(N -M)
K = <N -flJ/K
I = M
WIN = MI
NIN = NI
0 = SIN) -SIM)
IF (D*(SI(NI> -SI(MIM) 11,13,13
11 MIN = NI
NIN = MI
281
19 KI = IABS(NIN -HIN)
IF ( K I ) 21.21,15
15 KI a (NIN -MIN)/KI
21 II = HIM .KI
C =0.
IF ( M O D E ) 31,31,23
24 C =1.
31 II = II +KI
SS = SKID
33 I = I +K
IF (I -N) 35.37,35
35 IF ( D * ( S ( I ) - S S ) ) 33*33,37
37 J = I
I = I -K
SS a SS - S < I )
FPPPP a C*(FPPP(J) -FPPP(IM/(S(j) -S(I»
FF = FPPP(I) +.25*SS*FPPPP
FF a- FPP(U +SS*FF/3.
FF a FP(I) +.5*SS*FF
FI(II) a F(I) +SS*FF
IF (II - N I N ) 3l,m,31
<tl R E T U R N
END
SX = 2.
TX = 3.5
ST = 2.75
OY « 8./NZ
IF [Link].O) 60 TO 1
c INITIALIZE: PLOTTER IF IPLOT LT o
CAUL PLOTSBLdOOO,2HHANTONY JAMESON X109H03)
C DEFINE ORIGIN
CALL PLOT([Link]-3)
in =i
C WRITE TITLE FOR DRAWING OF WINS
£NCOD£U2«2»Rd) )
2 FORNATd2HVl£W OF WINS)
CALL SYMBOL(2.».5..14,R,0.»12)
C READ FIRST THREE SLICES OF POTENTIAL ARRAY FROM DISC FILE
11 00 12 Lsl,3
BUFFER IN (Ni.l) (S(l«l«L)«G(MX,nY,L))
c GIVE UP IN EVENT OF DISC FAILURE
IF (UNIT(Nl).GT.O.) GO TO 101
12 CONTINUE
K =2
C INCREMENT Z
21 K = K +1
IF ([Link].KTE2) GO TO 61
C SHIFT SLICES OF POTENTIAL ARRAY
DO 22 J=1,MY •
00 22 1=1,MX
S(I.J.l) = G(I»J,2)
22 G(I,J,2) = 6(1,J,3)
C READ SLICE OF POTENTIAL ARRAY FROM DISC FILE
BUFFER IN (Nl,l) <Gd,l,3),G(MX,«Y,3))
IF (UNIT(Nl).GT.O. ) GO TO 101
c GIVE UP IN EVENT OF oisc FAILURE
IF ([Link]) GO TO 21
11 = ITEl(K)
12 = ITE2(K)
C CALCULATE SURFACE SPEEO [Link]* NUMBER SM,PRESSURE COEFFICIENT CP
C ANO COORDINATES X,Y OF WING SECTION
CALL VELO <K,2,[Link]"I,CP«X,Y)
IF ([Link]) GO TO HI
C WRITE TITLE ANO FLOW PARAMETERS
C BEFORE GENERATING PLOT AT FIRST SPAN STATION
ENCODERS,32,R<1)) TITLE
32 FORMATd2A4)
CALL SYMBOL!.5,0.•.14»*«0.,48)
ENCODE(4o,34,Rd)> FMACH,YA,AL
31* FORMAT(4HM = ,F6.3,3X,6HYAW = ,F5.2,3xtSHALP = ,F5.2)
CALL SYM80L(,5,-.25«.14,R,0.,40)
£NCOOE(tO«36,R(l)) VLO,CL,CO
36 FORMAT(6HL/0 = ,F6.2.3X«5HCL = ,F6.t,^[Link] = ,F6.4)
CALL SYMBOH ,5,-,5» .14, R, 0. ,i»0 )
C SCALE ANO TRANSLATE COORDINATES ANO PRESSURE COEFFICIENT
m XPIIN = x(Il)
00 42 1=11,12
t»2 XI«IIN = AMINKX(I) ,XMIN)
284
SCAUX = 2.5/CHOROO
SCALP B -1.25
DO 44 I=IltI2
X(I) = <X(I) -XMIN)*SCALX *SX
Y(I) = Y(I)*SCALX +ST
44 CPU) = SCALP*CP(I) +SY
C INCREMENT VERTICAL SHIFT FOR NEXT SPAN STATION
sr = ST +OT
IF (M.E8.2) 60 TO 51
.C IF H = 1 DRAW WINS SECTION
N =12 -II +1
CALL LINE<X(II),Y(II)»N.I.O,I,O.,I.,O.«I.)
60 TO 21
C IF H = 2' PLOT PRESSURE COtFFIClENT OVER UPPER AND LOHER SURFACES
51 N =12 -LX 4-1
c PLOT UPPER SURFACE COEFFICIENT AT LEFT SIDE OF PAGE
CALL LlNE(X<LX),CP(LX)tNtltO,[Link].,[Link].)
N = tX .11 4-1
c TRANSLATE x COORDINATES ro RIGHT
00 52 1=[Link]
52 X(I) = X(I) +TX
c PLOT LOWER SURFACE COEFFICIENT AT RIGHT SIDE OF PAGE
CALL LINE<X(II).cPdi>.N. 1.0,1,o.,i.,o.,i.)
60 TO 21
61 REWINO Nl
H = M +1
c SHIFT ORIGIN FOR NEXT PLOT
CALL pLOT<i2.,o.,-j)
IF ([Link].2) 60 TO 71
C RESET HORIZONTAL AND VERTICAL SHIFTS
SX = 0.
SY = 2.75
c WRITE TITLES FOR PRESSURE PLOTS
ENCODEU<»t62,R(l))
62 FORHATC24HUPPER SURFACE PRESSURE )
CALL [Link].m.R.o..24)
ENCOOE(24.64,R(1))
64 FOR«AT(24HLOWER SURFACE PRESSURE )
CALL SYMBOL(3.5t.5t.l4tRt0.t24)
60 TO 11
C SET 10 TO INDICATE SUCCESSFUL COMPLETION
71 10 =1
RETURN
C SET 10 TO INDICATE DISC' FAILURE
101 10 =0
RETURN
END
285
the exact form because its computation time is about forty percent
longer than the listed option. The option is based on a centered
C ****TO USE THIS OPTION REPLACE THE SUBROUTINE 1UR«AN FOUND ****
c ****ON PAGES ^ THRU 9 OF THE LISTING OF PROGRAM H BY THE ****
C ****FOLLOWlNG NEW VERSION. ****
SUBROUTINE NlURMAN
C SET UP COEFFICIENT ARRAYS FOR THE TRIDIAGONAL SYSTEM USEO FOR LINE
C RELAXATION aND COMPUTE THE UPDATED PHI ON THIS LINE
COMMON P H I ( 1 6 2 i 3 1 ) t F P ( 1 6 2 i 3 1 ) « A ( 3 1 ) • B< 31) , C < 31) 10( 31) .E (31)
1 t R P ( 3 1 ) ,RPP ( 31) , R ( 3 l ) , R S ( 3 l ) , R I ( 3 1 ) , A A ( 1 6 2 ) , B 9 ( 1 & 2 ) . C 0 ( l 6 2 )
? .SId62)«PHIRd&2>, X C d & 2 ) . Y C U & 2 ) , FMd&8)>ARCL(1&2).OSUMds2)
3 «ftNGOLD(16?),XOLD(162).YOLO(162).ARCOLD(l&2),OELOLO(1&2)
C0.1M3N /A/ [Link]«[Link]«[Link]»[Link],YR
1 *XA,YA,TE,nT.O*,DELTH,OELR,RAtDCN,DSN,RAtt,EpSIL,[Link],Ca
2 «C«»,C5.C6,C7«BET«BETA«FSY,"[Link](!*) [Link]
3 <lK,JK,IZ,lTYP,!U|Of}[Link]<NFCtNCYtNRNtl\|GtIOlntN2tN3«Nt>MTtIXX
H , NPTS,LL,i,LSEP,>H*
DIMENSION VU(35),RPO(35)
DATA RPO/35*0./
BETP = BETA+.25
C 00 THE BOUNDARY
RP(1) = 0.
RP(NN) = 0.
KK = 0
PHIO = PHI(I,2)-2.*DR*CO(I)
PHIyP= PHI(I,2)-PHI(I,1)
PHIYY = PHIYP+PHIO-PHIU.I)
pHIxx = PHI(I+l,l)+PHI(I-ltl,-PHl(I,l)-PHI(I.l)
PHIXM = PHI(I*ltl)-PHl(I-l.l)
PHIxP = PHI(I+1,2)-PHI(I-1,2)
C CHECK FOR THE TAIL POINT
IF <[Link]> GO TO 10
C(l) = (Cl+Cl)*RS(l)
A(l) = -C<1)+XA*C1-C1
OU) = Cl*(PHlXX+RSm*PHlTY+RAt*CO<I)-Ei(l»
GO TO 40
10 U = PHIXM»OELTH-SI(I)
BQ = U/FPU.l)
OS = U*BQ
CS = C1-C2*BS
9Q = BQ*QS*(FP(I-1.1).FP(I'H«1) )
X = RAi»*([Link])*CO(I>
C(l) = (CS+CS>*RS<1>
CMOS = CS-QS
0(1) = CMQS*PHIXX*CS*RS(1)*PHIYY + Rl(l)*BQ + X + RPPrl)
PHIXT = BETP*ABS(U)+A8S(CMaSi
IF ([Link].O.) GO TO 30
C FLOW IS SUPERSONIC. BACKWARD DIFFERENCES
KK = 1
PHIXT s PHIXT-CMOS
A<1) = -(C(i)+PHIXT)
RPP(l) i CM8S*PHIXX
0(1) •= Od)-PHIXT*E<l)-RPP(l)
GO' TO *0
c FLOW SUBCRITICAL, CENTRAL DIFFERENCES
287
50 Ad) = XA*CMQS-Cd)-PHIXT
0(1) = 0 ( 1 ) - P H I X T * E ( 1 )
RPPd) = 0.
00 NON-BOUNDARY POINTS
1*0 DO 60 J =[Link]
PHIXX = PHI(I«-1.J)+PHI(I-1, J)-PHI(I, J ) - P H I ( I . J )
CU = PHIXP
PrilXP = PHld + l«J+l)-PHI<I-ltJ+l)
PHIXY = PHIXP-PHIXI
PHIXM = OU
DU = DU*DELTH
PHIyv| = PHIYP
PHIYP = PHI( ItJ + D-PHl (I. J)
PHIYY = PHIYP-PHIYM
U = R( J)*OU-SHI>
DV = R(J)*(PMl(I,J+i)-PHI(I,J-l))*DELR
V = OV*R(J)-cO(I)
VV(J) = V
RAV = R(J)*RA*V
9Q = l./FPd.J)
BQU = BQ*U
us = sau*u
UV = (BQU*BOU)*V
US = 8Q*V»tf
us = us+vs
CS = C1-C2*QS
CMVS = CS-Vs
CMUs = CS-Us
COMPUTE CONTRIBUTION OF RISHr-HANO SIDE FRO* LOW ORDER TERHS
0(J) =RA t >*(<CMVS+US-VS>*DV-UV*OU)+RI< J)*QS*3Q*(U*(FP(T-l«J)-
UV = .5*BQU*RAV
C(J) = RS( J)*Cf1\/S
B(J) = C(J)
0(J) = 0(J)+C(J)*PHIYY-UV*PHIXY*CNUS*PHIXX+RPP(J)
CSQS = CS/QS
CVQS - CSQS-1.
PHIXT = BETP*ABS{U)
PHIYT = BETP*ABS(RftV)
IF ([Link].O.) GO TO 50
SUPERSONIC FLOW, USE BACKWARD DIFFERENCING
KK = KK+1
PHIxT = PHIXT-C«QS*(US+US+ABS(UV))*CSQS*VS
PHIYT r PHIYT-CNQS*(RS(J)*<VS+VS>+ABS<UV))
B(J) = RSI J)*CSQS*'JS
C(J) = B(U)+PHIYT '
A t J ) = -<C(J)+B(J)+PHIXT)
RPP(J) = CMQS*(US*PHIXX -f Utf*PHIXY)
RP(J) = CHQS*RS(J)*VS*PHIYY-f-RPO(J)
R P O ( j ) = CI*IGS*UV*PHIXY
GO TO 60
SUBSONIC FLOW. USE CENTRAL DIFFERENCES
50 C(J) = C(J)*PHIYT
PHIXT = PHIXT+CMUS
A(J) = XA*C!"IUS-B(J)-C<J)-PHIXT
288
RPP(J) = 0.
RP(J) = RPO<J)
RPO<J) = 0«
60 0(J) = D<J)-PHIXT*E'<J)-RPP<J>-RP<J>
00 70 J = 2.N
IF (VV(J).LT.O.) SO TO 72
D(J) = 0<J)+RP(J+l)
GO TO 70
72 9Q = B(J)
B(J) = C(J)
C(J) = BO
0(J) = 0(J)*RP(J-l)
70 CONTINUE
75 MSP = NSP*KK
SOLVE THE TRIOIASOMAL SYSTEH
CALL TRIO
RETURN
END
289
IF (IABS(NRN).GT,999) FAC = .5
PAGE: 156 INSERT AFTEK LINE 11 THE FOLLOWING
TE = AINAG jeej
PAGE I5b DELETC LINES 31 THRU 33 AND REPLACE BY THE FOLLOWING
NRN = ISIGN(I"IOO(IABS(NRN) ,1000) »NRN)
CALL CPLOT {(3.0,2.0),-3)
C SF WILL BE THE CHORD LENGTH iN INCHES
SF = 5.
C DEFAULT TAIL EXTENSION TO 1 IF [Link].l OTHERWISE TE=H..6*DY
IF ([Link].O.) TE = l.+.6*AHAXl([Link](CC(6)t1.00001-TR))
IF (([Link].1.).OR.(CC(6).LE.O.)) TE = -TE
CALL GRF (NN,NNX»TE)
PAGE 156 DELETE LINES 39 AND ifQ AND REPLACE BY THE FOLLOWING
Xl«IAX = 22.*FAC
CALL CPLOT (CflPLX(.5*XnAX,it.S),-3)
SIZE = .It
REWIND N3
READ < N 3 » 9 0 ) (PG(I),I = 1,6)
CALL CSYPIBL ((-S.O.-[Link].&O)
SIZE = .07
PAGE 15b DELETE LINE 13 AND REPLACE BY THE FOLLOWING
CALL XYAXES ((0. • 0. ),1. + 3.»XflAX/ll..1.*3.*XMAX/11.,1.4/XMAX)
PAGE 156 DELETE LINE "»5 AND KEPLACE BY THE FOLLOWING
CALL XYAXES «0.»0.),l. + XMAX/il.,l.+XMAX/ll.,<*.t/XMAX)
PAGE I5b DELETE LINE ifB AND REPLACE BY THE FOLLOWING
40 FORMAT (3I5,20X,3F5.3,5X,F5.3)
PAGE 160 DELETE LINES l»? AND 50 ANO REPLACE BY THE FOLLOWING
IA = IA3S(II)
00 10 J = [Link]
PAGE 160 DELETE LINE 55 AND REPLACE BY THE FOLLOWING
IF [Link].O) T = T*CSQRT(l.*dP*flP/(T*T)iON£)-6P
PAGE 161 DELETE LINE 1 AND REPLACE BY THE FOLLOWING
JJ = 15 + IA3S(NK)
PAGE i$i DELETE LINE 3 AND REPLACE BY THE FOLLOWING
WRITE (N2.100) (FF(J)iJ =1.64)
PASE 161 DELETE LINE 9 AND REPLACE BY THE FOLLOWING
40 FORMAT (lXA4,2Al,A3,lXA4t9F5.3.1XA4)
PAGE 161 DELETE LINES 14 ANO 15 AND REPLACE BY THE FOLLOWING
90 FORMAT (///38X,6HTAP£ 7///4XA4tAliA3«A4.F6.2«2F5.2.F6.2.2F6.3
1 ,F7.3,F6.3,F5.2,A4/4X,16A4/4X,16A4)
PAGE 161 DELETE LINE 20 AND REPLACE BY THE FOLLOWING
2 F5.3.3X3HDY=F5.3.3X,4HT/C=FS.3/////i5X«14HTAP£ b, PATH O/)
Vol. 77: A Auslender, Problemes de Minimax via I'Analyse Con- Vol. 103: D. E. Boyce, A Farhi, R. Weischedel, Optimal Subset
vexe et les Inegalites Variationelles: Theorie et Algorithmes. VII, Selection. Multiple Regression, Interdependence and Optimal
132 pages. 1972. DM 16,- Network Algorithms. XIII, 187 pages. 1974. DM 20,-
Vol. 78: GI-Gesellschaft fur Informatik e.V. 2. Jahrestagung, Karls- Vol. 104: S. Fujino, A Neo-Keynesian Theory of Inflation and
ruhe, 2.-4. Oktober 1972. Herausgegeben im Auftrag der Gesell- Economic Growth. V, 96 pages. 1974. DM 18,-
schaft fiir Informatik von P. Deussen. XI, 576 Seiten. 1973. Vol. 105: Optimal Control Theory and its Applications. Part I.
DM 36,- Proceedings of the Fourteenth Biennual Seminar of the Canadian
Mathematical Congress. University of Western Ontario, August
Vol. 79: A. Berman, Cones, Matrices and Mathematical Program- 12-25, 1973. Edited by B. J. Kirby. VI, 425 pages. 1974. DM 35,-
ming. V, 96 pages. 1973. DM 16,-
Vol. 106: Optimal Control Theory and its Applications. Part II.
Vol. 80: International Seminar on Trends in Mathematical Model- Proceedings of the Fourteenth Biennial Seminar of the Canadian
ling, Venice, 13-18 December 1971. Edited by N. Hawkes. VI, Mathematical Congress. University of Western Ontario, August
288 pages. 1973. DM 24,- 12-25, 1973. Edited by B. J. Kirby. VI, 403 pages. 1974. DM 35,-
Vol. 81: Advanced Course on Software Engineering. Edited by
Vol. 107: Control Theory, Numerical Methods and Computer
F. L. Bauer. XII, 545 pages. 1973. DM 32,- Systems Modelling. International Symposium, Rocquencourt,
Vol. 82: R. Saeks, Resolution Space, Operators and Systems. X, June 17-21, 1974. Edited by A. Bensoussan and J. L Lions. VIII,
267 pages. 1973. DM 22,- 757 pages. 1975. DM 53,-
Vol. 83: NTG/GI-Gesellschaft fur Informatik, Nachrichtentech- Vol. 108: F. Bauer et al., Supercritical Wing Sections II. A Hand-
nische Gesellschaft. Fachtagung ..Cognitive Verfahren und Sy- book. V, 296 pages. 1975. DM 28,-
steme", Hamburg, 11.-13. April 1973. Herausgegeben im Auftrag
der NTG/GI von Th. Einsele, W. Giloi und H.-H. Nagel. VIII, 373
Seiten. 1973. DM 28,-
Vol. 84: A. V. Balakrishnan, Stochastic Differential System I.
Filtering and Control. A Function Space Approach. V, 252 pages.
1973. DM22,-
Vol. 85: T. Page, Economics of Involuntary Transfers: A Unified
Approach to Pollution and Congestion Externalities. XI, 159 pages.
1973. DM 18,-
Okonometrie und Unternehmensforschung
Econometrics and Operations Research
Vol. I Nichtlineare Programmierung. Von H. P. Kiinzi und W. Krelle unter
Mitwirkung von W. Oettli. - Mit 18 Abbildungen. XV, 221 Seiten.
1962. Geb. DM38,-
Vol. II . Lineare Programmierung und Erweiterungen. Von G. B. Dantzig. Ins
Deutsche ubertragen und bearbeitet von A. Jaeger. - Mit 103 Ab-
bildungen. XVI, 712 Seiten. 1966. Geb. DM 68,-
Vol. Ill Stochastic Processes. By M. Girault. - With 35 figures. XII, 126
pages. 1966. Cloth DM 28,-
Vol. IV Methoden der Unternehmensforschung im Versicherungswesen. Von
K.-H. Wolff. - Mit 14-Diagrammed VIII, 266 Seiten. 1966. Geb.
DM 49,-
Vol. V The Theory of Max-Min and its Application 4to Weapons Allocation
Problems. By John M. Danskin. - With 6 figures. X, 126 pages. 1967.
Cloth DM 32,- • •' . ' .
Vol. VI Entscheidungskriterien bei Risikp. Von H. Schneeweiss. - Mit 35
Abbildungen. XII, 214 Seiten. 1967. Geb. DM48,- . . .
Vol. VII Boolean Methods in Operations Research and Related Areas. By P.
L. Hammer (Ivanescu) and S. Rudeanu. With a preface by R. Bellman. -
With,25 figures: XVI, 329 pages. 1968. Cloth DM46,1 -. „;
Vol. VIII Strategy for R & D: Studies in the Microeconomics of Development.
By Th. Marschak, Th. K. Glennan JR., and R. Summers.'- With 44
figures..XIV, 330 pages. 1967. Cloth DM 56,80
, Vol. IX Dynamic Programming of Economic Decisions. By M. J. Beckmann. -
With 9 figures XII, 143 pages. 1968. Cloth DM 28,-
Vol. X Input-Output-Analyse. Von J. Schumann. - Mit 12 Abbildungen. X,
311 Seiten. 1968. Geb. DM 58,-
Vol. XI Produktionstheorie. Von W. Wittmann. - Mit 54 Abbildungen. VIII,
177 Seiten. 1968. Geb. DM 42,- .
Vol. XII Sensivitatsanalysen und parametrische Programmierung. Von W. Din-
kelbach. - Mit 20 Abbildungen. XI, 190 Seiten. 1969. Geb. DM 48,-
Vol. XIII Graphentheoretische Methoden und ihre Anwendungen.' Von W.
•• Knodel. - Mit 24 Abbildungen. VIII, 111 Seiten. 1969. Geb. DM 38,-
Vol. XIV Praktische Studien zur Unternehmensforschung. Von E. Nievergelt,
O. Muller, F. E. Schlaepfer und W^RrLandis. - Mit 82 Abbildungen.
XII, 240 Seiten. Geb. DM 58,-. "~~~—-
Vol. XV Optimale Reihenfolgen. Von [Link]-Merbach.-Mit'43 Abbildungen.
IX, 225 Seiten. 1970. Geb. DM 60,-
Vol. XVI Preispolitik der Mehrproduktenunternehmung in der statischen Theo-
. He. Von R-. Seiten. - Mit 20 Abbildungen. VIII, 195 Seiten. 1970. Geb.
DM 64,- , .-
Vol. XVII Information Theory for Systems Engineers. By L. P. Hyvarinen. - With-
42 figures. VIII, 197 pages. 1970. Cloth DM 44,-
Vol. XVIII Unternehmensforschung im Bergbau. Von F. L. Wilk'e. - Mit 29 Ab-
. bildungen: VIII, 150 Seiten. 1972. Geb. DM 54,-
This series aims to report new developments in mathematical eco-
nomics, econometrics, operations research, and mathematical systems,
research and teaching - quickly, informally and at a high level. The
type of material considered for publication includes:
1. Preliminary drafts of original papers and monographs
2. Lectures on a new field, or presenting a new angle on a classical field
3. Seminar work-outs
4. Reports of meetings, provided they are
a) of exceptional interest and
b) devoted to a single topic.
Texts which are out of print but still in demand may also be considered
if they fall within these categories.
The timeliness of a manuscript is more important than its form, which
may be unfinished or tentative. Thus, in some instances, proofs may be
merely outlined and results presented which have been or will later be
published elsewhere. If possible, a subject index should be included.
Publication of Lecture Notes is intended as a service to the international
scientific community, in that a commercial publisher, Springer-Ver-
lag, can offer a wider distribution to documents which would otherwise
have a restricted readership. Once published and copyrighted, they can
be documented in the scientific literature.
Manuscripts
Manuscripts should comprise not less than 100 pages.
They are reproduced by a photographic process and therefore must be typed with extreme care. Symbols
not on the typewriter should be inserted by hand in indelible black ink. Corrections to the typescript
should be made by pasting the amended text over the old one, or by obliterating errors with white cor-
recting fluid. Authors receive 75 free copies and are free to use the material in other publications. The
typescript is reduced slightly in size during reproduction; best results will not be obtained unless the
text on any one page is kept within the overall limit of 18x26.5 cm (7x10Vz inches). The publishers
will be pleased to supply on request special stationery with the typing area outlined.
Manuscripts in English, German or French should be sent directly to Springer-Verlag New York or
Springer-Verlag Heidelberg
ISBN 3-540-07029-X
ISBN 0-387-07029-X