Abstract:
A method for designing a flow device, including the steps of defining the physical properties of the flow device, and deriving from the physical properties the corresponding geometry of the flow device using inverse flow equations. The flow device can be at least a part of a gas turbine, a steam turbine, or a combined gas/steam turbine-cycle plant.

Description:
FIELD OF THE INVENTION 
     The present invention relates a method for designing a flow device. 
     PRIOR ART 
     The documents cited in the References section of the present application describe related technologies. 
     1. Introduction and Discussion of Background 
     The present invention is related to the design and construction of flow devices and especially turbine and compressor blade channels. 
     In recent years substantial efforts have been made to introduce inverse design tools for airfoil, turbine blades and vane channel, and other flow devices for which the optimum combination of aerodynamic forces and losses plays a crucial role. For an adequate overview of the state of the art the reader is referred to work by Jameson et al. [3], Giles and Drela [4], Dang and Isgro [5], Dulikravich [6], Pandya and Baysal [7]. An overview of applications to current turbomachinery design is given by Jennions [8]. 
     It should be understood that this list of publications on inverse design methods is far from complete. Furthermore, it should be mentioned that formulations of the Euler equations using two stream functions have been previously used (see Turner and Giles [9], for example). However, with the exception of genuine inverse methods (Keller et al. [10], [11]) that have been proposed for certain two-dimensional flows, the “inverse methods” that are presently used to arrive at an optimum design of three-dimensional flow devices are not really inverse. In general “Evolution Strategy”, genetic algorithms or simply Newton&#39;s method are used to guide an extensive series of flow calculations for gradually varying geometry toward some kind of optimum. 
     SUMMARY OF THE INVENTION 
     According to certain aspects of the present invention, a method for designing a flow device comprises the steps of (1) defining the physical properties of said flow device, and (2) deriving from said physical properties the corresponding geometry of said flow device by means of inverse flow equations. 
     Still other objects, features, and attendant advantages of the present invention will become apparent to those skilled in the art from a reading of the following detailed description of embodiments constructed in accordance therewith, taken in conjunction with the accompanying drawings. 
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS 
     The invention of the present application will now be described in more detail with reference to preferred embodiments thereof, given only by way of example, and with reference to the accompanying drawings, in which: 
     FIG. 1 illustrates distributions of flow speed along the two walls of an annular diffuser. The solid lines mark curves of maximum adverse pressure gradient. Smooth transitions between these curves have been generated using fourth-order polynomials, in order to define suitable distributions of flow speed along the diffuser walls. 
     FIG. 2 illustrates a map of streamlines and curves of constant σ that corresponds to the boundary conditions defined by the distributions of flow speed shown in FIG.  1 . 
     FIG. 3 illustrates a map of streamlines and curves of constant σ of a turbine inlet flow. 
    
    
     DETAILED DESCRIPTION OF EXEMPLARY EMBODIMENTS 
     Referring to the drawing figures, like reference numerals designate identical or corresponding elements throughout the several figures. 
     The key advantage of using inverse equations instead, is associated with the fact that the geometry of a flow device does not have to be defined in advance. As an alternative it is possible to impose the distribution of static pressure, flow speed, temperature, heat flux or any other physical quantity, as a boundary condition and ask for the geometry (as a part of the solution) that generates the desired distribution of static pressure, etc. Furthermore, if the original domain of definition is simply connected, a one-to-one mapping into a cubic domain of definition with Cartesian coordinates can be achieved. As a consequence, the generally tedious task of grid generation for numerical calculations can be avoided. 
     2. Stream Functions and Natural Coordinate 
     To define boundary conditions and to embed boundary layers in an inviscid bulk flow, for example, it may be useful to introduce the “natural coordinate” σ that is related to the arc lengths along stream lines. Following the ideas of Ref. [2] we introduce a pair of stream functions to define a vector potential for the velocity field. Without essential loss of generality we restrict our considerations to solenoidal vector fields for the mass flux. In order to extend the discussion presented in Ref. [2] to compressible flows the mass flux (rather than the velocity) needs to be related to the stream functions. Hence, introducing the two stream functions ψ and χ, we can write 
     
       
         ρ u =grad(ψ)×grad(χ),  (1)  
       
     
     where  u  is the velocity vector and ρ the density. The natural coordinate can be defined by                h                   grad        (   σ   )         =           grad        (   ψ   )       ×     grad        (   χ   )                  grad        (   ψ   )       ×     grad        (   χ   )                .             (   2   )                                
     where σ is the dimensionless natural coordinate and h=h(σ,ψ,χ) is a nonlinear stretching coefficient, such that hdσ denotes the actual length increment along the streamlines. As usual the nonlinear stretching needs to be included, in order to reach compatibility between the orthogonality requirement (2) and the uniform definition of hdσ as the length increment along the streamlines. 
     Furthermore, in inverse coordinates we can write                      ∂     x   _         ∂   σ       ·       ∂     x   _         ∂   σ         =     h   2       ,       x   _     =       x   _          (     σ   ,   ψ   ,   χ     )         ,           (   3   )                                
     wherein  x  is the coordinate vector in the original space. 
     First we need to express the time-independent Euler equations in terms of the two stream functions ψ and χ. 
     3. The Governing Equations in a Rotating Frame of Reference 
     An area of major interest to use inverse equations is the design of optimized turbine and compressor blade channels. For this reason it seems to be appropriate to further extend the governing equations to rotating frames of reference. In a frame of reference rotating steadily at the angular velocity  Ω  the Euler equation can be written in the form                      (       u   _     ·   grad     )          u   _       +       1   ρ        grad                 p       =         -   2          Ω   _     ×     u   _       -       Ω   _     ×     (       Ω   _     ×     x   _       )           ,           (   4   )                                
     wherein P refers to the static pressure. The terms on the right-hand side of (4) contribute to the fictitious volume force. The first term accounts for the Coriolis force and the second for the centrifugal force. Making use of the identities                  (       u   _     ·   grad     )          u   _       =       grad        (       1   2            u   _     ·     u   _         )       -       u   _     ×     ω   _                 (   5   )             and                               -     Ω   _       ×     (       Ω   _     ×     x   _       )       =     grad        (         1   2          [       Ω   _     ×     x   _       ]       ·     [       Ω   _     ×     x   _       ]       )         ,           (   6   )                                
     where  ω  refers to the vorticity,                      ω   _     =                  rot        (     u   _     )       =     rot        (       1   ρ          grad        (   ψ   )       ×     grad        (   χ   )         )                     =                    1   ρ          rot        (       grad        (   ψ   )       ×     grad        (   χ   )         )         +       grad        (     1   ρ     )       ×     [       grad        (   ψ   )       ×     grad        (   χ   )         ]                       (   7   )                                
     we can rewrite (4) as                  grad        (         1   2            u   _     ·     u   _         -         1   2          [       Ω   _     ×     x   _       ]       ·     [       Ω   _     ×     ]         )       -       u   _     ×     [       ω   _     +     2        Ω   _         ]       +       1   ρ        grad                 p       =   0.           (   8   )                                
     The first law of thermodynamics can be expressed as                Tds   =         c   v        dT     -       p   ρ            d                 ρ     ρ           ,           (   9   )                                
     where T is the temperature, s denotes the entropy per unit mass and c v  refers to the specific heat per unit mass at constant volume. The equation of state is                  p   ρ     =       (       c   p     -     c   v       )        T       ,           (   10   )                                
     and c p  refers to the specific heat per unit mass at constant pressure. Combining (8), (9) and (10) yields                  grad        (         c   p        T     +       1   2            u   _     ·     u   _         -         1   2          [       Ω   _     ×     x   _       ]       ·     [       Ω   _     ×     x   _       ]     ·     [       Ω   _     ×     x   _       ]         )       -       u   _     ×     [       ω   _     +     2        Ω   _         ]       -     T                 grad                 s       =   0.           (   11   )                                
     Introducing the total enthalpy              H   =         c   p        T     +       1   2            u   _     ·     u   _         -         1   2          [       Ω   _     ×     x   _       ]       ·     [       Ω   _     ×     x   _       ]                 (   12   )                                
     equation (11) can now be expressed in the form of Crocco&#39;s theorem in the rotating frame: 
     
       
         gradH− u ×[ ω +2 Ω ]−Tgrad s=0.  (13)  
       
     
     By virtue of (1) the condition of mass conservation is identically satisfied. Moreover, noting that the validity of the Euler equations is restricted to adiabatic flows, it is obvious that the s is constant on streamlines. Hence, the entropy is an integral of motion (in any frame of reference), 
     
       
         s=s(ψ,χ) u ·grads=0.  (14)  
       
     
     Furthermore, it is obvious that 
     
       
         ( u ×[ ω +2 Ω ])· u =0.  (15)  
       
     
     By virtue of (14) and (15) the scalar product of the velocity vector with equation (13) yields 
     
       
           u ·gradH=0H=H(ψ,χ)  (16)  
       
     
     This is Bernoulli&#39;s theorem in the rotating frame of reference. It shows that the total enthalpy in the rotating frame is an integral of motion. Making use of (1), (14) and (16) equation (13) can be rewritten as                    (         ∂   H       ∂   ψ       -     T          ∂   s       ∂   ψ           )        grad                 ψ     +       (         ∂   H       ∂   χ       -     T          ∂   s       ∂   χ           )        grad                 χ     -         1   ρ          [     grad                 ψ   ×   grad                 χ     ]       ×     [       ω   _     ×   2        Ω   _       ]         =   0           (   17   )                                
     Furthermore, using the identity 
     
       
         −[gradψ×gradχ]×[ ω +2 Ω ]=([ ω +2 Ω ]·gradχ)gradψ−([ ω +2 Ω ]·gradψ)gradχ  (18)  
       
     
     the vector products of (17) with gradψ and gradχ, respectively, lead to the following generalized Euler equations in terms of two stream functions                    ∂   H       ∂   ψ       -     T          ∂   s       ∂   ψ         +           1   ρ          [       ω   _     +     2        Ω   _         ]       ·   grad                   χ       =   0           (   19   )             and                               ∂   H       ∂   χ       -     T          ∂   s       ∂   χ         -           1   ρ          [       ω   _     +     2        Ω   _         ]       ·   grad                   ψ       =   0.           (   20   )                                
     For the subsequent discussion it is useful to introduce the mixed product that is defined by 
     
       
         { abc }=[a×b]· c   (21)  
       
     
     An aim of the present paper is to derive the inverse forms of Euler&#39;s equations, using the two stream functions ψ and χ and the natural (streamline) coordinate σ as independent variables. 
     4. Mapping to Inverse Coordinates 
     In this section we consider the mapping 
     
       
           x ←→(σ,ψ,χ).  (22)  
       
     
     The derivatives in the two coordinate system are related by                    grad        (   ψ   )       ×     grad        (   χ   )         =       1   J            ∂     x   _         ∂   σ           ,           (   23   )                     grad        (   χ   )       ×     grad        (   σ   )         =       1   J            ∂     x   _         ∂   ψ           ,                                 grad        (   σ   )       ×     grad        (   ψ   )         =       1   J            ∂     x   _         ∂   χ           ,                                            
     where the mixed product              J   =     {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ              ∂     x   _         ∂   χ         }             (   24   )                                
     may be interpreted as Jacobi&#39;s determinant. Moreover it can be shown that                  grad        (   σ   )       =         1   J          [       ∂     x   _         ∂   ψ       ]       ×     [       ∂     x   _         ∂   χ       ]         ,           (   25   )                   grad        (   ψ   )       =         1   J          [       ∂     x   _         ∂   χ       ]       ×     [       ∂     x   _         ∂   σ       ]         ,                             grad        (   χ   )       =         1   J          [       ∂     x   _         ∂   σ       ]       ×       [       ∂     x   _         ∂   ψ       ]     .                                                
     For a function ƒ( x )=g(σ,ψ,χ) the conversion of derivatives can be expressed as                  grad        (   f   )       =       1   J          {         (       [       ∂     x   _         ∂   ψ       ]     ×     [       ∂     x   _         ∂   χ       ]       )            ∂   g       ∂   σ         +       (       [       ∂     x   _         ∂   χ       ]     ×     [       ∂     x   _         ∂   σ       ]       )            ∂   g       ∂   ψ         +       (       [       ∂     x   _         ∂   σ       ]     ×     [       ∂     x   _         ∂   ψ       ]       )            ∂   g       ∂   χ           }         ,           (   26   )                                
     and for a vector  a ( x )= b (σ,ψ,χ) we can write                div        (     a   _     )       =       1   J          (       (       [       ∂     x   _         ∂   ψ       ]     ×     [       ∂     x   _         ∂   χ       ]       )     -       ∂     b   _         ∂   σ       +       (       [       ∂     x   _         ∂   χ       ]     ×     [       ∂     x   _         ∂   σ       ]       )     ·       ∂     b   _         ∂   ψ         +       (       [       ∂     x   _         ∂   σ       ]     ×     [       ∂     x   _         ∂   ψ       ]       )     ·       ∂     b   _         ∂   χ           )               (   27   )             and                           rot        (     a   _     )       =       1   J            (         (       [       ∂     x   _         ∂   ψ       ]     ×     [       ∂     x   _         ∂   χ       ]       )     ×       ∂     b   _         ∂   σ         +       (       [       ∂     x   _         ∂   χ       ]     ×     [       ∂     x   _         ∂   σ       ]       )     ×       ∂     b   _         ∂   ψ         +       (       [       ∂     x   _         ∂   σ       ]     ×     [       ∂     x   _         ∂   ψ       ]       )     ·       ∂     b   _         ∂   χ           )     .               (   28   )                                
     Using (25), (26) and (28) it may be shown that                  grad        (   χ   )       ·     rot        (       grad        (   ψ   )       ×     grad        (   χ   )         )         =              [       1   J            ∂     x   _         ∂   ψ         ]     ·     (       ∂     ∂   σ            [       1   J            ∂     x   _         ∂   σ         ]       )       -       [       1   J            ∂     x   _         ∂   σ         ]     ·     (       ∂     ∂   ψ            [       1   J            ∂     x   _         ∂   σ         ]       )                   (   29   )                                
     Similarly (or by analogy) we obtain                  grad        (   ψ   )       ·     rot        (       grad        (   ψ   )       ×     grad        (   χ   )         )         =              [       1   J            ∂     x   _         ∂   ψ         ]     ·     (       ∂     ∂   σ            [       1   J            ∂     x   _         ∂   σ         ]       )       +       [       1   J            ∂     x   _         ∂   σ         ]     ·     (       ∂     ∂   χ            [       1   J            ∂     x   _         ∂   σ         ]       )                   (   30   )                                
     Furthermore, using (25) it may be shown that                  Ω   _     ·     grad        (   χ   )         =         1   J          {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ            Ω   _       }       =       {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ            Ω   _       }       {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ              ∂     x   _         ∂   χ         }                 (   31   )             and                             Ω   _     ·     grad        (   ψ   )         =         1   J          {         ∂     x   _         ∂   σ            Ω   _            ∂     x   _         ∂   χ         }       =       {         ∂     x   _         ∂   σ            Ω   _            ∂     x   _         ∂   χ         }       {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ              ∂     x   _         ∂   χ         }                 (   32   )                                
     Making use of (7), the identity 
     
       
           a ×( b × c )=( a · c ) b −( a · b ) c   (33)  
       
     
     and (29) to (32) the equations of motion, (19) and (20), can be written as                    ∂   H       ∂   ψ       -       ∂     ∂   ψ            (         1   2          [       1     ρ                 J              ∂     x   _         ∂   σ         ]       ·     [       1     ρ                 J              ∂     x   _         ∂   σ         ]       )       -     T          ∂   s       ∂   ψ         +       [       1     ρ                 J              ∂     x   _         ∂   ψ         ]     ·     (       ∂     ∂   σ            [       1     ρ                 J              ∂     x   _         ∂   σ         ]       )       +       2     ρ                 J            {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ            Ω   _       }         =   0           (   34   )             and                               ∂   H       ∂   χ       -       ∂     ∂   χ            (         1   2          [       1     ρ                 J              ∂     x   _         ∂   σ         ]       ·     [       1     ρ                 J              ∂     x   _         ∂   σ         ]       )       -     T          ∂   s       ∂   χ         +       [       1     ρ                 J              ∂     x   _         ∂   χ         ]     ·     (       ∂     ∂   σ            [       1     ρ                 J              ∂     x   _         ∂   σ         ]       )       -       2     ρ                 J            {         ∂     x   _         ∂   σ            Ω   _            ∂     x   _         ∂   χ         }         =   0           (   35   )                                
     Note that the mapping described in this section is thus far valid for any choice of the variable σ. No use has been made yet of the fact that σ should ultimately be interpreted as the natural coordinate as defined by (2) and (3). 
     The subsequent part of the analysis shall now account for the special significance of σ as a natural coordinate. From (2), (3), (23) and (25) we obtain                  ∂     x   _         ∂   σ       =           h   2     J          [       ∂     x   _         ∂   ψ       ]       ×       [       ∂     x   _         ∂   χ       ]     .               (   36   )                                
     Hence, (36) implies that the covariant and contravariant forms of ∂ x /∂σ differ by the factor h 2 . Furthermore, (36) implies                      ∂     x   _         ∂   σ       ·       ∂     x   _         ∂   σ         =     h   2       ,           ∂     x   _         ∂   σ       ·       ∂     x   _         ∂   ψ         =   0     ,           ∂     x   _         ∂   σ       ·       ∂     x   _         ∂   χ         =   0.             (   37   )                                
     The first relation (37) is the same as (3). Using (37) the equations (34) and (35) can be simplified to                    ∂   H       ∂   ψ       -       1     ρ                 J            ∂     ∂   ψ            (       h   2       ρ                 J       )       -     T          ∂   s       ∂   ψ         +       2     ρ                 J            {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ            Ω   _       }         =   0           (   38   )             and                               ∂   H       ∂   χ       -       1     ρ                 J            ∂     ∂   χ            (       h   2       ρ                 J       )       -     T          ∂   s       ∂   χ         -       2     ρ                 J            {         ∂     x   _         ∂   σ            Ω   _            ∂     x   _         ∂   χ         }         =   0           (   39   )                                
     It is quite remarkable that for irrotational isentropic flow, 
     
       
         dH=0, ds=0,  Ω = 0 .  (40)  
       
     
     equations (38) and (39) can be integrated once. 
     5. Thermodynamic Relations 
     Making use of (12) and the definitions 
     
       
         H≡c p T 0 , u≡{square root over ( u )}· u ,  W ≡ Ω ×χ, W≡{square root over ( W )}· W , γ≡c p /c v , c 2 ≡(c p −c v )γT, M≡u/c,  (41)  
       
     
     we can write                  c   p          T   o       =         c   p        T     +       1   2          u   2       -       1   2          W   2                 (   42   )             or                           T     T   o       =         1   +       W   2       2        c   p          T   o             1   +         γ   -   1     2          M   1           .             (   43   )                                
     T 0  is the stagnation temperature in the absolute system, c is the sound speed and M the Mach number in the rotating frame. Hence, 
     
       
         T 0 =T 0 (ψ,χ).  (44)  
       
     
     Combining (14) with the first law of thermodynamics and the equation of state we can shown that the isentropic relations                  ρ     ρ   o       =       [     T     T   o       ]       1     γ   -   1           ,       p     p   o       =       [     T     T   o       ]       γ     γ   -   1           ,       c     c   o       =       [     T     T   o       ]       1   2         ,           (   45   )                                
     hold along streamlines. In these expressions ρ 0 , P 0  and c 0  should be interpreted as the stagnation density, the stagnation pressure and the stagnation sound speed, respectively, in the absolute system and, by virtue of (14), 
     
       
         ρ 0 =ρ 0 (ψ,χ), P 0 =P 0 (ψ,χ), c 0 =C 0 (ψ,χ).  (46)  
       
     
     It is important to note that with the help of (14) and 44) to 46) all thermodynamic quantities can be related to integrals of motion! For the mass flux we obtain (see (3) and the first equation (23))                        m   2     ≡                  ρ   2          u   2         =       [       grad        (   ψ   )       ×     grad        (   χ   )         ]     ·     [       grad        (   ψ   )       ×     grad        (   χ   )         ]                   =                    h   2       J   2       =       ρ   o   2              c   o   2          [     1   +       W   2       2        c   p          T   o           ]           γ   +   1       γ   -   1                  M   2         [     1   +         γ   -   1     2          M   2         ]         γ   +   1       γ   -   1           .                       (   47   )                                
     We may now note that the equations of motion (38), (39), together with the condition (36), that defines the natural coordinate σ, and the thermodynamic relations constitute a complete set of equations in inverse coordinates.. 
     6. Cylindrical Polar Coordinates as Dependent Variables 
     For many applications, including flows through turbines, it is convenient to express the coordinate vector  x  in terms of cylindrical polar coordinates. Introducing cylindrical polar coordinates (x,r,θ) it may be shown that Jacobi&#39;s determinant assumes the form                    J   =                             ∂   x       ∂   σ               ∂   x       ∂   ψ               ∂   x       ∂   χ                   ∂   r       ∂   σ               ∂   r       ∂   ψ               ∂   r       ∂   χ                 r          ∂   θ       ∂   σ               r          ∂   θ       ∂   ψ               r          ∂   θ       ∂   χ                                =                      ∂   x       ∂   σ              ∂   r       ∂   ψ          r          ∂   θ       ∂   χ         +         ∂   x       ∂   ψ              ∂   r       ∂   χ          r          ∂   θ       ∂   σ         +         ∂   x       ∂   χ              ∂   r       ∂   σ          r          ∂   θ       ∂   ψ         -                                    ∂   x       ∂   σ              ∂   r       ∂   χ          r          ∂   θ       ∂   ψ         -         ∂   x       ∂   ψ              ∂   r       ∂   σ          r          ∂   θ       ∂   χ         -         ∂   x       ∂   χ              ∂   r       ∂   ψ          r          ∂   θ       ∂   σ                         (   48   )                                
     Assuming that the angular velocity vector  Ω  is parallel to the axis and points in its positive direction, we find                {         ∂     x   _         ∂   σ              ∂     x   _         ∂   ψ            Ω   _       }     =       Ω                 r          ∂   r       ∂   σ              ∂   θ       ∂   ψ         -     Ω                 r          ∂   r       ∂   ψ              ∂   θ       ∂   σ                   (   49   )             and                           {         ∂     x   _         ∂   σ            Ω   _            ∂     x   _         ∂   χ         }     =       Ω                 r          ∂   r       ∂   χ              ∂   θ       ∂   σ         -     Ω                 r          ∂   r       ∂   σ                ∂   θ       ∂   χ       .                 (   50   )                                
     The relation (37), that defines the natural coordinate σ, may be expressed as                      (       ∂   x       ∂   σ       )     2     +       (       ∂   r       ∂   σ       )     2     +       (     r          ∂   θ       ∂   σ         )     2       =     h   2       ,           (   51   )                         ∂   x       ∂   σ              ∂   x       ∂   ψ         +         ∂   r       ∂   σ              ∂   r       ∂   ψ         +       (     r          ∂   θ       ∂   σ         )          (     r          ∂   θ       ∂   ψ         )         =   0     ,                                   ∂   x       ∂   σ              ∂   x       ∂   χ         +         ∂   r       ∂   σ              ∂   r       ∂   χ         +       (     r          ∂   θ       ∂   σ         )          (     r          ∂   θ       ∂   χ         )         =   0.                                            
     Inserting (48) to (50) in (38) and (39) leads to a complete set of equations in cylindrical polar coordinates, if the equations (51), that define the natural coordinate, and the thermodynamic relations are added. 
     7. Cartesian Coordinates as Dependent Variables 
     In most cases, including inverse design of wings, for example, the coordinate vector  x  is preferably expressed in terms of Cartesian coordinates and the frame of reference is not rotating. Hence, we can set  Ω = 0 . Introducing Cartesian coordinates (x,y,z) Jacobi&#39;s determinant assumes the form              J   =           ∂   x       ∂   σ              ∂   y       ∂   ψ              ∂   z       ∂   χ         +         ∂   x       ∂   ψ              ∂   y       ∂   χ              ∂   z       ∂   σ         +         ∂   x       ∂   χ              ∂   y       ∂   σ              ∂   z       ∂   ψ         -         ∂   x       ∂   σ              ∂   y       ∂   χ              ∂   z       ∂   ψ         -         ∂   x       ∂   ψ              ∂   y       ∂   σ              ∂   z       ∂   χ         -         ∂   x       ∂   χ              ∂   y       ∂   ψ              ∂   z       ∂   σ                   (   52   )                                
     The relation (37), that defines the natural coordinate σ, may be expressed as                      (       ∂   x       ∂   σ       )     2     +       (       ∂   y       ∂   σ       )     2     +       (       ∂   z       ∂   σ       )     2       =     h   2       ,           (   53   )                         ∂   x       ∂   σ              ∂   x       ∂   ψ         +         ∂   y       ∂   σ              ∂   y       ∂   ψ         +         ∂   z       ∂   σ              ∂   z       ∂   ψ           =   0     ,                                   ∂   x       ∂   σ              ∂   x       ∂   χ         +         ∂   y       ∂   σ              ∂   y       ∂   χ         +         ∂   z       ∂   σ              ∂   z       ∂   χ           =   0.                                            
     8. Isentropic Flows 
     As pointed out earlier, the equations (38) and (39) can be substantially simplified if we restrict the consideration to isentropic flows. In this case the conditions (40) lead to                                ∂     ∂   ψ            (       h   2       ρ                 J       )       =       0   ⇒       h   2       ρ                 J         =       f   1          (     σ   ,   χ     )                         ∂     ∂   χ            (       h   2       ρ                 J       )       =       0   ⇒       h   2       ρ                 J         =       f   2          (     σ   ,   ψ     )                 }     ⇒       h   2       ρ                 J         =     f        (   σ   )         ,           (   54   )                                
     where ƒ 1  and ƒ 2  are functions of integration. Moreover, as we have not yet defined the scaling of the natural coordinate σ the function ƒ(σ) can be “absorbed” by re-scaling the natural coordinate, 
     
       
         ƒ(σ)dσ→dσ.  (55)  
       
     
     For convenience we may also assume that the normalization of the density is included in this stretching of the natural coordinate and from now on consider ρ as the density that is normalized with the stagnation density. Hence we can write 
     
       
         h 2 =ρJ.  (56)  
       
     
     and              ρ   =       ρ        (   M   )       =         [     1   +         γ   +   1     2          M   2         ]         -   1       γ   -   1         .               (   57   )                                
     On the other hand, by suitable scaling of the stream functions, the stagnation sound speed can also be absorbed and (47) may be normalized accordingly,                  h   2       J   2       =         m   2          (   M   )       =         M   2         [     1   +         γ   -   1     2          M   2         ]         γ   +   1       γ   -   1           .               (   58   )                                
     Combining (36) and (56) yields                  1   ρ            ∂     x   _         ∂   σ         =       [       ∂     x   _         ∂   ψ       ]     ×       [       ∂     x   _         ∂   χ       ]     .               (   59   )                                
     Eliminating J and h from (56), (58) and the first equation (37) leads to                    ∂     x   _         ∂   σ       ·       ∂     x   _         ∂   σ         =         (     ρ   m     )     2     .             (   60   )                                
     We may now regard (59) as the equation for  x = x (σ,ψ,χ), where the density ρ is related to |∂ x /∂σ| by (57), (58) and (60). In other words, (57) and (58) may be considered as a parametric representation (in terms of the Mach number) of the function m=m(ρ) that can be inserted in (60) to obtain the desired expression for the density. 
     The structure of the equations becomes still more transparent by considering the special case of plane two-dimensional flows. This simplification is achieved, for example, by introducing Cartesian coordinates,  x =(x,y,z) and defining the stream function χ by 
     
       
         χ=z.  (61)  
       
     
     In this case (59) can be reduced to                    ∂   x       ∂   σ       =     ρ          ∂   y       ∂   ψ           ,         ∂   y       ∂   σ       =       -   ρ            ∂   x       ∂   ψ                   (   62   )                                
     and (57), (58) and (60) can be put into the form                      (       ∂   x       ∂   σ       )     2     +       (       ∂   y       ∂   σ       )     2       =       (     ρ   m     )     2       ,     ρ   =     ρ        (   M   )         ,     m   =       m        (   M   )       .               (   63   )                                
     The Jacobian reduces to              J   =           ∂   x       ∂   σ              ∂   y       ∂   ψ         -         ∂   x       ∂   ψ                ∂   y       ∂   σ       .                 (   64   )                                
     Combining (56), (58), (62) and (64) it may be shown that                    (       ∂   x       ∂   ψ       )     2     +       (       ∂   y       ∂   ψ       )     2       =       1     m   2       .             (   65   )                                
     Considering now an incremental stream tube it is easy to shown that the equations (62), (63) and (65) are in agreement with the equations of motion. Eliminating ρdσ from (62) and (63) we obtain the correct relation between massflow vector and velocity vector. On the other hand (65) relates the massflow to the cross-sectional area of the incremental stream tube. The massflow m reaches its maximum at M=1. Hence, the term on the right-hand side of (65) reaches its minimum at M=1. Considering an incremental Laval nozzle with slowly varying width dy along a line y=constant, for example, the convergent-divergent character of the cross-section is immediately apparent from (65). Nevertheless, it is surprising that the system of equations (62) is generally elliptic. It appears that the hyperbolic character of supersonic flow is hidden in the combination of (62) and (63) with ∂m/∂M&lt;0. 
     Finally, for the special case of incompressible flows (62) reduces to Laplace&#39;s equation,                    ∂   x       ∂   σ       =       ∂   y       ∂   ψ         ,         ∂   y       ∂   σ       =         -       ∂   x       ∂   ψ         ⇒           ∂   2        x       ∂     σ   2         +         ∂   2        x       ∂     ψ   2             =   0       ,             ∂   2        y       ∂     σ   2         +         ∂   2        y       ∂     ψ   2           =   0     ,           (   66   )                                
     which is in agreement with the properties of a special type of conformal mapping that is usually called hodograph transformation. 
     9. General Limitations with Respect to Supersonic Flow 
     An important limitation that applies to any coordinate inversion and is relevant for supersonic flow appears if the Jacobian vanishes, 
     
       
         J=0  (67)  
       
     
     Flows that satisfy the condition (67) are not included in the solution space of the inverse equations. The reader is referred to the corresponding discussion by Courant and Friedrichs [1]. Flows that satisfy (67) are called simple waves. An expansion fan in an isentropic flow, for example, is not a solution of the equations (57) to (60). To obtain a complete set of solutions we need to add the solutions of (67) to the solution space of any set of inverse equations. In Cartesian coordinates, for example, the general solution of (67) can be expressed as 
     
       
         G(x,y,z)=0.  (68)  
       
     
     The validity of (68) is easily confirmed by introducing (52) in (67) and inserting (68). 
     10. Axisymmetric Flows 
     For the special case of axisymmetric flow in the absolute frame of reference the equations (19) and (20) reduce to a generalized form of the Bragg-Hawthorne equation. Without giving the details of deriving the equation we simply state the results. For axisymmetric flows we introduce a single stream function ψ and the specific angular momentum C as an additional integral of motion:                  ρ                 u     =       1   r            ∂   ψ       ∂   r           ,       ρ                 v     =       -     1   r              ∂   ψ       ∂   x           ,     C   =       C        (   ψ   )       =   rw       ,           (   69   )                                
     where x and r are the axial and radial coordinates, and u, v and w denote the axial, radial and azimuthal velocity components, respectively. In this case the Euler equation is found to be                      1     ρ                 r              ∂     ∂   r            [       1   pr            ∂   ψ       ∂   r         ]         +       1     pr   2              ∂     ∂   x            [       1   ρ            ∂   ψ       ∂   x         ]           =            H          ψ       -     T             s          ψ         -       C     r   2                 C          ψ             ,           (   70   )                                
     where the total enthalpy 
     
       
         H=H(ψ)=c p T+ ½[u   2 +v 2 +w 2 ]  (71)  
       
     
     and the entropy s=s(ψ) are again integrals of motion. Equation (47) reduces to                        ρ   2          [       u   2     +     v   2     +     w   2       ]       =                    [       1   r            ∂   ψ       ∂   x         ]     2     +       [       1   r            ∂   ψ       ∂   r         ]     2     +                                  C   2       r   2              ρ   0   2         [     1   +         γ   -   1     2          M   2         ]       2     γ   -   1                         =                    ρ   o   2          c   o   2          M   2           [     1   +         γ   -   1     2          M   2         ]         γ   +   1       γ   -   1                         (   72   )                                
     The remaining equations are the same as in section 5. We may note that, regarding the density ρ and the entropy s as constants, replacing the total enthalpy (71) by the total head,                H   =       H        (   ψ   )       =       p   ρ     +       1   2          [       u   2     +     v   2     +     w   2       ]             ,           (   73   )                                
     and eliminating the density from the definition of the stream function,                u   =       1   r            ∂   ψ       ∂   r           ,     v   =       -     1   r              ∂   ψ       ∂   x           ,           (   74   )                                
     (70) reduces to the well-known Bragg-Hawthorne equation for incompressible axisymmetric flow:                      ∂   2        ψ       ∂     r   2         -       1   r            ∂   ψ       ∂   r         +         ∂   2        ψ       ∂     x   2           =         r   2               H          ψ         -     C               C          ψ       .                 (   75   )                                
     At this point we introduce the mapping 
     
       
         (x,r)←→(σ,ψ)  (76)  
       
     
     to inverse coordinates and, for simplification, we restrict the consideration to incompressible flow. In this case the streamline coordinate σ is defined by                  h          ∂   σ       ∂   x         =         ∂   ψ       ∂   r               (       ∂   ψ       ∂   x       )     2     +       (       ∂   ψ       ∂   r       )     2             ,           (   77   )                   h          ∂   σ       ∂   r         =       -       ∂   ψ       ∂   x                 (       ∂   ψ       ∂   x       )     2     +       (       ∂   ψ       ∂   r       )     2             ,                                            
     where h=h(σ,ψ) is the stretching function that provides the basis for compatibility between orthogonality of the curves σ(x,r)=constant and ψ(x,r)=constant and the uniform meaning of hdσ.                  h                 d                 σ     =               (   dx   )     2     +       (   dr   )     2         ⇒         (       ∂   x       ∂   σ       )     2     +       (       ∂   r       ∂   σ       )     2         =     h   2         ,           (   78   )                                
     as the length increment along streamlines. The condition for the curvilinear coordinates to be orthogonal is again                      ∂   x       ∂   σ              ∂   x       ∂   ψ         +         ∂   r       ∂   σ              ∂   r       ∂   ψ           =   0.           (   79   )                                
     Introducing Jacobi&#39;s determinant,                J   ≡           ∂   x       ∂   σ              ∂   r       ∂   ψ         -         ∂   r       ∂   σ              ∂   x       ∂   ψ             ,           (   80   )                                
     the partial derivatives of σ(x,r) and ψ(x,r) can be expressed as                    ∂   σ       ∂   x       =       1   J            ∂   r       ∂   ψ           ,         ∂   σ       ∂   r       =       -     1   J              ∂   x       ∂   ψ           ,         ∂   ψ       ∂   x       =       -     1   J              ∂   r       ∂   σ           ,         ∂   ψ       ∂   r       =       1   J              ∂   x       ∂   σ       .                 (   81   )                                
     For a function ƒ(x,r)=g(σ,ψ) the conversion of derivatives can be expressed as                    ∂   f       ∂   x       =       1   J          [           ∂   r       ∂   ψ              ∂   g       ∂   σ         -         ∂   r       ∂   σ              ∂   g       ∂   ψ           ]         ,         ∂   f       ∂   r       =     -         1   J          [           ∂   x       ∂   ψ              ∂   g       ∂   σ         -         ∂   x       ∂   σ              ∂   g       ∂   ψ           ]       .                 (   82   )                                
     Inserting (81) in (82) we obtain the conversion of the second derivatives:                      ∂   2        ψ       ∂     x   2         =       1   J          [           ∂   r       ∂   σ            ∂     ∂   ψ            (       1   J            ∂   r       ∂   σ         )       -         ∂   r       ∂   ψ            ∂     ∂   σ            (       1   J            ∂   r       ∂   σ         )         ]         ,           (   83   )                     ∂   2        ψ       ∂     r   2         =         1   J          [           ∂   x       ∂   σ            ∂     ∂   ψ            (       1   J            ∂   x       ∂   σ         )       -         ∂   x       ∂   ψ            ∂     ∂   σ            (       1   J            ∂   x       ∂   σ         )         ]       .                                              
     Adding the second derivatives according to (83), and making use of (78), (79) and (80), integration by parts yields                      ∂   2        ψ       ∂     x   2         +         ∂   2        ψ       ∂     r   2           =       1   J          ∂     ∂   ψ              (       h   2     J     )     .               (   84   )                                
     Furthermore, from (78), (79) and (80) it follows that                  ∂   x       ∂   σ       =             h   2     J            ∂   r       ∂   ψ         ⇒       1   r            ∂   ψ       ∂   r           =         1   rJ            ∂   x       ∂   σ         =         h   2       rJ   2                ∂   r       ∂   ψ       .                   (   85   )                                
     Combining (75) with (84) and (85) we finally obtain the inverse Bragg-Hawthorne equation                  r   J          ∂     ∂   ψ            (       h   2     rJ     )       =         r   2               H          ψ         -     C               C          ψ       .                 (   86   )                                
     A complete system of equations is obtained by added (78) and (79). The equation (86) is again remarkably simple. To check the previous results a somewhat more elaborate analysis has been chosen in this case to derive the inverse equation. It is tedious, but rather straightforward to eliminate the Jacobian from (86) in favor of partial derivatives of the radius. Making use of (78), (79) and (80) the equation (86) can be written in the form                      (       ∂   r       ∂   ψ       )     2              ∂   2        r       ∂     σ   2           +         ∂   r       ∂   σ              ∂   r       ∂   ψ                ∂   2        r         ∂   σ          ∂   ψ           +       [       h   2     -       (       ∂   r       ∂   σ       )     2       ]          [           ∂   2        r       ∂     ψ   2         +       1   r            (       ∂   r       ∂   ψ       )     2         ]       -       1   h            ∂   r       ∂   σ                ∂   r       ∂   ψ            [           ∂   r       ∂   σ              ∂   h       ∂   ψ         +         ∂   r       ∂   ψ              ∂   h       ∂   σ           ]         +           h   2          (       ∂   r       ∂   ψ       )       3          {         r   2               H          ψ         -     C             C          ψ           }         =   0           (   87   )                                
     and the derivative of the coordinate stretching factor with respect to the stream function can be expressed as                  ∂   h       ∂   ψ       =             ∂   r       ∂   ψ           h   2     -       (     ∂     ∂   σ       )     2              [           ∂   h       ∂   σ              ∂   r       ∂   σ         -     h            ∂   2        r       ∂     σ   2             ]       .             (   88   )                                
     The idea is now to solve (87) by means of a Gauss-Seidel iteration scheme and determine the coordinate stretching factor h using a stepwise integration of (88). From (74), (78), (81) and (85) we obtain                        u   2     +     v   2       =                  1     r   2            [         (       ∂   ψ       ∂   x       )     2     +       (       ∂   ψ       ∂   r       )     2       ]                   =                    1       r   2          J   2              [         (       ∂   x       ∂   σ       )     2     +       (       ∂   r       ∂   σ       )     2       ]       =         h   2         r   2          J   2         =         h   2     -       (       ∂   r       ∂   σ       )     2           h   2              r   2          (       ∂   r       ∂   ψ       )       2                           (   89   )                                
     The flow speed can now be expressed as                    u   2     +     v   2     +     w   2         =           u   2     +     v   2     +       (     C   r     )     2         =               h   2     -       (     ∂     ∂   σ       )     2           h   2              r   2          (       ∂   r       ∂   ψ       )       2         +       (     C   r     )     2         .               (   90   )                                
     11. Difficulties Encountered with Applications of Inverse Equations 
     For the transform of the equation to their inverse form it is essential to choose the domain of definition such that a one-to-one mapping is generated. For simply connected domains of definition, such as turbine blade channels, for example, this difficulty would generally not appear. However, in the case of flow fields that are not simply connected the difficulty arises that the entire flow field must be suitably cut to obtain domains for which a one-to-one mapping exists. A further point that should be addressed concerns the chosen domain of definition for the stream functions. In this context the reader is referred to Keller [12] for more details. If a wake appears within the domain of definition there may be a problem with streamlines, originating on the downstream side, for which the integrals of motion have not been defined. For this reason it is a good idea to include a check in a numerical solution procedure that warns the user in cases where integrals of motion have been taken beyond their intended domain of definition. 
     12. Two Simple Examples 
     To demonstrate the capabilities of inverse equations we consider two simple axisymmetric problems. The design of a compressor diffuser of a gas turbine is often strongly restricted by the space requirements of turbine and combustor. As a consequence there is usually the necessity of choosing a curved diffuser. Although there is broad empirical knowledge on the performance of straight annular diffusers (see Kline, Abbott and Fox [13], for example), it may often be difficult to apply this information to arbitrary curved diffusers. However, knowing the relation between the largest acceptable adverse pressure gradient and the ratio of (local) length and inlet width of a diffuser we can specify restrictions for the static pressure distributions along the inner and outer walls of a curved annular diffuser. As indicated in FIG. 1, a possible strategy to design a curved diffuser may be defined as follows. Near the diffuser inlet the flow speeds along both walls first follow a curve the corresponds to the largest acceptable pressure gradient. After a certain arc length a smooth transition to a different curve of maximum pressure increase is chosen for the outer diffuser wall. For a certain distance the flow speed distributions along the two walls follow two different curves of largest pressure increase. This domain will generate the curved portion of the diffuser. Toward the exit of the diffuser a smooth transition of the flow speed along the inner wall is introduced, in order to obtain a velocity profile at the diffuser exit that is as uniform as possible. The second transition beings the flow speed along the inner wall back to a common curve of maximum pressure increase. Knowing the flow speed distributions along both diffuser walls enables us to express ∂r/∂σ in terms of ∂r/∂ψ, making use of equation (90). At the inlet to the diffuser we set ∂r/∂σ=0, and at the exit of the diffuser we choose a linear transition for ∂r/∂σ between the corresponding values at the endpoints of the two walls. Choosing any kind of distributions for the total head and the circulation at the diffuser inlet does not cause additional difficulties, because H=H(ψ) and C=C(ψ) appear as integrals of motion in the governing equations (87) and (88). It is now easy to find a numerical solution with the help of the Gauss-Seidel iteration scheme. For more details the reader is referred to Keller et al. [10]. For the simple example presented it has been assumed that the flow is incompressible, the total head is uniform and the flow is free of swirl. FIG. 2 illustrates a streamline map that corresponds to the boundary conditions defined by the flow speed distributions shown in FIG.  1 . Furthermore, FIG. 2 shows lines of constant σ. It should be pointed out that, according to the definition of the natural coordinate, lines of constant σ are always orthogonal to streamlines (i.e., lines of constant ψ and χ). Thus we may observe that one of the remarkable spin-offs of using inverse equations is the automatic generation of an orthogonal contour-fitted coordinate system in the original space. Although the actual computation is done in a cubic domain of definition with Cartesian coordinates in inverse space, a curvilinear orthogonal grid is simultaneously generated together with the solution in the original space. Furthermore, the density of grid lines gets automatically high in domains of large gradients. For these reasons the accuracy of the method is remarkably high (for a fixed number of grid points) and the convergence is very fast. 
     The second example (see FIG. 3) is a turbine inlet flow. In this case the geometry of an annular duct is imposed as a boundary condition, combined with a linear variation of the flow angles at the inlet and exit, respectively. According to ABB&#39;s convention the flow is from right to left. This is one of the examples that has been used to check the equations. To include and discuss three-dimensional examples of inverse design would exceed the scope of this paper. Such results will be presented elsewhere. 
     13. Concluding Remarks 
     With the transformation of the Euler equation to their inverse form, using a pair of stream functions and the natural coordinate as independent variables, a theoretical basis has been introduced for a broad class of extremely powerful inverse computation methods. Instead of searching for a geometry of a flow device that will lead to certain physical properties, the physical properties can be defined in advance and the solution of the inverse equations directly leads to the corresponding geometry. It has become easier to arrive at an optimized design. However, the difficulty of defining an optimum remains. 
     From a practical point of view it may be convenient to proceed as follows. As a preparation for the optimum design of a certain flow device we introduce a first guess for its geometry, that may satisfy certain basic constraints. Secondly we use the inverse method for a “backward computation”. In other words, we specify the geometry according to the first guess and determine the static pressure distribution along the surface. In a second step we modify the static pressure distribution along the surface according to the requirements of the optimization and introduce the modified static pressure distribution as the boundary condition for the “forward computation”. It should be noted that the inverse equations provide an excellent basis for using both physical quantities or the geometry of a flow device to define the boundary conditions. In fact it is even easy to use combinations of geometrical and physical boundary conditions. 
     Hence, it seems that there is justified hope that computation methods that are based on inverse equations will become extremely efficient tools to arrive at optimized flow devices. 
     Although it might appear that the existence of integrals of motion is essential to derive a set of antisymmetric equations of the form (19), (29). However, it is interesting that the extension to non-adiabatic flows of both antisymmetric equations and transformation to inverse coordinates is relatively simple, including the case of the general time-independent Navier-Stokes equations. 
     References 
     [1] R. Courant and K. O. Friedrichs, Supersonic flow and shock waves, Interscience, New York (1967). 
     [2] J. J. Keller, A pair of stream functions for three-dimensional vortex flows, Z. angew. Math. Phys. 47 (1996), 821. 
     [3] A. Jameson, N. A. Pierce, and L. Martinelli, Optimum aerodynamic design using the Navier-Stokes equations, AIAA paper 97-0101, 35th Aerospace Sciences Meeting, Reno, Nev. 1997. 
     [4] M. B. Giles and M. Drela, Two-dimensional transonic aerodynamic design method, AIAA Journal 25 (1987), 1199 
     [5] T. Dang and V. Isgro, Inverse method for turbomachine blades using existing time-marching techniques, ASME paper 94-GT-20, International Gas Turbine and Aeroengine Congress, The Hague, Netherlands, 1994. 
     [6] G. S. Dulikravich, Shape inverse design and optimization for three-dimensional aerodynamics, AIAA paper 95-0695, 33rd Aerospace Sciences Meeting, Reno, Nev. 1995. 
     [7] M. J. Pandya and O. Baysal, Gradient-based aerodynamic shape optimization using ADI method for large-scale problems, AIAA paper 96-0091, 34th Aerospace Sciences Meeting, Reno, Nev. 1996. 
     [8] I. K. Jennions, The role of CFD in the design process, AGARD Paper, 7 Rue Ancelle, 92200 Neuilly sur Seine, France, Lecture Series 195, May to June 1994. 
     [9] M. G. Turner and M. B. Giles, Design and analysis of internal flow fields using a two stream function formulation, ASME Winter Annual Meeting 1990, FED Vol. 103, 203. 
     [10] J. J. Keller, W. Egli, and J. Exley, Force- and loss-free transitions between flow states, Z. angew. Math. Phys. 36 (1985), 854. 
     [11] J. J. Keller, W. Egli, and R. Althaus, Vortex breakdown as a fundamental element of vortex dynamics, Z. angew. Math. Phys. 39 (1988), 404. 
     [12] J. J. Keller, On the interpretation of vortex breakdown, Phys Fluids 7 (7), (1995), 1696. 
     [13] S. J. Kline, D. E. Abbott and R. W. Fox, Optimum design of straight-walled diffusers, J. Basic Eng., Trans. ASME 81, Series D (September 1959). 
     While the invention has been described in detail with reference to preferred embodiments thereof, it will be apparent to one skilled in the art that various changes can be made, and equivalents employed, without departing from the scope of the invention.