This work presents a high-fidelity in-house Fluid Structure Interaction (FSI) solver devel- oped by combining discrete forcing Immersed Boundary Method (IBM) with a RK-4 based structural solver. Classification of the grid points as fluid, solid and IB points in the IBM framework and the solution of the pressure correction equations are the two most computationally expensive section in the numerical solver. These computational efforts can be significantly reduced by implementing OpenMP techniques. However, the successive over-relaxation (SOR) iterative method used in the serial code is not suitable for OpenMP parallelization as it shows data dependencies from previous iterations. Therefore, the Red-Black (RB) SOR is implemented to avoid the data dependencies

Majumdar, D

Sarkar, S

Shah, CL

English

Spiral - Imperial College Digital Repository

21st Annual CFD SymposiumAugust 8− 9, 2019,Bangalore.Performance Enhancement of an Immersed Boundary MethodBased FSI Solver Using OpenMPChhote Lal Shah1,∗, Dipanjan Majumdar1, and Sunetra Sarkar1Department of Aerospace Engineering, IIT Madras, Chennia-600036, India* Communicating Author, E-Mail: lal.chhote.shah@gmail.comAbstract. This work presents a high-fidelity in-house Fluid Structure Interaction (FSI) solver devel-oped by combining discrete forcing Immersed Boundary Method (IBM) with a RK-4 based structuralsolver. Classification of the grid points as fluid, solid and IB points in the IBM framework and thesolution of the pressure correction equations are the two most computationally expensive section in thenumerical solver. These computational efforts can be significantly reduced by implementing OpenMPtechniques. However, the successive over-relaxation (SOR) iterative method used in the serial code is notsuitable for OpenMP parallelization as it shows data dependencies from previous iterations. Therefore,the Red-Black (RB) SOR is implemented to avoid the data dependencies.Keywords: Fluid-structure interaction, Immersed Boundary Method, OpenMP, Red-Black SOR, parallelisa-tion1 IntroductionFluid–structure interaction (FSI) is of great relevance to a wide range of real life engineering applications.One of the main difficulties experienced by researchers in numerically simulating these FSI problems isassociated with the handling of complex geometry and time dependent movement of the solid boundary. Inthe conventional approaches, Arbitrary Langrangian Eulerian (ALE) method based flow solver is coupledwith a structural solver to study the combined FSI effects [1]. However, these ALE based approaches are notsuitable as significant effort is required to regenerate the mesh at every time-step as the solution is advanced.One alternative method to solve these FSI problems efficiently is to incorporate the Immersed BoundaryMethod (IBM) [2] as flow solver. Instead of direct application of physical boundary conditions at the gridpoints, in IBM, the presence of a solid boundary inside a fluid is incorporated by imposing constraints tothe model governing equations. Therefore, in IBM framework, even a complex geometry problem involvingmoving boundaries can easily be handled using simple structured grids without the requirement of re-meshingat every time step. The primary advantage of IBM over body-fitted ALE method is that the task of gridgeneration is greatly simplified which results in optimized CPU cost and ease of implementation.In the present work the authors have developed an in-house FSI solver in the C++ environment by com-bining discrete forcing IBM with RK−4 based structural solver. Although the mesh generation process istrivial in IBM framework, the process of classification of the grid points at every time step for a movingboundary problem is computationally intensive process. The solution of the pressure correction equations isalso one of the most expensive computational section. However, parallelization of all these computationallyexpensive processes using OpenMP (Open Multi-Processing) techniques can significantly reduces the com-putational time. In the Parallelized solver the red-black Successive Over-Relaxation (SOR) [3, 4] algorithmis implemented for solving the set of algebraic equations to avoid data dependencies from previous iterationwhich is otherwise present in the conventional SOR algorithm. Passive heaving of an elliptical foil with anactive pitching actuation in a free stream flow is simulated to test the efficacy of the parallelized code. Simu-lations are performed with the following iterative algorithms: (a) the traditional SOR, (b) RB-SOR and (c)RB-SOR with OpenMP on a Linux cluster with sixteen dual processor 2.60GHz 64 bit Intel Xeons, totalling32 processors. Details of the structural and flow solver and the parallelization process is discussed next inSection 2; Computational time taken by the three numerical approaches are compared in Section 3; Finally,Section 4 summarizes the present paper with the concluding remarks.122 Numerical methodology2.1 Structural solverAn elliptic cylinder is considered to be elastically mounted in the plunging degree of freedom and is allowedto oscillate in a free stream flow. The equation of motion of the center of the elliptic section is given in thenon-dimensional form below [5]:¨̄y +(cDMU∞)˙̄y +(2πfNcU∞)2ȳ =2πARρρbCL, (1)where,M is mass per unit span of the elliptic cylinder; ȳ denotes the vertical displacement non-dimensionalizedby the length of the major axis (D); c is the damping coefficient, U∞ is the free stream velocity, fN is thestructural natural frequency, AR is the ratio of minor to major axis of elliptic cylinder, ρ/ρb is the ratio offluid to solid density. CL is the coefficient of lift force acting at the center of the elliptic body and is suppliedby a flow solver. An active control is also implemented via a pitching motion as: θ = Aθsin (2πfDτ); where, θis the angular displacement, Aθ is the pitching amplitude, fD is the driving frequency and τ(= tU∞/D) is thenon-dimensional time. The structural governing equation is solved using the explicit 4th order Runge-Kutta(RK−4) method.2.2 Flow solverThe flow equations are solved broadly following the discrete forcing IBM proposed by Kim et al. [2]. Inthis framework of IBM, the boundary conditions are satisfied at the solid boundaries by adding a additionalforcing term to the momentum conservation equation. A mass source/sink term is also added to the continuityequation to ensure the mass conservation across the immersed boundary. The flow equations are solved on abackground Eulerian mesh and the location of the solid boundary is tracked by a set of Lagrangian markers.The Eulerian mesh nodes are needed to be classified as fluid points and solid points at every time step. Themomentum forcing and the mass source/sink term are applied throughout inside the entire solid domain toreduce the Spurious Force Oscillation (SFO). The flow around the oscillating body is primarily governed bythe incompressible Navier-Stokes Equations. The flow governing equations are given in the non-dimensionalform below:∂u∂τ+∇. (uu) = −∇p+ 1Re∇2u + f, (2)∇.u− q = 0, (3)where u = (u, v) denotes the flow velocity non-dimensionalised by the free stream velocity U∞, p is thepressure non-dimensionalised with ρU2∞ where ρ is the fluid density. Re = U∞D/ν is the Reynolds numberwhere ν is kinematic viscosity of the fluid. Here, f = (fx, fy) is the momentum forcing applied to enforce theno-slip boundary condition at the solid boundary immersed in the fluid.The flow governing Eqs. (2) and (3) are solved using Fractional Step Method (FSM) on a Cartesianstaggered grid arrangement. The convective term of the momentum equation is discretized using Adams-Basforth discretization whereas the diffusion term is discretized using the Crank-Nicolson technique. Thevelocity is corrected using a pseudo pressure correction term so that the continuity is satisfied at every timestep. The discretised form of Eqs. (2) and (3) are as follows:ûk − uk−1∆τ= −∇pk−1 − 32H(uk−1)+12H(uk−2)+12ReL(uk−1)+12ReL(uk−2)+ fk, (4)∇2φk = 1∆τ(∇.ûk − q̂k), (5)uk = ûk −∆τ∇φk, (6)3pk = pk−1 + φk − ∆τRe∇2φk, (7)where L(∗) is the discretised form of the Laplace operator and H(∗) denotes the discretised form of theconvective term, ûk is the intermediate velocity at kth time step, φ is the pseudo pressure correction termand∆τ is the time increment. All the spatial derivatives are evaluated using the second-order central differencescheme.In the present FSI framework, the flow equations are solved to get the flow-field around the body and theaerodynamic loads acting on the body are evaluated at every time step. Then this lift force computed by theflow solver is supplied to the structural equation to compute the position of the body in the next time step.Thus, at every time step, the flow solver and structural solver exchange informations in a staggered mannerresulting in a weak coupling FSI solver.2.3 ParallelizationOpenMP. The most computationally costly sections in the sequential code are: (i) classifications of thegrid points as fluid, solid and IB points in the IBM framework and (ii) solution of the pressure correctionequations. OpenMP is implemented in these sections along with the other possible places to enhance theoverall performance of the solver. OpenMP is an Application Programming Interface (API) that serves aconvenient, scalable model to the developers of shared memory parallel applications [6]. It uses explicit directmulti-threaded and shared memory parallelism. OpenMP uses the fork-join model of parallel execution aspresented in Fig. 1. OpenMP employs multiple processors to execute the same code and thus acceleratesthe numerical process, hence decreases computational time. It is necessary to check data dependencies, dataconflicts, race conditions and deadlock situations before implementing OpenMP as it may cause erroneousresults, if not taken care properly. OpenMP is called in the loop executable using “#pragma omp parallel”.OpenMP is also implemented for solving discretized momentum equations and for reading and writing severaldata files parallelly.Fig. 1. Fork-Join model for OpenMP parallelisationRed-Black SOR. The traditional SOR algorithm used in the sequential code suffers from data dependencyon previous iteration and therefore cannot be parallelized using OpenMP. The RB-SOR method [3,4], on theother hand, eliminates this issue and takes advantage of OpenMP parallelization. The RB ordering techniqueand its implementation are schematically shown in Fig. 2. The red and black nodes are located alternativelyas shown in the figure. In RB-SOR algorithm, the values at the black nodes are updated using the values atthe red nodes and the values at the red nodes are updated using the values at the black nodes. Therefore, itcan be partitioned into independent tasks that can be performed by multiple threads simultaneously.4(a) (b)Fig. 2. Red-Black ordering technique and implementation of the SOR algorithm: (a) updating the values of the blacknodes using red nodes and (b) updating the values of the red nodes using black nodes.3 Results and discussion3.1 Simulation ResultsThe flow domain and the mesh grid used in the present study are shown schematically in Fig. 3. At the initialtime, the centre of the elliptic foil lies at the origin. A Dirichlet condition is applied at the left-hand sideboundary as inflow, and the slip boundary conditions are applied at upper and lower boundaries. A Neumann(a) (b)Fig. 3. (a) Flow domain and corresponding boundary conditions and (b) Cartesian grid.type boundary condition is employed at the velocity outlet. The boundary conditions are shown in Fig. 3(a).5The Non-uniform cartesian grid is considered in this study as shown in Fig. 3(b). For the pressure correctionequation, Neumann type boundary condition is applied at all the domain boundaries.To test the efficacy of the present FSI solver, passive heaving with active pitching actuation of an ellipticcylinder is simulated with parameter values as: cD/MU∞ = 0.503, fNc/U∞ = 0.2, AR = 1/6 and ρ/ρb = 0.2and the results are compared with that of Griffth et al. [5]. The simulations are performed at a low Reynoldsnumber ofRe = 200. The lengths of the computational domain along the stream-wise and transverse directionsare 32.5D and 25D, respectively. The flow domain is discretized using a structured Cartesian mesh with718 × 1017 grid points along the horizontal and vertical directions, respectively. A constant time step of10−4 is considered to perform the simulations. The computational domain, grid spacing and time step sizeare chosen after performing suitable convergence test. Simulation results of SOR, RB-SOR and RB-SORwith OpenMP solvers match very closely with each other. Here, we present the comparison of instantaneousFig. 4. Comparison of flow-field: (a) Griffith et al. [5] and (b) present work using OpenMP parallelisation(a) (b)470 472 474 476 478 480-2.5-2-1.5-1-0.500.511.522.5470 472 474 476 478 480-2.5-2-1.5-1-0.500.511.522.5Fig. 5. Comparison of (a) CL, and (b) vertical displacement y obtained from the present work using OpenMP paral-lelisation with that of Griffith et al. [5].6vorticity contour, evolution lift coefficient (CL) and vertical displacement (y) from the present work usingOpenMP parallelisation with that of Griffith et al. [5] in Fig. 4 and 5, respectively. Very good agreement isobserved between the results from the present OpenMP FSI solver and that of Griffith et al. [5].3.2 Performance ResultsThe performance enhancement of the numerical solvers is evaluated in terms of computational time taken bythe solvers to perform 1000 time step iterations i.e. τ = 0 to τ = 0.1. The solver speed (S) is defined asS =TSOR − TxTSOR× 100%, (8)where, Tx indicates the real clock time taken by RB-SOR or RB-SOR with OpenMP. The number of threadsconsidered for the execution of code in parallel are sixteen for all the computations. The performance ofthe iterative solvers on a Linux cluster with sixteen dual processor 2.60GHz 64 bit Intel Xeons, totaling 32processors is compared in Table 1. The RB-SOR with OpenMP shows significant speedup as compared tosequential FSI solver.Table 1. Time taken to compute 1000 time steps (non-dimensional time (τ))Solver Time taken (T/sec.) Speedup (S)SOR 12920 −Red-Black SOR 11780 8.77%Red-Black SOR withOpenMP3740 70.61%4 ConclusionIn this work, a high-fidelity in-house Fluid Structure Interaction (FSI) solver developed by combining dis-crete forcing Immersed Boundary Method with a RK-4 based structural solver is presented. The RB-SORmethod is successfully implemented to get rid of data dependencies from previous iterations and the solveris parallelize using OpenMP technique. The computational time taken to execute 1000 time-steps have beencompared between SOR, RB-SOR and RB-SOR with OpenMP parallelization. The RB-SOR with OpenMPparallelization is seen to run significantly faster and has shown the notable speed up as compared to serialcode.References1. H. T. Ahn and Y. Kallinderis, “Strongly coupled flow/structure interactions with a geometrically conservative alescheme on general hybrid meshes,” Journal of Computational Physics, vol. 219, no. 2, pp. 671–696, 2006.2. J. Kim, D. Kim, and H. Choi, “An immersed-boundary finite-volume method for simulations of flow in complexgeometries,” Journal of Computational Physics, vol. 171, no. 1, pp. 132–150, 2001.3. C. Zhang, H. Lan, Y. Ye, and B. D. Estrade, “Parallel sor iterative algorithms and performance evaluation on alinux cluster,” Naval Research Lab Stennis Space Center MS Oceanography Div, Tech. Rep., 2005.4. I. Yavneh, “On red-black sor smoothing in multigrid,” SIAM Journal on Scientific Computing, vol. 17, no. 1, pp.180–192, 1996.5. M. D. Griffith, D. L. Jacono, J. Sheridan, and J. S. Leontini, “Passive heaving of elliptical cylinders with activepitching–from cylinders towards flapping foils,” Journal of Fluids and Structures, vol. 67, pp. 124–141, 2016.6. B. Chapman, G. Jost, and R. Van Der Pas, Using OpenMP: portable shared memory parallel programming. MITpress, 2008, vol. 10.

Performance enhancement of an immersed boundary method based FSI solver using OpenMP

http://spiral.imperial.ac.uk/bitstream/10044/1/101664/6/P27-Performance%20Enhancement%20of%20an%20Immersed%20Boundary%20Method.pdf

Performance enhancement of an immersed boundary method based FSI solver using OpenMP

Abstract

Similar works

Full text

Available Versions

Spiral - Imperial College Digital Repository