O. Hassan, K. Morgan and N. P. Weatherill School of Engineering, University of Wales

1. A Method for Time Accurate Turbulent Compressible Fluid Flow Simulation with Moving Boundary Components Employing Local Remeshing O. Hassan, K. Morganand N. P. Weatherill School of Engineering, University of Wales Swansea, United Kingdom Workshop on Mesh Refinement of Unsteady Flows 7 December 2005 - Oxford University

2. Outline Problem of Interest and Adopted Approach Solution Algorithm Unstructured Mesh Generation Techniques Mesh Adaptation Techniques Mesh Adaptation for Unsteady Flow Problem with Moving Boundary Parallel Implementation Conclusion

3. Industrial End-User Configurations DASA - F16 Airbus - A340 Dassault - Falcon Airbus - A380

4. Adopted Approach • The governing equations are the Navier-Stokes equations • The application requires the ability to model complex geometries • The simulation of turbulence affords a real challenge • For many applications, the Euler equations are appropriate • Computational requirements for realistic geometries can be expected to be large Unstructured grid technology provides the required flexibility for these (and other) applications

5. Governing Equations • The Favre Averaged Navier Stokes Equations Where and • Turbulent is modelled by adding the one equation model of Spalart and Allmaras

6. Solution Algorithm • Edge Based Data Structure: Typical interior node I • The ALE term for an interior node I is: Where: should lead to a numerical ALE flux that satisfies GCL • Resulting Equation

7. Solution Algorithm • The turbulent viscosity equation is discretised in a similar fashion • Stabilisation achieved by replacing the actual flux function by JST flux function • Discontinuity capturing achieved by the addition of a switched artificial diffusion • For steady state Runge-Kutta relaxation and local timestepping is utilised • Convergence acceleration is achieved by using the Full Approximation Storage (FAS) Multigrid scheme • Coarse grids are achieved by agglomeration • Volume weighted operator is used for restriction • Injection is used for prolongation • Parallelimplementation allows agglomeration across partitions

8. AIAA test case: Drag Prediction Workshop 2001 M = 0.75 Re = 3 x 106 1.6 million points 35 viscous layers 5 grids levels Steady Turbulent Flow Pressure distribution Cl vs a Cl vs Cd

9. Time Discretisation • For unsteady problems, the second order approximation is adopted • An implicit formulation is employed and this removes the stability constraints associated with explicit schemes • At each time step, the equation is solved by explicit relaxation with multi-grid acceleration • This approach can be thought of as converging the set of steady state equations with the addition of the time source for every physical timestep • No significant memory penalties compared to explicit procedures

10. Mesh Generation • The surface is defined as a set of: • Surface Components: bi-cubic patches, NURBS • Curve components: cubic splines, NURBS • Mesh control: • Background Mesh • Point, Line, circular and planar sources • Curvature Controlled e: the gap between the element and the surface k1, k2 : the two principle curvatures

11. Mesh Generation • Surfacemesh generation: Advancing Front • Volumemesh generation: Delaunay Triangulation with automatic point insertion (Requires 100Mb/106 elements) • Boundary layer generation: Hybrid meshes by the Advancing Layer method

12. Improved Volume Mesh Quality • 3D Edge Swap • 3D Edge Collapse • 3D Nodal Smoothing • 3D Local Re-generation

13. Improved Surface Mesh Quality • Super surfaces eliminate small and distorted patches generated by the CAD systems • Merge neighboring surfaces based on continuity of the normal • Patch Independent Remeshing • Starting from any triangulation, re-triangulate using edge splitting, edge collapse and edge swapping.

14. Mesh Adaptation • Error Analysis • Error indicator based on posterior error analysis are also possible, but not very practical for unsteady flow. • A Error indicator based upon interpolation theory and accounts for directionality is employed. • Assuming exact nodal values, estimate the local error for each elements as: • Equidistribution of the error results in a mesh spacing d for the new mesh: • In 2D/3D: Apply the 1D criterion separately to each principal direction of the Hessian Matrix

15. Mesh Adaptation • Mesh Enrichment • Advantages: • Simple and quick to implement • Trivial interpolation • Disadvantages: • Multiple refinement results in large meshes • De-refinement require excessive storage • Incorporating stretching results is distorted elements • Not suitable for unsteady flow with moving components

16. Mesh Adaptation • Mesh Enrichment • Special care is required in 3D to ensure compatibility of adjacent elements • Special care is also needed to ensure the validity of the mesh after projecting the added points onto the surface Geometry Initial Grid Adapted Grid

17. Mesh Adaptation • Mesh Movement • Replace the sides of the mesh by spring • Spring stiffness depends on the flow properties • Where • Move the nodes until nodal equilibrium • Solve by iteration

18. Mesh Adaptation • Mesh Movement • Advantages: • Simple and quick to implement • Can handle moving components • Disadvantages: • Expensive interpolation • Initial mesh may lack the required resolution to resolve all the flow features • Hard to control the quality of the moved mesh • Coupling of mesh movement and mesh enrichment can over come most of the restrictions.

19. Mesh Adaptation • Adaptive Remeshing • Use the current mesh as a background mesh • At each node compute the mesh parameters using the equidistribution criterion • Use the geometry definition to regenerate the mesh

20. Mesh Adaptation • Adaptive Remeshing Using a Background Mesh • Advantages: • Simple and quick to implement • Can handle moving components • High quality meshes • Disadvantages: • Expensive to regenerate the complete mesh

21. No coarsening beyond • the initial mesh • Multiple refinement can • generate generate distorted • elements • Memory intensive in 3D Useful for steady state • Reliable moving methods • are expensive • Not easy to guarantee • valid mesh • May not have enough • initial points Efficient for unsteady flow with small moving boundaries Essential for unsteady flow with large moving boundaries Time consuming in 3D Adaptation Mesh Enrichemnet Mesh Movement Remeshing For unsteady flow with moving boundary components it is essential to develop a scheme which utilise the advantages of the various methods

22. Adaptation for Unsteady Flow • Generate initial mesh • Calculate initial quality measure, qo and spacing, do Loop over timesteps • Update coordinates of moving nodes • Apply deforming mesh algorithm • Compute new quality measure, qn • Compute new spacing required based on the equidistribution criterion, dn Mark element to be deleted Form holes from marked elements T F Remesh hole Interpolate solution Compute qo, do

23. Interaction between a strong shock and an object of complex shape 1228 elements 650 points 3780 elements 1946 points 8554 elements 4348 points 11726 elements 5956 points Density contours

24. Confined Blast Wave From Rupturing Cylindrical Pressure Vessel Experimental Apparatus 1499 < Number of elements < 21022 Comparison of pressure history at various transducers

25. Surface Adaptation The Geometry definition is utilised for the regeneration of the surface portion of the hole

26. Surface Adaptation

27. Unsteady Simulation B60 Configuration M = 0.801 a0 = 2.780 am = 10 Zm = 1m 745198 Elements 135760 Points

28. Unsteady Inviscid Flow Store Separation Simulation ainit = zero degrees M= 0.5 Degree Container motion computed 2.7 million tetrahedra 15 time steps 50 multigrid cycle per time step Geometry for a complete F16 Configuration 8h Wall clock time Solver: 16 R14000 CPUs Preprocessing and adaption : 1 CPU

29. Store Separation Simulation Surface Pressure Distribution

30. Update coordinates of moving nodes and the boundary layer nodes Apply deforming mesh algorithm Calculate new dihedral angles (Din) F Adaptation for Unsteady Turbulent Flow Loop over physical time-steps • Generate initial mesh • Calculate minimum dihedral angle (Dio) • Store the layer number for all nodes in the boundary layer Mark elements to be deleted Form holes from marked elements Din<0 • Check elements intersection • Check if any boundary layer node can grow further • Remesh hole T Interpolate solution & Compute Dio

31. Unsteady Turbulent Flow • NACA64A010 Aerofoil: • Prescribed sinusoidal oscillation • Amplitude 1.01 degrees • ainit = zero degrees • St = 15.567 • 3D stacked hybrid mesh • 173720 nodes • 300 multigrid cycle/time step • 32time steps per cycle 16 R14000 processors One movement cycle 5 h clock time Estimated speed up 175

32. Unsteady Turbulent Flow k = 0. Lift Polar k = 90 k = 270

33. Unsteady Turbulent Flow Shuttle Booster Separation Simulation ainit = zero degrees M= 0.85 Degree Re = 3 * 10 6 Prescribed Shuttle movement Initial mesh: 2.9 million elements Final mesh: 3.3 million elements 20 time steps 300 multigrid cycle per time step 36h Wall clock time Solver: 16 R14000 CPUs Preprocessing and adaption : 1 CPU

34. Unsteady Turbulent Flow Shuttle Booster Separation Simulation Cut through the volume mesh Meshes on the symmetry plane

35. Unsteady Turbulent Flow Shuttle Booster Separation Simulation Meshes of the symmetry plane after remeshing Cut through the volume mesh after remeshing

36. Unsteady Turbulent Flow Shuttle Booster Separation Simulation

37. Parallel Implementation • Elements are select to be remeshed in each domain separately • Selection Based on Deviation from Prescribed Spacing • Selection Based on Element Quality • Selection Based on Intersection Tests • To determine intersection of elements due to moving geometries, one ovelapping ghost layer of elements is used. • If intersection with the ghost element has occurred, the search will also take place in to the domain which own the ghost cell.

38. Domain 2 Domain 1 3 1 2 1 2 1 1 Interfaces Domain 3 Domain 4 2 Parallel Implementation • Constraint Repartitioning is employed to ensure that each region to be remeshed will • be contained completely on one process

39. Unsteady Inviscid Flow Geometry for a complete F18 Configuration Store Separation Simulation ainit = 0.46 degrees M= 0.96 Container motion computed 2.1-2.3 million Nodes 12.1-13.4 million Elements 40 Physical timesteps with sub-cycling 10.4 h on 24 Processors

40. Unsteady Inviscid Flow Store Separation Simulation

41. Unsteady Inviscid Flow Store Separation Simulation • CFD Solution 40% • Motion Application 57% • Mesh Deformation 10.3% • Volume mesh Analysis 3.7% • Volume remeshing 37% • Re-partitioning 6% • I/O 3%

42. Conclusions • A hybrid unstructured finite volume method for aerodynamic flows has been presented • Turbulent flows are treated via the one equation Spalart and Allmaras model • Computational performance is enhanced by the use of multigrid acceleration and parallelisation • Transient moving boundary flows are treated by an ALE approach • Mesh movement and adaptive remeshing have been employed to handle the deformation due to the moving components • Adaptive remeshing was extended to meshes with stretched elements in the boundary layers • Parallel implementation of the adapted remeshing has been completed • A number of challenging problems have been simulated and the agreement with available experimental observations is good