On aerodynamic force computation in fluid–structure interaction problems — Comparison of different approaches

https://doi.org/10.1016/j.cam.2023.115208Get rights and content

Abstract

This paper is interested in the numerical approximation of a two-dimensional fluid–structure interaction problem. A special attention is paid to the numerical realization of aerodynamic force (AF) calculations. Three different AF computations within finite element method framework are described and compared. First, the extrapolation of fluid stress tensor components from the interior of the fluid domain is the most straightforward method; however, this approach is of lower theoretical accuracy due to the explicit computation of velocity gradient. Second, the local reconstruction method applied on the velocity gradient is based idea of the least square fit of the gradients evaluated at finite element centroids from local patch. Third, the weak reformulation technique using the weak formulation of the fluid flow problem leads to the expression of AF in an integral form. As it avoids direct velocity gradient computation it is of higher accuracy compared to previous methods.

The experimental convergence of all mentioned methods for the drag and lift coefficients is tested on the benchmark of a static cylinder in a cross-flow and further extended to the case with prescribed motion of the cylinder. Further, two FSI simulations compare the transferred energy from the fluid to the structure for the case of structure vibrations submerged in the viscous fluid and also for the flow-induced structural vibrations. The critical inlet velocity of flutter instability is determined and it is shown that for a given numerical setting the choice of AF computation can lead to the different FSI behavior.

Introduction

The fluid–structure interaction (FSI) is a very interesting multi-physical problem with many technical applications as aeroelasticity, hydroelasticity, operating life of compressor blades, see [1], [2], as well as many biomechanical applications, e.g. hemodynamics or human phonation, [3], [4]. The understanding of the FSI phenomena and their further subsequent optimization depend on an accurate solution of these FSI problems. Since the FSI problems have very rarely an analytical solution, reliable and fast numerical simulations are of great interest. However, numerical simulation of FSI is a computational challenge since the governing equations for solid and fluid regions are different and the algorithm needs to solve the locations of solid–fluid interfaces simultaneously with the dynamics in both regions where the kinematic (e.g., non-slipping) and dynamic (e.g., stress matched along the normal to solid–fluid interface) boundary conditions are imposed at the interface.

The FSI problem is the non-linear problem consisting of two subproblems, the fluid flow problem and the elastic structure deformation problem, coupled by the boundary conditions at the common interface. The coupling between the two subproblems often leads to a highly nonlinear problem, which needs to be addressed by a robust and efficient solution technique, see [5]. The numerical treatment of FSI problems is complicated due to the need for mixed description — while the spatial (Eulerian) description is suitable for the fluid, the natural description for the solid part is the material (Lagrangian) description. Thus usually a mixed arbitrary Lagrangian–Eulerian (ALE) description is employed introducing an additional nonlinearity into the resulting equations. In general, the numerical schemes for solving the FSI problems may be classified into two approaches: the monolithic method, see e.g. [6], [7], [8], [9], [10] or the segregated method, see e.g. [11], [12], [13], [14]. The monolithic solvers are usually considered to be more robust, whereas computationally more demanding. Their advantage is that interface conditions can be automatically satisfied as these are included in the fully coupled (weak) formulation, see e.g. [9]. In a monolithic approach, the system of nonlinear algebraic equations arising from the discretizations is solved typically using a variant of Newton’s method, see e.g. [6], [8]. In [7] a block preconditioner is presented for the efficient solution of the linear systems that arise when employing Newton’s method.

Even though the monolithic – or fully coupled – approach is usually considered to be more robust, still the segregated approach, or partitioned scheme, is widely used for multi-physics problems due to its modularity: the specialized well-tuned solvers for the subproblems can be coupled only by the exchange of the data at the interface. In practice, usually the Dirichlet–Neumann (DN) partitioned scheme is used, see e.g. [11], where a simple and robust fixed-point FSI solver with dynamic relaxation was presented suitable for a wide range of applications. On the other hand the DN partitioned schemes are sensitive to the artificially added-mass effect which can cause the coupling iteration to diverge, see e.g. [15], [16]. In the last decade, the interface quasi-Newton method was developed in order to control the added-mass effect and simultaneously to speed up convergence, see e.g. [17]. Although such approach has several advantages, still in many applications the standard DN coupling is used, see e.g.[18], where the numerical investigation of free dynamic response of a rigid, symmetric airfoil was performed using standard commercial software. It was confirmed in [18], in particular, that the choice of the numerical parameters – such as domain size, grid size, and time step – can have a big influence on the accuracy of the numerical results.

In all the above mentioned approaches, the numerical procedure used for evaluation of the aerodynamic quantities seems to be somehow overlooked. In many papers the procedure of the evaluation of the aerodynamic forces (AF) is not specified even though this is one of the key components of the algorithm and there exist several possibilities of how the AF can be numerically approximated, [19]. The AF is mathematically defined as the surface integral of the AF density represented by the fluid stress tensor. From the mathematical point of view the evaluation of the surface integral is not well established for the standard choice of the functional spaces L2 and H1 for pressure and velocity, respectively. From the numerical evaluation point of view, if we omit the possible discontinuities of the pressure, the viscous part of the tensor determined by the symmetric part of the velocity gradient needs to be approximated e.g. by numerical differentiation formula and it is called the derived quantity, i.e. the quantity which can be obtained by postprocessing of the primary variables, see [19], [20]. However the numerical differentiation, or more precisely the straightforward extrapolation of the stress tensor components from the interior domain on its boundary, can significantly pollute the solution by numerical errors and it possibly reduces the theoretical accuracy of the method.

In literature there are known two approaches applicable for the AF evaluation within the finite element method (FEM) framework of how to overcome the above mentioned possible loss of accuracy. The first approach is motivated by the superconvergence phenomenon observed for FE approximations of a scalar linear problem, see e.g. [21], [22], [23]. Here, it was shown that derivatives of the approximate solution at appropriately chosen points of the elements (e.g. in triangle centroids) have a higher accuracy compared to the standard error estimates. This phenomenon is called superconvergence. One way how to exploit is to compute quantities of interest at these points [23] or to use the averaging of derived quantities directly computed from numerical solution over local patch of surrounding elements, see [24]. A technique followed here is to use the least-square method instead of pure averaging, see e.g. [25]. These procedures are common components of the FEM commercial softwares like Comsol or Abaqus, etc. [26], [27], and they improve theoretical accuracy for sufficiently fine meshes, see [23]. The second alternative approach of the AF evaluation is provided by their weak reformulation, which is sometimes called the consistent force formulation, see [28]. In this case the evaluation of the forces is rewritten with the aid of the Green theorem in combination with a suitable test function in the form using only volume integrals consistent with the weak form of the problem. Therefore, this method retains the accuracy of FE velocity approximation, see error analysis in [20]. The application of this method for a simplified FSI problem was published e.g. in [15], [29]. This approach seems to be also available in the Comsol toolbox, see [26], where however the detailed mathematical description is missing.

There are several FSI problems motivated by biomechanical applications e.g. by the flow-induced vibrations of vocal folds, see e.g. [30], [31]. Such a problem can be modeled as a FSI problem formed by the elastic subproblem described using the linear elasticity theory and by the fluid flow modeled by the incompressible Navier–Stokes equations in the ALE formulation allowing to incorporate the effects of the time dependent fluid domain, see [32]. The FEM is used for the numerical approximation of both FSI subproblems and the strongly coupled partitioned approach is implemented in order to increase accuracy and robustness of the numerical method especially for the simulation of the flutter onset, see [15], [31], [33]. Let us mention that in this case the classical Theodorsen method is not applicable due to complicated geometry of vocal folds, but the linear stability theory, which in the discrete version is based on the eigenvalues determination of associated system matrix, was applied e.g. in [30], [34] using several simplifying assumptions (as simplified vocal folds geometry, the potential flow model and a simplified AF). This simplified approach was compared with the FSI simulation based on the continuum flow model discretized by the FEM in [35], [36] with a good yet partial correspondence. The determination of aeroelastic stability boundary represented by the flutter velocity is an important task also in engineering areas like e.g. airfoils and turbomachines, [1], [2].

In this paper we pay the main attention to the problem of numerical approximation of the AF namely in the FSI problems, where the effect of the accuracy of approximation of AF can be expected to be substantially stronger due to the substantial dependence of the interface position on the acting AF and vice versa. Here, the three possibilities of AF calculations are described motivated by the above mentioned techniques, and these numerical procedures are tested on several benchmark problems similar to those in [37]. First, the classical DFG benchmark problem is considered, see [38]. However, for the FSI problems the evaluation of the AF can be strongly influenced by the substantial dependence of the interface position on the acting AF and vice versa, which brings motivation for the second test case of harmonically oscillating structure. Further, the AF evaluation techniques are applied for solution of a relevant FSI problem motivated by [36] and the influence of the applied techniques on the flutter velocity is tested. The main novelty of the present paper lies in significantly extended comparison of AF calculations studied for the benchmark fluid problems in [37] to the FSI problems motivated by the flow-induced vibrations of vocal folds, see e.g. [30], [31], [34].

The paper is structured as follows: The mathematical description of the FSI problem is given in the first section of this article. The second section contains brief description of the numerical scheme with a special attention paid to the three considered techniques of aerodynamic force calculation (within the FEM). The numerical results presented in the third section compare different AF calculations on four test cases. The first two cases consider static and dynamic benchmark problems of a cylinder in a cross-flow, while the next two FSI simulations investigate the structural energy dissipation and accumulation including determination of the critical inlet velocity for the considered problem setting.

Section snippets

Mathematical model

Let us consider a two-dimensional FSI problem consisting of an elastic body represented by domain Ωts and flow domain represented by domain Ωtf and their interface ΓWt, see Fig. 1, where the time dependences are highlighted by lower index t. At t=0, both the domains are given as part of initial conditions of the FSI problem and they are called reference domains, i.e. Ωrefs=Ω0s, Ωreff=Ω0f. Since the Lagrange coordinates are used for description of the elastic body, the time index is further

Numerical model

For the spatial discretization of the subproblems (1), (5) the FEM is utilized. For the purpose of time discretization the time interval [0,T] is divided into N equidistant time steps tn=nΔt with Δt=TN and the Newmark and backward differentiation formula (BDF) is used.

Numerical results

The described aerodynamic forces computation procedures are tested on four tests and the results are compared. First, these techniques are applied for computation of the drag and the lift coefficients for the flow around the cylinder and their convergence on refined meshes is investigated. Furthermore, a test with the oscillating cylinder is realized testing the methodology on moving meshes. This test is followed by two FSI cases motivated by simulation of human vocal folds. In the first FSI

Conclusion

This paper presented numerical techniques applicable for evaluation of aerodynamic forces in the context of a FSI problem. The mathematical model of a fluid–structure interaction (FSI) problem was presented and its numerical solution based on the finite element method was described, where the so-called minielement, i.e. P1-bubble/P1 element, was used for approximation of the velocity and pressure. Special attention was paid to different numerical realizations of aerodynamic force (AF)

CRediT authorship contribution statement

Jan Valášek: Formal analysis, Investigation, Methodology, Software, Validation, Visualization, Writing – original draft, Writing – review & editing. Petr Sváček: Conceptualization, Funding acquisition, Methodology, Supervision, Writing – review & editing.

Acknowledgments

Authors acknowledge support from the ESIF, Czech Republic, EU Operational Programme Research, Development and Education, Czech Republic, and from the Center of Advanced Aerospace Technology, Czech Republic (CZ.02.1.01/0.0/0.0/16_019/0000826), Faculty of Mechanical Engineering, Czech Technical University in Prague as well as support provided by Premium Academiae of Šárka Nečasová provided by the Czech Academy of Sciences .

References (48)

  • BodnárT. et al.

    Fluid–Structure Interaction and Biomedical Applications

    (2014)
  • LasotaM. et al.

    Anisotropic minimum dissipation subgrid-scale model in hybrid aeroacoustic simulations of human phonation

    J. Acoust. Soc. Am.

    (2023)
  • HeilM. et al.

    Solvers for large-displacement fluid–structure interaction problems: Segregated versus monolithic approaches

    Comput. Mech.

    (2008)
  • HronJ. et al.

    A monolithic FEM/multigrid solver for an ALE formulation of fluid–structure interaction with applications in biomechanics

  • WangY. et al.

    An optimal control method for time-dependent fluid–structure interaction problems

    Struct. Multidiscip. Optim.

    (2021)
  • KüttlerU. et al.

    Fixed-point fluid–structure interaction solvers with dynamic relaxation

    Comput. Mech.

    (2008)
  • BazilevsY. et al.

    Isogeometric fluid–structure interaction: theory, algorithms, and computations

    Comput. Mech.

    (2008)
  • ValášekJ. et al.

    On the application of acoustic analogies in the numerical simulation of human phonation process

    Flow Turbul. Combust.

    (2019)
  • FösterC.

    Robust methods for fluid–structure interaction with stabilised finite elements

    (2007)
  • E.H. van Brummelen, Added mass effects of compressible and incompressible flows in fluid–structure interaction, J....
  • SpenkeT. et al.

    A robin-neumann scheme with quasi-newton acceleration for partitioned fluid–structure interaction

    Internat. J. Numer. Methods Engrg.

    (2022)
  • EbrahemM. et al.

    Numerical investigation of transient response of a coupled two-degrees-of-freedom symmetric airfoil before flutter

    Int. J. Aeroacoust.

    (2018)
  • GreshoP.M. et al.

    Incompressible Flow and the Finite Element Method

    (1998)
  • V. John, M. Tabata, L. Tobiska, Error estimates for nonconforming finite element approximations of drag and lift in...
  • View full text