II. PROGRESS IN GRAND CHALLENGE APPLICATIONS

II.1. Coupled Fields

This component considers the parallel simulation of coupled field problems in structural and fluid mechanics. Three Grand Challenge focus application problems were listed in the original proposal:

  1. Aeroelasticity of a complete aircraft
  2. Distributed control of flexible structures
  3. Coupling of electrothermomagnetic and superconducting phase-change behavior

The project focused initially on the aeroelastic problem, which had been studied in our Center since 1990. Work in the two other problems started on April 1993. The third problem was terminated in 1995 and resources focused on the other two application areas and their eventual combination. The original goal of the aeroelastic project: aeroelastic analysis of a complete aircraft in steady motion, was first attained in late 1994. Work in that area is shifting of the simulation of arbitrary maneuvers of a complete aircraft and the development of local turbulence models. The original goal of the distributed-control application was attained in early 1995, at which point it was decided to link up the first two focus problem so as to target the first coupled simulation of aeroelasticity and vibration control of a complete flexible aircraft.

Past and projected timelines for the various tasks are as follows:



A more detailed account of our overall approach to coupled fields simulation and our progress in the focus application areas follows.


a. General Approach

Our approach to the parallel simulation of coupled problems is based on a two-level problem decomposition:

  1. Field Decomposition divides the overall problem into subproblems generically called "partitions". These are usually physically distinct entities, such as the structure and the fluid in aeroelastic problems. Sometimes it is useful to designate a non-physical computational device as a separate partition; for instance the fictitious spring network that models a dynamic fluid mesh.
  2. Domain Decomposition divides the discrete spatial model of each partition into subdomains that are assigned to individual processors. The decomposition is carried out through mesh decomposition techniques.

Partitions are treated by separate programs. For example in the aeroelasticity application there are three partitions: structure, fluid and dynamic fluid mesh; which are processed by three individual programs. The discretization and time integration techniques for each program are tailored to the subproblem. Thus in the aeroelastic problem, it is natural to treat the structure by finite element methods and implicit time stepping. However the exterior fluid is better handled through finite volume discretization and a choice of explicit or implicit time stepping, depending on the type and goal of the simulation.

A distributed-memory, message-passing paradigm was adopted for all applications following early experiments on hybrid machines such as the KSR-1, as discussed in II.2.b. This computational model is used also in shared-memory machines such as Crays and SGIs. Programs treating subproblems are mapped to available processors as per expected load balance. For example, if 128 processors are available to an aerolastic simulation, the fluid, dynamic mesh and structure could be allocated 80, 28 and 20 processors, respectively. Each of those processors runs its program copy on the different data structures predefined by the domain decomposition.

The time stepping is handled by staggered integration techniques originally developed in the 1970s by Park and Felippa on sequential machines. Programs run with their own time steps, exchanging information at synchronization points. For example, if in an aeroelastic smulation the structure advances implicitly whereas the fluid is advanced explicitly (a common situation), the timestep for the fluid might be 100 or even 1000 times smaller. Such techniques are collectivelly called subcycling methods.

Information transfer between programs is of two types. Intrapartition transfers (e.g. structure to structure or fluid to fluid) are handled by standard message passing techniques. Interpartition transfers (e.g. structure to fluid or control to structure) are also handled by message passing but may require a more elaborate mesh interpolation procedure as further discussed in II.2.b.


b. Aeroelasticity of Complete Aircraft

This GC focus problem addresses the parallel simulation of the dynamic interaction between a complete flexible aircraft and the surrounding fluid traversed in the transonic through hypersonic regimes. This project has been directed by Professor C. Farhat, with present participation by Cai, Felippa, Mandel, McBryan, Djellouli, Lesoinne, Koobus, Manzouri, Sarkis, Chen, Degand, McArthur and Pierson, and past participation by Crivelli, Laney, Lanteri, Maman, Manzouri, Rixen, Roux, Stern and Haugen.

At the physical level the problem decomposes naturally into two partitions: structure and fluid, which are treated by Lagrangian and Eulerian descriptions, respectively. There is also a third computational partition formed by the fluid mesh near the aircraft that is displaced by vibratory structural motions. The three fields are handled by different codes, identified as the "fluid", "structure" and "dynamic mesh" in the sequel. The structure code is based on finite element discretization in space. The fluid code is based on a combination of finite volume and finite element discretization techniques on unstructured meshes. The dynamic mesh code uses a network discretization consisting of fictitious springs and masses placed along the edges and on nodes of the fluid mesh elements.

For parallel processing each of the preceding components is treated by mesh and domain decomposition techniques that separate each discrete model into subdomains. Subdomain data is assigned to one processor (CPU), which contains a copy of the appropriate program. The simulation is advanced in time through staggered solution procedures. Each partition runs with its own time step, exchanging information at synchronization stations. Because the fluid and structure meshes are generally nonconforming, inter-field information transfer must be done through interpolation techniques on the fluid-structure interface.

The project started from the foundation provided by two-dimensional aeroelastic codes developed on the Cray XMP and Intel 860 Hypercubes during 1990-1992. Beginning in 1992, these were ported to third-generation parallel machines: the KSR-1 at CU and Cornell, the Paragon XP/S at NOAA, and the CM-5 at NCSA. Much of the initial development work was done on the KSR-1. Two implementations, one based on a Cray-like shared memory model and the other on a message-passing (MP) distributed memory model, were tested on the KSR. The latter was found to be significantly more efficient. On the basis of these results, it was decided in 1993 that subsequent developments for heterogeneous operation were to follow the MP model. Early implementations followed either PVM (if available) as well as vendor-supplied MP libraries such as TCGMSG. Starting in 1995, we decided to stardardize on MPI on all major platforms.

The first three-dimensional simulations, carried out in 1993, dealt with portions of the target problem. As the project developed and more complex simulations were undertaken, there emerged a set of issues interweaving physics modeling, numerical analysis, and computer sciences. Major developments in the 1993-95 period are noted below.

Scalable Iterative Structure Solvers. A new solution approach called Finite Element Tearing and Interconnecting (FETI) was developed beginning in 1990, and progressively refined during this project. This approach solves the domain interface force equations iteratively through a projected conjugate gradient technique, while the internal subdomain equations are processed by direct solvers. Solvers based on this approach were shown to scale into the MPP regime (that is, with 64 or more processors) if appropriate subdomain aspect ratios are maintained and if plate bending and shell elements are not used in the finite element discretization. A 1995 extension of the original method, called FETI2, removed the second scalability constraint. Using FETI2, a plate problem with almost one million degrees of freedom was solved in 45 iterations and 105 seconds CPU on a 64-processor IBM SP2.

Memory FETI Solver. A strategy to "reuse" previouly found conjugate directions was demonstrated to significantly reduce the solution time for multiple correlated right hand sides as they occur, for example, in the linear dynamic analysis of the structure model.

Fluid-Structure Mesh Information Transfer. A parallel MATCHER program that links separately constructed, unstructured meshes for the fluid and structure partitions was developed and used in the 3D computations.

Heterogeneous Processing Experiments. The first heterogeneous fluid-structure computations were initially demonstrated on an iSPC Hypercube. In these simulations the fluid and structure partitions run concurrently on different computers. Experiments were also performed on the KSR-1, Paragon XP/S, and IBM SP2 machines. Binary incompatibility between different platforms was found to be a more significant communication bottleneck than network bandwidth.

Parallel Staggered Procedures. A new class of partly-augmented, subcycled staggered solution procedures aimed to reduce structure-fluid communication overhead while maintaining satisfactory stability and accuracy was developed and theoretically studied in 1D model problems. This study considered the effect of subcycling and fluid mesh motion. The best among a set of alternative implementations was then implemented and tested on realistic 2D and 3D problems.

Geometric Conservation Laws. It was discovered that the presence of a moving fluid mesh driven by structure deformations can induce serious violations in the conservation of various geometric quantities unless special care is taken in the time integration procedure to minimize such deviations.

Fluid Time Integration. For expedience, standard explicit time-stepping algorithms such as MUSCL were used in the initial phases of this effort. Beginning in 1995 a more sophisticated family of implicit time integrators were tested in 2D problems and then implemented in 3D simulations as alternative to the explicit schemes. Cai and Sarkis of Computer Sciences were instrumental in the development of the necessary theory and algorithmic framework.

The combination of these interdisciplinary efforts was highlighted in late 1994 with the first aeroelastic simulation of a complete aircraft: a Falcon jet driven through a smooth parabolic trajectory. Since late 1995, the focus on the aeroelastic project has been redirected toward the simulation of arbitrary manuevers, including emergency ones such as those required to avoid colisions. Such simulations are of high practical interest because they cannot be pretested in wind tunnel models. Achieving this ambitious target requires the integration of additional modeling and computational capabilities:

  1. Large deflection effects in the structure model.
  2. Turbulence model in dynamic unstructured fluid meshes
  3. Interaction of propulsion forces and servo control.

Work during Period 4 has focused on the first two items, and is further discussed under item e below.


c. Distributed Control Structure Interaction

The control-structure interaction (CSI) focus problem considers the parallel simulation of a control system with a flexible structure. The project is supervised by K. C. Park with present participation by Felippa, Manzouri, Farhat, Lesoinne and Chen, and past participation by Menon, Hwang and Stern.

The project was initiated in 1993. Emphasis was initially placed on the control of vibrations of space structures. The key development in the period 1993-1994 was a method was developed to design and place a large number (hundreds or thousands) of actuators and sensors on complex structures. The method analyzes the structure under representative loading conditions and identifies the spatial energy concentration pattern and the regions of high-rate energy transmissions. Rate feedback controllers (or active dampers) are placed on the high energy concentration regions and inertial controllers (or force actuators) are placed on the high-energy transmission regions, respectively. The method, developed and refined by Menon and Park, was implemented and tested on the Paragon and KSR parallel machines.

The merging of this activity with the aeroelasticity problem is discussed in item f below.


d. Coupling of Electrothermomagnetics and Phase Change

This focus problem considered the interaction between mechanical, thermal and electromagnetic fields to develop and test a model capable of representing co-existing normal and superconductivity phases of a conductor. The project was carried out by Felippa and Schuler.

This project was initiated in 1993. The problem was found to be extremely nonlinear, with interaction of quantum-mechanics and macromechanic effects varying by up to 30 orders of magnitude. Much of the initial development was algorithmic, with focus on overcoming such difficulties. The first simulation of a transition state was achieved on the KSR-1 in 1994. These involved the solution of the coupled thermomechanical and Ginsburg-Landau equations on a 128x128 axisymmetric mesh.

In early 1995 a decision was made to conclude this focus application to concentrate resources on the other two focus problems (aeroelasticity and control-structure interaction) and their eventual combination.


e. Aeroelasticity of a Maneuvering Aircraft

As mentioned under item b, achieving success in maneuvering aeroelastic simulations requires the consideration of nonlinear structural finite elements, turbulence models for dynamic meshes, and interaction of propulsion forces with trajectory control.

The use of nonlinear finite elements is required because displacements and rotations in the Lagrangian description of the aircraft structure cannot longer be assumed small. Modeling such motions calls for the incorporation of geometrically nonlinear plate, shell and beam elements in the structure code. Preliminary work in this area was conducted by Felippa and Manzouri from October 1995 through March 1996. It was decided to chose a finite element corotational description of large motions since that simplifies the reuse of existing linear elements under large deflection "wrapper" routines. The ability of such a description to simulate dynamics of trajectory-specified large rotations (for example an airplane turning 90 deg in a few seconds) was tested in simple test problems. The next stage, scheduled to start by January 1997, is to incorporate corotational elements in a rewritten version of the structure code presently under development.

Another key component in maneuver simulation is the computation of turbulent flows in 3D moving fluid meshes, a topic not previously addressed by CFD researchers. Koobus and Farhat started by developing a methodology for simulating turbulent flows on 2D unstructured dynamic grids. Turbulence effects are modeled through a wall-law k-epsilon relation. By integrating the fluid and turbulence equations up to a small distance from the wall, the use of a highly refined mesh near the wall is avoided. Hence the model results in less expensive numerical simulations. These savings can be essential for unsteady 3D calculations with fluid meshes that can contain millions of tetrahedra.

The wall-laws of Reichardt are used, in combination with the generalized boundary conditions for k-epsilon proposed by Jaeger and Dhatt. This model allows the use of computational domains which extend not too far from the wall. It improves the prediction of unsteady separated flows and vortex shedding, which may appear in rapid maneuvers. For the numerical solutions we have used we use a Finite Element / Finite Volume formulation for the space discretisation of the fluid and turbulence equations. The numerical method is implicit and second order accurate both in space and time and uses a Defect Correction approach. Dynamic meshes are handled by a specific treatment of the convective and diffusive fluxes.

A solver for simulating low Mach (< 0.01) flow (laminar and turbulent) has also been developed for fixed and moving grids. This capability is important for aircraft maneuvers at very low speed and dynamic stall analysis. This methodology has been applied to the simulation of turbulent flow problems around oscillating airfoils at different Reynolds, Mach and angle of attack. We have also considered the problem of the vortex shedding with the simulation of turbulent flows over a square and a H-profile bridge. The results are in good agreement with the experiments when available. The next stage calls for extending this methodology to the simulation of turbulent 3D flow problems on unstructured moving grids.

The coupling of propulsion forces and control surfaces with trajectory control has not been yet addressed.


f. Aeroservoelasticity of Complete Aircraft

In the summer of 1995 it was decided to gradually merge control-structure interaction and aeroelasticity applications, having as ultimate goal the simulation of the fluid and structure interacting with a control system for trajectory, vibration and flutter control. This three-discipline coupling, known as aeroservoelasticity, may be said to represent a "Ultra Grand Challenge" since two-discipline couplings are still viewed as Grand Challenges.

The initial merge step attempted to interfacing the "dry structure-control" interaction code with the parallel aeroelasticity code for simulating the vibration and flutter control of a vibrating airfoil. This coupling run into numerical stability difficulties that precluded the use of a reasonable common time step for the CSI coupling. To circumvent that difficulty, a new implicit-implicit staggered procedure was theoretically derived, implemented and tested. This new algorithm was interfaced to the parallel aeroelastic code over the period April through October 1996.

After many adjustments, the first proof-of-concept simulation of the vibration control of a flat plate submerged in a subsonic flow was successfully carried out in November 1996 by Manzouri on a Power Onyx. This is probably the first time in which three totally separate programs: fluid with dynamic mesh, structures and control, have run concurrently in a message-passing environment. Attempts to control the flutter that develops at supersonic speeds were not successful, however, since the initial controller design did not take into account the significant changes in the structural eigenpectrum as the flutter speed is appraoched. Overcoming that limitation will require the use of adaptive control laws.