A noniterative method for the numerical generation of orthogonal curvilinear coordinates for plane annular regions between two arbitrary smooth closed curves was developed. The basic generating equation is the Gaussian equation for an Euclidean space which is solved analytically. The method is applied in many cases and these test results demonstrate that the proposed method can be readily applied to a wide variety of problems. The method can also be used for simply connected regions only by obtaining the solution of the linear equation under the changed boundary conditions

Thompson, J. F.

Warsi, Z. U. A.

English

NASA Technical Reports Server

. • ~~N 81 -
Numerical Generation of Two-Dimensional Orthogonal 
Curvilinear Coordinates in an Euclidean Space* 
z. U. A. Warsi and J. F. Thompson 
Department of Aerospace Engineering 
Mississippi State University 
Mississippi State, MS 39762 
Abstract 
14726 
In this paper a non-iterative method for the numerical generation 
of orthogonal curvilinear coordinates for plane annular regions between 
two arbitrary smooth closed curves has been developed. The basic 
generating equation is the Gaussian equation for an Euclidean space 
which has been solved analytically. The method has been applied in 
many cases and these test results demonstrate that the proposed method 
can be readily applied to a wide variety of problems. The method can 
also be used for simply connected regions only by obtaining the solution 
of the linear equation (19) under the changed boundary conditions. 
Details on the work reported in this paper are available in Reference [1]. 
*This work has been supported in part by the Air Force Office of 
Scientific Research, under Grant AFOSR No. 76-2922 and AFOSR No. 80-0185. 
519 
https://ntrs.nasa.gov/search.jsp?R=19810006212 2020-03-21T14:25:41+00:00Z
520 
Fundamental Ideas of the Method 
All methods of numerical coordinate generation in a two-dimensional 
plane and classified under the method of "elliptic equations" (Refs. 
(2]-(10]), have depended invariably on the solution of Poisson equations 
for the curvilinear coordinates ~(x,y) and n(x,y): 
n2c ~ - l (g rl 2g rl + g rl) 
v S g 11 22 - 12 12 22 11 
g22 
= - - p (~ 11) g , CIa) 
gll 
- -g- Q(t;, 11) (lb) 
~here P(t;,I1), Q(~,I1) are arbitrarily specified control functions, the gij 
i 
are the fundamental metric coefficients, the r jk are the Christoffel 
symbols of the second kind 
(2) 
and 
Implicitly equation (1) implies t~o things: (i) that the coordinates 
for the same domain can also be obtained by solving the Laplace equations 
(3) 
~-
i 
and (ii) since the r jk have first partial derivatives of gij in them, 
equation (1) can also be interpreted as providing a set of constraints 
or relations among the g .. and their first partial derivatives. 
1.J 
In this paper we present another method based on elliptic 
equations and state the problem as follows. 
The three functions gIl' g12' g22 of the curvilinear coordinates 
~,n define an element of length ds in a plane if the Gaussian 
equation with zero curvature 
a Ig ri2 
a~ ( g ) 
11 
o (4) 
holds for every point in the plane, and then the Cartesian coordinates 
are given as 
x = x(~,n), y = y(~,n) 
Equation (4) is identically satisfied by a function a(~,n) defined as 
Specifically, a is the angle of inclination with respect to the x-axis 
of the tangent to the coordinate line n = canst. directed in the sense 
of increasing values of the parameter~. The first partial derivatives 
of x and yare 
x~ I gil cos a, y~ -/gn sin a 
= - I (;g cos a - gl2 sin 
Ig
u 
521 
522 
Then 
The 
x = J[/g11 cos a d~ + _1_ (g12 cos a + ;g sin o.)dn] I g11 
Y =-j[/g11 sin a d~ _1_ (Ii cos a - g12 sin o.)dn] I
g11 
. 
a = -f Ii (r 2 d~ 
g11 11 + ri2 
dn) 
inverse relations of (5) are 
~x = (Ii cos a - g12 sin a) / Iggu 
~y = -(g12 cos Ci. + ig sin a) /lgg11 
For the case of orthogonal coordinates, the coefficient 
g12 = 0, i.e., 
o 
~hich is satisfied by the equations 
~here F > 0 is a continuous function of gIl and g22 [ll] . 
Referring to Figure 1, let the boundary r 2 of a bounded region 
in an Euclidean t~o-dimensional space be a simple curve x = xoo(~)' 
y = y (~), ~ith a uniformaly turning tangent. In the region Q, let 
00 
(6) 
(7) 
(8) 
(9) 
Lx 
I)-I).., 
.., 
r 
If 
n=T) 
o 
I) 
r-- r 
r 
2 
r* 
2 
* r 
1 
~ -
-
Transformed Plane 
(Natural Coordinates) 
I 
J 
I 
J 
I 
I 
: 
i 
Figure 1.- Physical and transformed planes. 
* ra 
r 
3 
r 
It 
523 
Q
s 
be an annular subregion bounded by the inner boundary r l and the 
outer boundary r 2 . The region Q is to be mapped onto a rectangular s 
region R in the sn-plane by a transformation so as to have 
x = x(s,n) } 
where IlS and Iloo are the actual parametric values associated with the 
boundaries fl and f2' respectively, and x,y are periodic in the ~-argument. 
Substituting g12 = 0 in the fundamental equation (4) we get 
a [_l_~ (p2gll ) ] +~ 
1 agu 0 a~ Fgu aE; all 
[pg-- --an-] 
11 
(10) 
where 
g22 = F2g } 11 g == (FgU ) 2 (11) 
Before we solve the problem of orthogonal coordinate generation 
based on the elliptic equation (10), we digress and state the following 
results: Following Potter and Tuttle [ 6] we assume that the E;-curves 
in the physical xy-plane are free from sources and sinks. This 
condition establishes a unique correspondence between the E;-points on 
each pair of Il = const. lines. In the absence of sources and sinks, 
we have 
div[grad ~(Il)] = 0 (12) 
where wen) is an arbitrary differentiable function of Il and grad w(n) 
is oriented along the normal to the curve Il == canst. Carrying out the 
differential operation in (12) and using the expressions 
Igrad III = 1/lg22 
524 
in (12), we obtain 
a d2 ljJ d~ an Un Ig ll / g22 )= - -/-d 
'I dn 2 11 
Writing ~ = liven) and denoting the arbitrary function of ~ due to dn 
integration by ~n~(~), we obtain the result 
= l/F 
This result shows that for the case of orthogonal coordinates the 
ratio gll/g22 is a product of the positive functions ~(~) and v(n). 
The result in (13) also provides the condition for the two distinct 
families of orthogonal curves 
~ = const., 11 = const. 
(13) 
to divide the physical plane in infinitensimal squares. (See Cohen [12)). 
We now introduce new coordinates ~'(~) and n'(n) as 
~' = J~(~)d~ n' = J ~ 
, v (n) (14) 
Thus 
so that 
Defining 
525 
526 
it can be shown that 
2 ag 
a 2 p , a P' \I ( n) {a 1 a 2 )] +..l.. [_l_~]} 
-- + -- = -(t") 3F [Fg
ll 
~ (F gl1 an Fg
ll 
an a~,2 an'2 )J <, <, 
(15) 
Using (15) in (10), we get a much simpler equation 
Another important result can be obtained based on (13). Using 
the orthogonality condition g12 = 0 in (7), we have 
n = ~ IF n = - ~ IF Y x' x y 
so that 
a (~/F) + -!- (~ IF) = 0 a~ x oy y (17) 
Carrying out the transformation (14) in (17), we get 
(18) 
Equation (18) provides the uniqueness condition for the solution of 
equation (16). 
Based on the preceeding analysis we can state that if an exact 
analytic solution of equation (10) can be obtained for F = 1, i.e., 
g22 = gIl' then the solution for any other coordinate system ~ and n, 
where ~ = ~ (~) and n = fen), can simply be obtained by the substitution 
of ~ and f in place ~ and n respectively. With this scheme in mind, 
we solve the equation 
(19) 
where 
P = 2n gll (20) 
under the boundary conditions 
= P (c) at n = n 
00 s 00 
I o ..::. f,; < 21T (21) 
where the subscripts Sand 00 denote the inner and outer boundaries, 
respectively. The periodicity requirement is that 
P(f,;,n) = P(f,;+21T,n) (22) 
Further, the f,;-coordinate must be such that the equation 
(23) 
is always satisfied. 
A general analytic solution of equation (19) under the boundary 
conditions (21) and the periodicity condition (22) is 
00 
+ L 1 sinh nn (c cos nf,; + d sin nO / sinh nno:> (24) 
n= n n 
where 
K = (c -a )/n 
00 00 
(25) 
and 
tThere is no loss of generality in setting the parametric value nS = o. 
The value n must be interpreted as the difference between the actual 
00 
~' values of n at the outer and inner boundaries. 
527 
27T I fa Ps(O cos ns ds, a n 7T 
I 27T 
c fa P (s) cos ns ds, n 7T 00 
For orthogonal Coordinates 
therefore for g22 = gIl 
and consequently 
1 ap 
a =-- a 
s 2 an' n 
I 27T PS(O b = - fa n 7T 
I 27T 
d =-J P (0 n 7T a 00 
1 ag 22 
:;:: - ------ZIg a~ 
cosh n(n -n) 
00 00 
sin ns ds 
sin ns ds 
a(Cn) = a(E;,o) + ~=l 2 sinh n n (b cos nE; - a sin nO n n 
00 
00 cosh nn 
+ ~=l 2 sinh n n (c sin nE; - d cos ns) n n 
00 
00 
- [ 
n=l 2 sinh n v 
cosh n noo 
(b cos nE; - a sin ns) 
n n 
00 
00 1 [ (c sin ns - d cos nE;) 
n=l 2 sinh n noo n n 
} 
Having determined gIl and a, ~e can find the Cartesian coordinates 
x(s,n) = x(s,o) + J; I g22 sin a dn 
y(~,n) = y(s,o) + f: I g22 cos a dn } 
The preceding solution is for the case ~hen g22 = gIl' i.e., 
(26) 
(28) 
F = 1. Ho~ever, as stated earlier, the solution for any other coordinate 
system E;,n in ~hich g22 ¥ gIl can be obtained by replacing E; and n 
in (24), (27) and (28) by ~(~) and fen). respectively. This feature can ~ 
528 
be used to redistribute the coordinate lines in the desired regions. 
Since the functions ~ and f are at our disposal, they play the same role 
as P and Q in equation (1). Further, since the Fourier coefficients 
(26) are invariant to a coordinate transformation ~ = ~(~), where 
~(~o) = 0, ¢(~m) = 2n and ~o corresponds to ~ = 0, ~m corresponds to 
2n, these coefficients need not be recalculated. 
The procedure of transformation from ~,n to s,n is as follows. 
On transformation from (~,n) to (~,n), the covariant metric 
coefficients transform through the equation 
so that, on using the relations g22 = gIl and g12 0, we have 
[(~) 2 an 2 ) + (--=) ] g11 ae.; ae.; [ (a~) 2 an 2 (29) + e--=) ] g11 
an an 
We now introduce the transformation 
n fen) } (30) 
where the functions ¢ and f are continuously differentiable and satisfy 
the conditions 
where e.; = 0 corresponds to e.; = sand n 
a 
Defining 
o 
o corresponds to n = nS. 
529 
530 
~e obtain from (30) 
1 Comparing ~ith (13), ~e find ~ ~ A, V =-8 (31) 
The salient feature of the preceding analysis is that the solution 
for the case g22 = gIl can be used to obtain the solution for the case 
g22 I gIl by coordinate transformation. 
Before solving any specific problem, it is important first to 
establish an orthogonal correspondence between unique points of the 
inner and outer boundary curves. This condition is satisfied if we 
choose the ~-curves satisfying the equation 
(32) 
The inner and outer boundary curves are available either in tabular 
or functional form as 
y(x ) 
00 
(33) 
For equation (32) to be satisfied, ~e can take s as the angle traced 
out in a clockwise sense by the common radius of the concentric circles 
in a conformal representation of the inner and outer boundary curves. 
If a and A, respectively, are the radii of the inner and outer circles 
in the transformed conformal plane, then 
s]ds I 
s]ds (34) 
1 2TI 
A = -- f [x (s) cos S - y (s) sin 2n 0 00 00 
" . ",- -,-
n~··:' .. 
r .1 ' 
As is ~ell kno~n, the preceding scheme is an iterative numerical scheme. 
In lieu of this, ~e have developed a method ~hich is fast and direct, 
and is equivalent to satisfying equation (32). 
We circumscribe circles around the inner and outer boundary curves. 
T~o cases arise depending on whether the circles are concentric or 
nonconcentric. 
Case I: If the circumscribed circles are concentric (Fig. 2(a), then we 
select those values of the ordinates ~hich correspond to the abscissae 
(35) 
~here rs and r L are the radii of the circumscribed circles in the 
physical plane. 
~ Case II: If the circumscribed circles are nonconcentric (Fig. 2(b», 
then ~e first use the formula for the conformal transformation of non-
concentric to concentric circles, Kober [13] , and choose the ordinates 
corresponding to abscissae given by the formula 
~here 
xes) = [(1 - cy cos s) {xL(l - cy cos s) + cYYL sin ~ 
- r (c sin ~ + y sin(~-~»}J L 
rL,r
s 
= radii of outer and inner circumscribed circles 
(36) 
531 
Lx 
532 
(a) Concentric circumscribed circles Cl 
and C2 of radii rs and rL' respectively, 
with center at the origin. 
(b) Nonconcentric circumscribed circles Cl 
and C2 of radii rs and rL and centers 
at Zs and zL' respectively. 
Figure 2.- Circumscribed circles. 
---
c [(d 2 + r2 L 
y = 1 for the outer boundary 
y ~ rLjd-tj for the inner boundary 
r t 
s 
t = cr L 
Having determined the appropriate sets (xs(~)' Ys(~» and 
(x (~), y (~», we use (34) to obtain the values of a and A. The 
00 00 
parametric difference n is connected in some manner with the "modulus" 
00 
of the domain which, however, by itself is a separate problem (see 
Burbea [14] and Gaier [15]). In this work we have defined n based 
00 
on the knowledge of a and A as discussed above by the formula 
A Q,n(-) 
a 
(37) 
For Figures 3 to 8, we have used the following functional forms 
of <P and f: 
fen) 
21T(~-t: ) 
o 
n -n 
00 6 
t: -~ 
m 0 
K(n-n@) 
K(nro-n S) 
533 
534 
so that 
A = 21T 
[1 + 
where K > 1 is an arbitrary constant, and ~ = ~ , n 
m 
1100 correspond, 
respectively, to ~ = 21T and n = noo' We treat ~ and n as integers so 
that ~ = 1 ~ 
o 'm lMAX, 118 = 1, and noo = JMAX. Since 11 is known 00 
from (37), hence by specifying the numerical values to K and JMAX 
we can create the desired mesh control in the direction of T). The 
value of K between 1.05 and 1.1 is quite sufficient [16] to have a 
fine grid near the inner boundary. 
The number of terms to be retained in the series (24) is usually ~ 
small for convex inner and outer boundaries, though we have retained 
(IMAX-l)/2 number of coefficients in each computation. This number 
is the optimum number of terms in a discrete Fourier series [17] 
having lMAX number of points in one period. 
Figure 3 shows the classic case of confocal ellipses with coordinate 
contraction in n. The value of K is 1.05. The orthogonal correspondence 
between ~-points of the inner and outer boundary has been established 
by using Case I, Eq. (35). 
Figure 3.- Confocal ellipses. Semimajor axes 1.48, 5.0, 
and semiminor axes 0.5, 4.802, respectively. Only 
38 n = Const. lines shown for detail. 
535 
Figure 4 presents orthogonal coordinates for a blunt body with elliptical '--
outer boundary. Here K = 1.01. For orthogonal correspondence between ~-points, 
Eq. (35) has been used. 
Figure 4.- A blunt body section with elliptical outer boundary. 
536 
~ Orthogonal coordinates for nonconcentric circles are presented in 
Figure 5. Here K = 1.01. For orthogonal correspondence between F,-points 
between the inner and outer boundary, Eq. (36) has been used. Data shown on 
the figure. 
Figure 5.- Nonconcentric circles: rs 
Zs = (0,0), zL = (1,0). 
2.5, 
537 
Orthogonal coordinates for a Joukowsky's airfoil with slightly rounded '--
trailing edge are shown in Figure 6. Eq. (35) is used for orthogonal 
correspondence. Here K = 1.02. 
538 
Figure 6.- Joukowsky's airfoil with slightly 
rounded trailing edge. 
Figure 7 presents orthogonal coordinates for nonconcentric ellipses. 
Centers of the inner and outer ellipses are at (0,0) and (1,0), respectively. 
Here K = 1.01. For orthogonal correspondence Eq. (36) has been used. 
Figure 7.- Nonconcentric ellipses. Size data same as in 
Figure 3. Zs = (0,0), zL ~ (1,0). 
539 
Orthogonal coordinates for an arbitrarily deformed upper part of 
Figure 4 are shown in Figure B. The placement of outer boundary is limited 
to avoid intersecting normals (Eiseman [lBJ). This figure shows that we need 
some attraction near those sections of the outer boundary which face the 
concave side. 
540 
Figure B.- Generated coordinates for body having convex, 
concave and straight portions. Placement of outer 
boundary is decided by the radius of the osculating 
circles of the concave portions. 
Summary of Numerical Experimentation 
In the course of this investigation a number of cases of inner 
and outer boundary shapes and orientations have been tested through 
the developed computer program. The main conclusions are listed below: 
(i) The method ~orks very effectively for smooth and convex boundaries 
of any shape and orientation. 
(ii) For concave boundaries a method similar to that of Eiseman has 
to be used in the placement of the outer boundary to avoid intersecting 
normals. Another remedy ~ould be to introduce some type of 
attraction near the outer boundary facing the concave side of the 
inner boundary. 
(iii) Sharp turns and corners are not admissible and have to be rounded 
to avoid singularities in the metric data. 
541 
References 
1. Z. U. A. Warsi, R. A. Weed, and J. F. Thompson, "Numerical Generation 
of Two-Dimensional Orthogonal Coordinate in an Euclidean Space", 
Engineering and Experimental Research Station, Mississippi State 
University, Rep. No. MSSU-EIRS-ASE-80-3, June 1980. 
2. A. J. Winslow, "Numerical Solution of the Quasi-Linear Equation in a 
Non-Uniform Triangular Mesh", Journal of Computational Physics, 2, 
149 (1966). 
3. W. D. Barfield, "An Optimal Mesh Generator for Lagrangian Hydrodynamic 
Calculations in Two Space Dimensions", Journal of Computational Physics, 
6, 417 (1970). 
4. W. H. Chu, "Development of a General Finite Difference Approximation 
for a General Domain, Part I: Machine Transformation", Journal of 
Computational Physics, 8, 392 (1971). 
5. A. A. Amsden and C. W. Hirt, "A Simple Scheme for Generating General 
Curvilinear Grids", Journal of Computational Physics, 11, 348 (1973). 
8. S. B. Pope, "The Calculation of Turbulent Recirculating Flows in 
General Orthogonal Coordinates", Journal of Computational Physics, 
26, 197 (1978). 
9. J. F. }-liddlecoff and P. D. Thomas, "Direct Control of the Grid Point 
Distribution in Meshes Generated by Elliptic Equations", AlAA 
Computational Fluid Dynamics Conference, Paper No. 79-1462 (1979). 
10. R. L. Sorenson and J. L. Steger, "Simplified Clustering of Nonorthogonal 
Grids Generated by Elliptic Partial Differential Equations", NASA-TM-
73252, August 1977. 
11. G. Starius, "Constructing Orthogonal Curvilinear Meshes by Solving 
Initial-Value Problem", Numerische Mathematik, 28, 25 (1977). 
12. A. Cohen, "An Introduction to the Lie Theory of One-Parameter Groups", 
G. E. Stechert and Co., New York (1931), p. 2-5. 
13. H. Kober, "Dictionary of Conformal Representations", Dover Publications, 
Inc. (1952). 
542 
.........,..--
14. J. Burbea, "A Numerical Determination of the Modulus of Doubly 
Connected Domains by Using the Bergman Curvature", Mathematics 
of Computation, 25, 743 (1971). 
15. D. Gaier, "Determination of Conformal Modulus of Ring Domains and 
Quadrilaterals", In-Lecture Notes in Mathematics, No. 399, Springer-
Verlag, Berlin (1974). 
16. Z. U. A. I\farsi and J. F. Thompson, "Machine Solutions of Partial 
Differential Equations in the Numerically Generated Coordinate 
Systems", Engineering and Experimental Research Stateion, Mississippi 
State University, Rep. No. MSSU-EIRS-ASE-77-1, August (1976). 
17. F. B. Hildebrand, "Introduction to Numerical Analysis", McGraw-Hill 
Book Company, Inc., (1956). 
18. P. R. Eiseman, "A Coordinate System for a Viscous Transonic Cascade 
Analysis", Journal of Computational Physics, 26, 307 (1978) . 
543 


Numerical generation of two-dimensional orthogonal curvilinear coordinates in an Euclidean space

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server