Supercomputer Simulations of Fluid-Structure Interaction Problems Using an Immersed Boundary Method

The paper describes a supercomputer application in simulations of fluid-structure interaction problems. A compressible flow solver based on a high-accuracy scheme for unstructured hybrid meshes is considered. It combines an immersed boundary method with a dynamic mesh adaptation method in order to represent motion of solid objects in a turbulent flow. The use of immersed boundaries allows you to dynamically adapt the mesh resolution near moving solid surfaces without changing the mesh topology. Multilevel MPI + OpenMP parallelization of these components fits well with the architecture of modern cluster systems. The proposed implementation can engage thousands of CPU cores in one simulation efficiently. An example application is presented in which a high-speed turbulent flow around a cavity with a deflector is simulated.


Introduction
Mathematical modeling of the fluid-structure interaction (FSI) problems is crucial for fluid mechanics engineering, as well as for design of structures and reliability control. At the same time scale-resolving simulation of turbulent flows in dynamic geometrically complex configurations imposes high computing demands.
An immersed boundary condition (IBC) method is used in combination with a dynamic mesh adaptation method in order to reduce computing costs and increase efficiency of simulation as compared to traditional body-fitted approaches. Originally IBC methods were created for modeling of flows around arbitrary-shaped obstacles using structured cartesian grids. In the present work the Brinkman penalization method [1] is used within a high-accuracy compressible flow solver on unstructured hybrid meshes. The dynamic mesh adaptation technique based on variational principles [2] significantly reduces the required number of mesh nodes and improves the accuracy of the solid surface representation.
In contrast to a body-fitted approach, when solid surfaces are represented by exterior mesh faces, the use of IBC allows to easily track a solid surface by adding an external force field. Thanks to the ability of nodes to pass through solid objects, the mesh adaptation becomes much easier, affecting only coordinates of the nodes while the mesh topology remains constant. This makes recalculation of control volumes much cheaper. Furthermore, in this case there is no workload imbalance, which would otherwise appear when adding and removing mesh nodes.

Mathematical Model and Numerical method
The mathematical model is based on the compressible Navier-Stokes (NS) equations. At the fluid-solid interface the no-slip condition is imposed using the Brinkman penalization method [1], which does not require matching of the mesh nodes to the solid boundary. The special penalty functions are added as source terms into the NS system. These functions differ from zero only at the mesh nodes inside obstacle. In order to model the FSI problems, the NS system is coupled with the obstacle motion equation.
The NS system is descretized on an unstructured hybrid mesh using the high-accuracy edgebased scheme EBR [3]. This scheme provides high accuracy (up to 5-th order on translationally invariant meshes) at about the same computing cost as a basic low-order scheme.
The time integration is performed using an implicit second-order scheme based on the Newton linearization. The corresponding Jacobi system of linear equations is solved using a preconditioned Bi-CGSTAB solver [4].
The algorithm of the time integration step consists of several stages. Firstly, the motion equation is solved in order to update the coordinates and velocities of the mass centers of the solid bodies. Then the positions of solid surfaces are updated and tracked by the mesh adaptation method. Finally, new penalty functions are determined, and the penalized NS system is solved.

Parallel Implementation
The proposed method was implemented in the CFD code NOISEtte [5]. It has MPI+OpenMP parallelization for supercomputers made of multi-and manycore processors. A multilevel mesh partitioning is used for workload distribution among supercomputer nodes and among CPU cores inside nodes. A detailed description of the parallel algorithm can be found in [5].
The immersed boundary penalization method was implemented following the same parallelization approach. The calculation of penalty functions produces no significant load imbalance since it is not computing intensive.
The dynamic mesh adaptation follows MPI+OpenMP parallelization as well, but it introduces notably more complexity to the parallel algorithm. It involves a different parallel iterative solver [6] for a linear system of equations represented with another sparse matrix format.

Verification and Validation
Various test cases with available experimental and numerical data were considered in order to validate the accuracy of FSI modeling using Brinkman penalization method coupled with dynamic mesh adaptation.
The oscillations of 2D and 3D obstacles, such as a cylinder or a sphere, caused either by a given external force or by vortex-induced forces were investigated. The computational results are in a good agreement with the reference data.
A flow around a flapping foil represents more complex geometry and motion law. In this test case the NACA0012 foil oscillates with plunging and pitching mechanisms. The instantaneous vorticity magnitude for different phase angles is shown in Fig. 1. An example of the automatic mesh adaptation that tracks the surface of the foil is shown in Fig. 2. A good agreement between numerical and experimental results was observed in comparison of the measured and calculated mean thrust and power coefficients. Further details on comparison with experimental data and body-fitted numerical results can be found in [7].

Industry-Oriented Applications
The immersed boundary method can be efficiently combined with a body-fitted approach. Static objects that require accurate boundary layer resolution can be represented with a bodyfitted mesh. At the same time, moving objects or static objects that allow less accurate resolution can be represented with IBC. For instance, if we study a helicopter main rotor using a bodyfitted mesh and want to account for effect of the helicopter fuselage on the flow, then we can simulate the fuselage using IBC. Deployment of spoilers on a wing, opening doors of a gear bay, deployment of landing gears, release of a payload, etc., can be simulated with this combined approach.
The proposed method was applied in a simulation of a flow over a deflector mounted at the upstream side of a cavity (Fig. 3). The flow parameters correspond to a flight regime of a high-speed aerial transport vehicle. The aim of the study was to reduce acoustic and vibration load on the surface of the cavity by improving the shape of the deflector. Multiple configurations of the complex shape deflector were numerically studied using IBC. The use of IBC allowed us to easily change the shape with a given surface parametrization, and to use one spatial mesh for all the configurations. The computational results demonstrated that the deflector considerably reduces the mean pressure coefficient and the overall pressure fluctuation level at the cavity bottom and on its rear wall. Further details on this study can be found in [8]. IBC applicability tests were performed for static meshes with up to one billion elements using up to 10,240 CPU cores of Lomonosov supercomputer with parallel efficiency of doubling the number of cores above 90%.

Conclusions
The developed approach can be applied for supercomputer modeling of FSI problems on unstructured meshes. The immersed boundary method allows to efficiently handle moving obstacles. Parallel dynamic mesh adaptation provides it with necessary resolution near solid surfaces. All the components fit well with MPI+OpenMP parallelization and can be easily implemented in the existing parallel CFD codes.