Three dimensional flow simulations with the finite element technique over a multi-stage rocket

TECHNICAL PAPERS

Three dimensional flow simulations with the finite element technique over a multi-stage rocket

L. C. Scalabrin^I; J. L. F. Azevedo^I; P. R. F. Teixeira^II; A. M. Awruch^III

^IInstituto de Aeronáutica e Espaço Centro Técnico Espacial CTA/IAE/ASE-N 12228-904 São José dos Campos, SP. Brazil [e-mail] [e-mail]

^IIDepartamento de Materiais e Construção Fundação Universidade do Rio Grande-FURG Av. Itália, Km 8, s/nº Campus Carreiros 96206-900 Rio Grande, RS. Brazil [e-mail]

^IIIPrograma de Pós-Graduação em Eng. Mecânica Universidade Federal do Rio Grande do Sul PROMEC/UFRGS Av. Osvaldo Aranha,99, 3º. Andar 90035-190 Porto Alegre, RS. Brazil [e-mail]

ABSTRACT

Aerodynamic flow simulations over the first Brazilian satellite launch vehicle, VLS, during its first-stage flight are presented. The three dimensional compressible flow is modeled by the Euler equations and a Taylor-Galerkin finite element method with artificial dissipation is used to obtain the numerical solution. Transonic and supersonic results for zero angle-of-attack are presented and compared to available experimental results. The influence of mesh refinement and artificial dissipation coeffcient on the transonic flow results are discussed. The results obtained for the supersonic simulations present good agreement with experimental data. The transonic simulation results capture the correct trends but they also indicate that this flight condition requires more refined meshes.

Keywords: Aeronautical engineering, computational fluid dynamics, finite element technique, Taylor-Galerkin method, transonic flows, supersonic flows

Introduction

In the present work, the results obtained for the simulation of aerodynamic flows over the first Brazilian Satellite Launcher, VLS, on its first-stage flight configuration are presented. The VLS, in this configuration, is composed by four strap-on booster arranged symmetrically around the central core. Three-dimensional grid generation for this configuration is complex and the computational fluid dynamics group at Centro Técnico Aeroespacial, CTA, decided to approach this problem in two different ways: with unstructured tetrahedral meshes and with multiblock hexahedral structured meshes grouped together through the Chimera technique. The group already obtained very good results with the Chimera procedure (Antunes, Basso and Azevedo, 2001; Basso, Antunes and Azevedo, 2000a, 2000b) and the present work shows the progress obtained with the unstructured grid approach. Furthermore, there was also interest in using this opportunity to develop the capability of performing such simulations using the finite element method, which was not yet explored in the group.

To rapidly reach this objective, a partnership with another Computational Fluid Dynamics, CFD, group in Brazil, from Universidade Federal do Rio Grande do Sul, UFRGS, was established and an existing finite element code was transfered to Instituto de Aeronáutica e Espaço, IAE, to start the project. This code, which had been developed at UFRGS as the graduate work of one of the authors (Teixeira, 1996), is intended to take full advantage of vector processing computers and this test over such complex geometry was an important opportunity for the continuation of its validation process.

The code is based on the Taylor-Galerkin finite element method, which is a solution proposed by Donea (1984) for the problems faced by the finite element community (Chung, 1978; Heinrich et al., 1977; Christies et al., 1976) to solve flow problems in a wide range of velocities (Löhner, Morgan and Zienkiewicz, 1985). This method is an extension of the Lax and Wendroff (1960; 1962; 1964) concepts and it is basically a Taylor series expansion of the solution where the temporal derivatives are approximated by the spatial derivatives according to the partial differential equations that govern the problem at hand. Despite the fact that the Lax-Wendroff method adds 4th difference artificial dissipation terms, it is necessary to add 2nd difference terms to stabilize the solution near shocks and, in this work, this is performed following the approach proposed on Peraire et al. (1988) and Morgan et al. (1991). The spatial derivatives are discretized by the Galerkin weighted residuals method (Zienkiewicz and Morgan, 1983) using linear and constant shape functions and the weak formulation. The boundary conditions are quite simple in this work and some diffculties appear only on subsonic flow entrance and exit boundaries, where the flux vector splitting technique is implemented (Warming, Beam and Hyett, 1975; Hughes and Tezduyar, 1984).

The paper first presents the theoretical formulation of the equations and the numerical method, followed by the boundary condition formulations. Next, the results and discussion are presented with experimental comparisons.

Nomenclature

a = speed of sound

A = Jacobian matrix

CC =artificial dissipation coeffcient

CFL = stability limit number

D = artificial dissipation term

e = total energy

E , F , G = inviscid flux vectors

l_k = unit normal vector components

L = diagonal matrix

M =mass matrix

N =linear shape function

p =static pressure

P = constant shape function

Q = conserved variables vector

R , S = finite element method integrals

t = time

u , v , w = flow velocity components

x , y , z = cartesian coordinates

X = eigenvector matrix

Greek Symbols

r = density

g = ratio of specific heats

Dt = time step

W = volume

G = face, face area n numerical pressure sensor

l = eigenvalues

Subscript

¥ = freestream variable

0 = reference variable

e = element variable

G = face variable

L = lumped

Superscript

n = time step counter

Theoretical Formulation

A first approach to solve the aerodynamic problem over the VLS is to assume the flow to be compressible and inviscid, i.e., it can be modeled by the Euler equations. The dimensionless conservative form of the Euler equations is

where is the dimensionless time, are the dimensionless Cartesian coordinates, x₀ is the reference length, and is the vector of dimensionless conserved variables, defined as

In these equations, is the dimensionless density, and are the dimensionless Cartesian velocity components, is the total energy per unit of mass, a is the speed of sound and the subscript ¥ refers to freestream quantities.

The , and are the dimensionless inviscid flux vectors, which can be written as

In the previous expression, the dimensionless pressure can be obtained from the equation of state for a perfect gas as

where g is the fluid ratio of specific heats. In the interest of simplifying the nomenclature, the bar that indicates a dimensionless variable will be dropped in the forthcoming equations with the understanding that, unless otherwise stated, all variables were made dimensionless as described.

Finite Element Approximation

In this work, the Taylor-Galerkin finite element method was used to obtain the numerical solution of the system of equations previously presented. The Taylor-Galerkin method can be viewed as a Lax-Wendroff method (Lax and Wendroff, 1960, 1962, 1964) for the time march with the spatial derivatives discretized with the finite element concepts. The Lax-Wendroff method for the Euler equations can be written as a predictor-corrector, time-marching procedure as

where

is the spatial part of the vectorial form of Eq. (1). These terms have to be discretized and this was accomplished by adopting a tetrahedral grid and using the finite element technique.

In the first stage of the predictor-corrector process the variables are discretized with constant element shape functions, in order to keep this stage explicit. The terms in the beginning of the stage are discretized with linear node shape functions and, consequently, the result of the second stage of the process has to be discretized with this type of shape function too. The spatial discretization expressions for the first stage of the method are given by

where nnm and nem are, respectively, the total number of nodes and elements in the grid, N_j is the linear shape function for node j and P_e is the constant shape function for element e. The result of the first stage of the predictor-corrector method, after using the spatial discretization, is given by

where the j-th summation is taken over the 4 nodes that define the E-th tetrahedral element and the Dt is the lowest time interval in the mesh, calculated using the CFL number (Strauss and Azevedo, 1999) and the characteristic length of the tetrahedral element (Scalabrin, 2002).

In the second stage, the necessary spatially discretized terms are

where N_j and P_e are as described before. Using this discretization and the weak formulation usually taken in finite elements for this class of problem (Donea, 1984), one can obtain the expression for the second stage of the predictor-corrector Taylor-Galerkin method in its matrix form as (Teixeira, 1996)

where

and are the components of the unit normal vector to the l' _e face. It is clear that this stage of the method involves the M matrix inversion to obtain the Qⁿ⁺¹ vector. In this work, the Donea and Giuliani (1981) approach was utilized to solve for Qⁿ⁺¹ and an iterative form for the solution of Eq. (12) was established as follows

In the above equation, p is the iteration counter and M_L Hinton lumped mass matrix which is calculated as

where the symbol d_ei indicates that the summation had to be taken only on the elements which contain the node i.

The flux terms in the R integral are calculated using

and the similar terms in the S integral are approximated to the surface using the ideas proposed by Argyris, Doltsinis and Friz (1989), i.e.,

The term is obtained from the node values and node shape functions on the face, is an arithmetic mean of flux node values over the nodes that define the e-th element. The term , for elements that do not have faces on flow entrances or exits is calculated according to Eq. (16). The calculation of this term for elements on entrances and exits will be discussed in the boundary condition section.

It is well known that the Lax-Wendroff method implicitly adds 4th-difference dissipation terms to the problem, but this dissipation is not enough to overcome the instability generated in the vicinity of shocks. In this work, explicit 2nd-difference dissipation was added after the second stage of the predictor-corrector process to stabilize the numerical method, as suggested by Morgan et al. (1991) and Peraire et al. (1988). This 2nd-order term is turned on only where it is necessary, i.e., near the shocks, which is achieved by the use of a pressure sensor. The expression for the pressure sensor is given by

where the maximum function is taken over the 4 nodes that define the e-th tetrahedral element. The damping vector itself, which is responsible for adding the artificial dissipation to the problem, is

where Dt is the lowest time step found in the grid, Dt_e is the maximum possible time step for the e-th element that still keeps the stability of the numerical procedure and CC is a constant defined by the program user. This damping vector is added to the solution after the second stage of the predictor-corretor process by

Initial and Boundary Conditions

A well posed mathematical problem needs boundary and initial conditions. In this work, freestream conditions are assumed as initial condition in the whole field, except on the VLS surface. Along the vehicle surface, velocity is forced to be tangent to the wall, whereas r and e also assume freestream values. The boundary conditions, however, demand more detailed explanations.

In the VLS problem, there are basically four types of boundary conditions: wall, symmetry and flow entrance and exit conditions. The symmetry boundary condition, in an inviscid flow, is established in the same fashion as the wall boundary condition. The wall boundary condition, in an Euler context, states that the flow must be tangent to the wall, or mathematically

The wall and symmetry boundary conditions are applied to the nodes at the end of the second stage of the method. This can be achieved by decomposing the velocity vector in a local coordinate system with the x_l axis normal to the boundary node, setting the u_x component in this local system equal to zero and using the other components to return to the global coordinate system. The node normal axis is calculated as an arithmetic average of the normal axis of the boundary faces which share the node. The calculation of the node normal axis for nodes which are in symmetry axes has some details which can be found in Teixeira (1996).

The flow entrance and exit boundary conditions in this work were modeled using the flux vector splitting technique (Warming, Beam and Hyett, 1975; Steger and Warming, 1981) that takes advantage of the homogeneous property of the Euler flux vectors. This property states that

where dF_k/dQ clearly represents the flux Jacobian matrix, A_k, in the k-th direction. As the main problem here concerns the boundary conditions, it is more adequate to work with the projection of the flux vector on the direction normal to the face, i.e.,

Here, are the components of the unit vector normal to the face.

The matrix has a complete set of eigenvalues and, according to Warming, Beam and Hyett (1975) and to Hughes and Tezduyar (1984), it can be diagonalized by

where . Furthermore, the eigenvector X matrix can be formed as

where

and

The eigenvalues of the matrix can be viewed as a sum of positive and negative eigenvalues which define the directions of propagation of the characteristics. Therefore, one can decompose the diagonal L matrix as

where

and . Hence, one can rearrange the matrix in terms of the positive and negative set of eigenvalues and, using the homogeneous property of the flux vectors of the Euler equations, the flux vectors can also be divided in two sets according to the sign of the eigenvalues (Steger and Warming, 1981). Therefore, one can write

where

After the definition of the split flux vectors, one can set the boundary conditions over flow entrance and exit boundaries in accordance with the appropriate characteristic directions. The boundary face is assumed subsonic if the flow in the element which contains this face is subsonic. The calculation of the Mach number in the element uses the properties based on elements, i.e., the properties obtained after the first stage of the Taylor-Galerkin method. In subsonic flow entrances, the flux on the inlet surface is calculated using

where Q is the vector of conserved variables for the element at the boundary and Qµ is the vector of conserved variables for the freestream. Similarly, for subsonic flow exit boundaries, one has

The above expressions are used to evaluate the in Eq. (17) for the elements which have faces on subsonic flow entrance or exit boundaries. If the flow in the element is supersonic, the boundary face is considered supersonic. For supersonic flow entrance and exit boundaries, these calculations are not necessary, because all characteristics propagate in one direction. Hence, for a supersonic entrance, one can write

and, for a supersonic exit,

Results and Discussion

The present work discusses results for the VLS first-stage flight configuration at zero angle-of-attack and freestream Mach numbers M =0.9 and M =2.0. These flight conditions were chosen considering the availability of experimental results, the dominant supersonic characteristic of the vehicle flight and the numerical difficulties that arise at transonic simulations. Two different tetrahedral unstructured grids were used in the simulations and both grids took advantage of the symmetry of the problem. This means that only one-quarter of the complete geometry has been included in the simulation, since only zero angle-of-attack cases are considered.

Mesh 1 has 581,000 elements and 106,000 nodes, and one longitudinal plane of it can be visualized in Fig. 1. The afterbody nozzles were omitted, as the flow over this part of the configuration cannot be correctly described by an Euler formulation (Strauss and Azevedo, 1999). A more detailed view of the nose region can be seen in Fig. 2 and a view of the booster nose cap is in Fig. 3. The grid near the booster nose cap is clearly not fine enough and this has created difficulties for an adequate capturing of all the aerodynamic phenomena that occur in this region, as it will be seen later. Figure 4 presents a schematic view of the outlet plane. In this figure, the azimuthal symmetry A plane is indicated. This is the computational plane along which comparisons between numerical and experimental results are performed.

The simulation of transonic conditions using mesh 1 did not reach satisfactory results. Therefore, a finer mesh was created, with 1.16 million elements and 212,000 nodes. This more refined mesh is compared to the original mesh in Fig. 5, which shows the two meshes in the region around the central body nose cap. One can see in the figure that, at least near the body wall, mesh 2 has almost three times more points than mesh 1. Clearly, the relation between the number of grid points does not remain constant throughout the flowfield, but Fig. 5 emphasizes that mesh 2 is considerably finer than mesh 1 in the more relevant regions of the flow.

Case 1: M=2.0 and Zero Angle-of-Attack

The supersonic test case was run only using mesh 1. Furthermore, the artificial dissipation coefficient was set to 10 for the present simulation. A general view of the pressure contours for the converged solution in this case is shown in Fig. 6. A more detailed view of the flow around the forebody portion of the vehicle central core is presented in Fig. 7. In this latter figure, one can clearly see the stagnation region and the rapid expansion over the spherical nose cap just downstream of the stagnation region. One can also see a slight recompression along the forebody conical section, followed by the flow expansion on the forebody cone - payload fairing intersection. The expansion in the boattail region is also evident in Fig. 7. This last expansion region is terminated by an oblique shock wave at the boattail-afterbody cylinder intersection.

A closer view of the flow around the forebody region of one of the boosters is presented in Fig. 8. As before, the main flow features in this region are readily identified in Fig. 8, including the booster stagnation point and the detached shock wave upstream of the booster. The relevant fact concerning this region is the interaction among the detached shock waves from the boosters and their reflections, which causes the formation of a high-pressure region between the boosters and the central core. This high-pressure region forces the flow to deflect outwards, as can be seen in Fig. 9, and it causes a low pressure region downstream along the afterbody. There is a strong expansion of the flow that crosses the high pressure region, as one can see in Fig. 10. This figure shows the pressure coefficient distribution along the central body surface in a azimuthal symmetry plane between two boosters, the A plane in Fig. 4

A more detailed view of the pressure coefficient distribution on the rocket nose surface can be seen in Fig. 11. This figure indicates a very good agreement between the numerical and the experimental results for the forebody region. As one can see in Fig. 10, the agreement in the afterbody portion of the vehicle is not as good. This should come as no surprise since the complexity of the phenomena and the various interactions present in the downstream portion of the flow account for a much more difficult simulation problem than the aerodynamically clean flowfield in the vehicle forebody. Nevertheless, it is also clear from Fig. 10 that the present simulation is able to qualitatively reproduce many of the features observed in the experimental data for this downstream region of the flow. Based on the experience in the simulation of a similar problem with a different numerical technique (Antunes, Basso and Azevedo, 2001; Basso, Antunes and Azevedo, 2000a, 2000b), the authors are clearly aware that there is a need for mesh refinement in the afterbody region in order to improve the correlation between numerical and experimental data. However, as this work is part of a larger effort to achieve numerical simulation capability over realistic configurations, the authors opted to invest the available resources to develop adaptive mesh refinement routines (Scalabrin, 2002).

Case 2: M=0.9 and Zero Angle-of-Attack Using Mesh 1

The first attempt at the solution of the transonic test case also used mesh 1 and a value of the artificial dissipation coefficient (CC) of 1.0. It should be observed that this coefficient value is ten times lower than the value used on the supersonic test case and this is justified by the fact that the shock waves in a transonic flow are weaker than the shock waves which appear in a supersonic flow. An overview of the pressure coefficient contours on the XZ plane for the converged solution in this case is shown in Fig. 12. Figure 13 presents a closer view of the central core forebody portion, where one can see the stagnation region and the subsequent expansion over the spherical nose cap. The more relevant feature in this region is the supersonic zone over the cylindrical payload fairing, which is created by the expansion on the transition from the conical section and which is terminated by a shock wave. The expansion in the upstream part of the boattail region is also clear from Fig. 13 and, as this is a subsonic flow, one can see the recompression along the boattail. A slight compression on the transition from the boattail to the afterbody cylinder can be seen as well.

A detailed view of the flow over one of the boosters is presented in Fig. 14. In this figure, one can see the stagnation region on the booster nose cap, the supersonic zone created on the cylindrical booster section after the transition from the booster cone and the high pressure region present between the booster and the vehicle central body. In Fig. 14, the poor definition of the contour lines is evident. This is mostly true in the region between the booster nose cap and the vehicle core. At first, the authors attributed this problem to a lack of artificial dissipation because of the decreased value of the CC coefficient described in the beginning of this sub-section.

However, tests performed with higher values of the artificial dissipation coefficient showed that an increase of the artificial dissipation, in this case, deteriorates the results. This is illustrated in Figs. 15 and 16, in which the numerical pressure coefficient distributions along the central body calculated using mesh 1 and CC coefficients equal to 1.0 and 7.5 are presented. This fact encouraged the authors to test the mesh refinement influence over the numerical results.

Case 3: M=0.9 and Zero Angle-of-Attack Using Mesh 2

This test case was studied in order to verify the influence of the mesh refinement on the solution of the transonic case. It was run using mesh 2, which is globally two times finer than mesh 1 and, in aerodynamic relevant regions, it is almost 4 times finer than mesh 1. In this case, the artificial dissipation coefficient was set to 1.0, as the previous test with the transonic case showed that an increase in this value can cause the solution to deteriorate.

One can see the improvement in solution quality obtained with the mesh refinement comparing Figs. 13 and ¹⁷. This latter figure presents more well defined features, especially the payload fairing shock wave and the stagnation region, and smoother contour lines. For the flow over the boosters, this improvement is even more evident as on can see in Figs. 18 and 19. Here, the perturbation zone caused by the booster is smaller than the one present in the results of Fig. 14. Furthermore, the aerodynamic phenomena in this region, i.e., the stagnation zone, the high pressure region between the booster and the vehicle central core and the booster cylindrical region shock wave, are better defined.

Figures 20 and 21 show a comparison between the solutions obtained for the central body surface on the A plane using mesh 1 and mesh 2. The numerical solution using the finer mesh 2 clearly provides a better approximation of the experimental pressure curve extremes and it also has a better fitting of some trends of the pressure coefficient distribution, as one can see in Fig. 21 at

Concluding Remarks

The paper has presented results for 3-D Euler simulations of the flow over the complete first-stage flight configuration of the first Brazilian satellite launcher, the VLS. To the authors' knowledge, this is the first time that such detailed simulations of the complete VLS system have been performed with an unstructured grid approach. The results here presented also indicate the current stage of development in the group of the capability of simulating complex, three-dimensional flows using unstructured grids. The code utilized is based on the Taylor-Galerkin finite element method and the spatial discretization uses tetrahedral elements with constant element shape functions and linear node shape functions. An artificial dissipation term is added to overcome the instability generated on strong shocks and the flux vector splitting technique is implemented to help define the boundary conditions on flow entrances and exits. The code was tested on supersonic and transonic flows over the VLS and studies concerning the artificial dissipation coefficient and mesh refinement were undertaken.

The quality of the supersonic results along the rocket nose was very good even for such a fairly coarse mesh. Some problems occurred on the vehicle afterbody, and the authors believe that most of this problems can be solved with a more refined mesh in this region. However, it cannot be discarded that non-modeled physics, due to the use of the Euler equations, can be playing an important influence on the results in the region in which there is strong aerodynamic interaction among the various bodies.

The transonic problem has required the use of finer grids in order to generate acceptable solutions. It was observed that the artificial dissipation coefficient has to be decreased in order to improve the solution quality when compared to the value used in the supersonic test case. This should come as no surprise if one remembers that this artificial dissipation term was added to overcome the instability generated in strong shock waves. Furthermore, it was noted that the mesh refinement has a stronger influence on the result quality than the artificial dissipation coefficient. Mesh refinement makes the aerodynamic phenomena better defined and generates smoother contour lines in the solution.

The results here presented indicate that further grid refinement is still required in order to adequately resolve some features of such complex flows. However, the maximum run-time memory available, in the computational facilities the authors can access, has already been reached by the grids here used. Furthermore, it is clear to the authors that, in order to provide a simulation tool with a high degree of flexibility to handle different problems, a more automatic form of grid refinement is required. Hence, efforts have already been initiated in the implementation of an adaptive refinement capability in the context of the present flow solver. These efforts, however, are beyond the scope of the present paper. The main goals of the paper were to demonstrate the capability presently available for flow simulation over such complex configurations and to report on the first detailed results for the complete VLS configuration using an unstructured grid code. The authors believe that these goals have been fully accomplished.

Acknowledgments

The authors would like to acknowledge Fundação de Amparo à Pesquisa do Estado de São Paulo, FAPESP, which provides a graduate scholarship for the first author, and Conselho Nacional de Desenvolvimento Científico e Tecnológico, CNPq, which partially supported the project under the Integrated Project Research Grants No. 522413/96-0 and 501200/2003-7. The authors are also indebted to Centro Nacional de Supercomputação, CESUP/UFRGS, and to Núcleo de Atendimento em Computação de Alto Desempenho, NACAD-COPPE/UFRJ, which have provided the computational resources used for the present simulations.

Paper accepted March, 2004. Technical Editor: Aristeu da Silveira Neto

• Antunes, A.P., Basso, E., and Azevedo, J.L.F., 2001, "Chimera Simulations of Viscous Flows over a Complex Satellite Launcher Configuration," AIAA Paper No. 2001-2475. Proceedings of the 19th AIAA Applied Aerodynamics Conference, Anaheim, CA, (publication in CD-ROM format without page numbering).
• Argyris, J., Doltsinis, I.S., and Friz, H., 1989, "Space Shuttle: Exploration of Reentry Aerodynamics," Computer Methods in Applied Mechanics and Engineering, Vol. 73, pp. 1-51.
• Basso, E., Antunes, A.P., and Azevedo, J.L.F., 2000a, "Three Dimensional Flow Simulations Over a Complete Satellite Launcher with a Cluster Configuration," AIAA Paper No. 2000-4514, Proceedings of the 18th AIAA Applied Aerodynamics Conference, Denver, CO, pp. 805-813.
• Basso, E., Antunes, A.P., and Azevedo, J.L.F., 2000b, "A Realistic Application of Computational Fluid Dynamics for the VLS Project," Mini Symposium on High Performance Computation, Proceedings of the 21st Iberian Latin American Congress on Computational Methods in Engineering -CILAMCE 2000, Rio de Janeiro, RJ, (publication in CD-ROM format without page numbering).
• Christies, I., Griffths, D.F., Mitchell, A.R., and Zienkiewicz, O.C., 1976, "Finite Element Methods for Second Order Differential Equations with Significant First Derivatives," International Journal for Numerical Methods Engineering, Vol. 10, pp. 1389-1396.
• Chung, T.J., 1978, Finite Element Analysis Fluid Dynamic, McGraw-Hill, New York.
• Donea, J., 1984, "A Taylor-Galerkin Method for Convective Transport Problem," International Journal for Numerical Methods in Engineering, Vol. 20, pp. 101-119.
• Donea, J., and Giuliani, S., 1981, "A Simple Method to Generate High-Order Accurate Convection Operators for Explicit Schemes Based on Linear Finite Elements," International Journal of Numerical Methods in Fluids, Vol. 1, pp. 63-79.
• Heinrich, J.C., Huyakorn, P.S., Zienkiewicz, O.C., and Mitchell, A.R., 1977, "An Upwind Finite Element Scheme for Two-Dimensional Convective Transport Equations," International Journal for Numerical Methods in Engineering, Vol. 11, pp 131-143.
• Hughes, T.J.R., and Tezduyar, T.E., 1984, "Finite Element Methods for First-Order Hyperbolic Systems with Particular Emphasis on the Compressible Euler Equations," Computer Methods in Applied Mechanics and Engineering, Vol. 45, pp. 217-284.
• Lax, P.D., and Wendroff, B., 1960, "Systems of Conservation Laws," Communications on Pure and Applied Mathematics, Vol. 13, p. 217.
• Lax, P.D., and Wendroff, B., 1962, "On the Stability of Difference Schemes," Communications on Pure and Applied Mathematics, Vol. 15, p. 363.
• Lax, P.D., and Wendroff, B., 1964, "Difference Schemes for Hyperbolic Equations with High Order of Accuracy," Communications on Pure and Applied Mathematics, Vol. 17, p. 381.
• Löhner, R., Morgan, K., and Zienkiewicz, O.C., 1985, "An Adaptive Finite Element Procedure for Compressible High Speed Flows," Computer Methods in Applied Mechanics and Engineering, Vol. 51, pp. 441-465.
• Morgan, K., Peraire, J., Peiro, J., and Hassan, O., 1991, "The Computation of Three Dimensional Flows Using Unstructured Grids," Computher Methods in Applied Mechanics and Engineering, Vol. 87, pp. 335-352.
• Peraire, J., Peiro, J., Formaggia, L., Morgan, K., and Zienkiewicz, O.C., 1988, "Finite Element Euler Computations in Three Dimensions," International Journal for Numerical Methods in Engineering, Vol. 26, pp. 2135-2159.
• Scalabrin, L.C., 2002, "Numerical Simulation of Three-Dimensional Flows over Aerospace Configurations," M.S. Dissertation, Instituto Tecnológico da Aeronáutica, São José dos Campos, SP.
• Steger, J.L. and Warming, R.F., 1981, "Flux Vector Splitting of the Inviscid Gasdynamic Equations with Application to Finite Difference Methods," Journal of Computational Physics , Vol. 40, No. 2, pp. 263-293.
• Strauss, D., and Azevedo, J.L.F., 1999, "A Numerical Study of Turbulent Afterbody Flows Including a Propulsive Jet," AIAA Paper No. 99-3190, Proceedings of the 17th AIAA Applied Aerodynamics Conference, Norfolk, VA, pp. 654-664.
• Teixeira, P.R.F., 1996, "Numerical Simulation of Compressible Three-Dimensional Flows Using Finite Elements," M.S. Dissertation, Civil Engineering Department, Universidade Federal do Rio Grande do Sul, Porto Alegre, RS, (in Portuguese, original title is "Simulação Numérica de Escoamentos Tridimensionais de Fluidos Compressíveis Aplicando Elementos Finitos").
• Warming, R.F., Beam, R.M., and Hyett, B.J., 1975, "Diagonalization and Simultaneous Symmetrization of the Gas-Dynamic Matrices," Math. Comp., Vol. 29, No. 132, pp. 1037-1045.
• Zienkiewicz, O.C., and Morgan, K., 1983, Finite Elements and Approximation, Wiley, New York.

Publication Dates

• Publication in this collection
  12 Aug 2004
• Date of issue
  June 2004

History

• Accepted
  Mar 2004