{'fraction_non_alphanumeric': 0.033331626644821054, 'fraction_numerical': 0.014515385796938202, 'mean_word_length': 4.69265520256485, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Application of foam-extend on turbulent fluid-structure interaction

Turbulent flow around flexible structures is likely to induce structural vibrations which may eventually lead to fatigue failure. In order to assess the fatigue life of these structures, it is necessary to take the action of the flow on the structure into account, but also the influence of the vibrating structure on the fluid flow. This is achieved by performing fluid-structure interaction (FSI) simulations. In this work, we have investigated the capability of a FSI toolkit for the finite volume computational fluid dynamics software foam-extend to simulate turbulence-induced vibrations of a flexible structure. A large-eddy simulation (LES) turbulence model has been implemented to a basic FSI problem of a flexible wall which is placed in a confined, turbulent flow. This problem was simulated for 2.32 seconds. This short simulation required over 200 computation hours, using 20 processor cores. Thereby, it has been shown that the simulation of FSI with LES is possible, but also computationally demanding. In order to make turbulent FSI simulations with foam-extend more applicable, more sophisticated turbulence models and/or faster FSI iteration schemes should be applied.


Introduction
Fatigue is responsible for the majority of failures in structural and mechanical components [1]. One of the many causes of fatigue loading is the action of a turbulent flow field on a structure, e.g. turbulent flow through water turbines [2], wind turbines [3], or around flexible structures like marine risers [4][5][6][7]. In order to assess the safety and reliability of these structures, it is necessary to estimate their fatigue life, which depends on the actual stress history experienced at their critical spots.
While being subjected to the turbulent flow, the structures will vibrate. If the structure is very stiff compared to the fluid, the influence of its vibration on the fluid flow will be negligible [8]. Then, the stress history may be found by first evaluating the flow field around or through a rigid model of the structure, and then imposing the resulting fluid load history to a structural dynamics analysis of the structure [9]. This is called an uncoupled simulation. However, slender structures like marine risers are flexible. Therefore, they will be subjected to significant deformations, thereby influencing the fluid flow. In these cases, it will be necessary to perform coupled fluid-structure interaction (FSI) simulations, in which the deformed shape of the structure is explicitly taken into account when calculating the fluid flow, in order to obtain realistic estimations for the stress history [10,11].
Because of their complexity, FSI problems are normally solved using computational techniques. These techniques may be divided into two groups: Monolithic approaches and partitioned approaches. In monolithic approaches [12,13], a unified mathematical model is formulated for the entire regime (fluid and structure), and solved as a single large system of equations [14] for each time step. Stability and convergence rate may be better in a monolithic approach, but the system matrices may become illconditioned, requiring appropriate preconditioners [12].
In partitioned approaches [14][15][16][17] on the other hand, the structure and fluid regimes are each solved independently by separate mathematical models, i.e. computational solvers. At the interface between the fluid and the structure, information is interchanged between the two solvers, thereby enabling the fluid-structure interaction to be taken into account. The FSI problem is solved sequentially, as shown in Figure 1. The partitioned approach makes it possible to use existing, well-developed solvers to compute each of the two regimes, e.g. to use a finite element solver to compute the structural vibrations and a finite volume solver to compute the fluid flow. Because of this flexibility, the partitioned approach is very popular, and is employed in the majority of the FSI research reported [15].
For the partitioned approaches, it is possible to employ either weak or strong coupling [12,15]. When weak coupling is employed, the fluid and structure regimes are each solved once at each time step. However, this scheme does not guarantee that the instantaneous equilibrium between the fluid and structure is fulfilled. In strong coupling, this is ensured by iteratively computing the fluid and structure fields until instantaneous equilibrium is satisfied, as indicated in Figure 1. According to Münsch and Breuer [18], weak coupling leads to severe stability problems in most applications, and strong coupling should therefore be employed, especially for flexible structures.
Fluid-structure interaction has been studied for a number of flexible structures. Holmes, et al. [5] have simulated the vortex-induced vibration of marine riser models using the AcuSolve finite element (FE) computational fluid dynamics (CFD) solver. In these simulations, the motions of the riser are computed by combining the CFD results with the eigenvalues and eigenvectors of the riser. Holmes, et al. propose to use a detached eddy simulation (DES) or an unsteady Reynolds averaged Navier-Stokes (URANS) turbulence model, but do not specify which turbulence model they have used. The same riser model has also been simulated by Kamble and Chen [6], using the finite volume (FV) CFD solver FANS3D, together with the chimera (overset) grid technique, in which a very fine mesh around the riser overlaps the mesh for the rest of the domain. Using a large eddy simulation (LES) model and about two millions cells, together with a finite difference scheme to simulate the curvature of the riser, good correlation with experimental results has been obtained.
Vortex-induced vibration of marine risers has also been simulated by Menter, et al. [19], by using the element-based finite volume software ANSYS-CFX to compute a laminar simulation of the flow, and the ANSYS structural FE software to compute the structural vibrations. Nguyen and Nguyen [7] used the FV CFD solver OpenFOAM to simulate vortex-induced vibration of rigid, elastically mounted circular cylinders, using a DES turbulence model. Other simulations of offshore structures include the uncoupled simulation of wave actions on a simplified platform deck structure by Marzban, et al. [20], and the strongly coupled simulation of the movement of a closed riser filled with mud by Paczkowski, et al. [21].
In the case of turbulent FSI, several test cases have been simulated by Münsch and Breuer [18], using strong coupling between a FV CFD solver and a FE structural solver, employing LES. The simulation results have only been compared to other computational results or laminar experimental data, however. Nilsson, et al. [22] have employed weak coupling and LES to simulate structural excitation of thin plate structures subjected to acoustics and gas flow in OpenFOAM, commenting on some short-comings. Recently, Cesur, et al. [16] studied the flow past a flexible cantilever beam, employing implicit large eddy simulation and strong FSI coupling, providing a thorough discussion of the results.
As indicated by the above literature review, the majority of FSI studies employ a finite element solver for the structural dynamics. However, recently a FSI toolkit was developed for foam-extend, a branch of the open source finite volume CFD code OpenFOAM [23], by Tuković, Cardiff, Karač, Jasak and Ivanković [24]. The toolkit employs the finite volume method (FVM) to solve both the fluid and structure regimes, in a partitioned approach, and is based on their previous works [25][26][27][28][29][30]. The toolkit supports parallel computation of FSI problems with weak or strong coupling between a laminar or turbulent incompressible flow and an elastic structure subjected to large deformations. Ponweiser, et al. [31] have investigated the scalability of this toolkit for multiple processors, and Šekutkovski, et al. [32] have used the toolkit together with a hybrid RANS-LES turbulence model, in order to simulate flow around aircraft wings.
The objective of this paper is to investigate whether the FSI toolkit for foam-extend 4.0 is readily applicable to predict the fatigue life of flexible structures subjected to a turbulent flow. Available test cases for FSI are limited. Turek and Hron's test case [33] is for laminar flow, while Pereira Gomes and Lienhart's turbulent test case [34] employs two structural materials of different densities, which is currently not supported in the FSI toolkit. Therefore, it is chosen to implement a LES turbulence model to Richter's laminar test case of a flexible wall in an enclosed channel flow [13], which also is included as a tutorial in the toolkit. This test case has also been simulated by Mehl, et al. [14] and Gillenbaart, et al. [17].
LES is chosen as the turbulence model, because Reynolds-averaged Navier-Stokes (RANS) models do not resolve the turbulent fluctuations which may be important for fatigue evaluation, while direct numerical simulation (DNS) has a high computational cost [16]. While previous FSI simulations with LES have often employed the Smagorinsky model [6,18,22], the one-equation k SGS eddy viscosity model is used in the current work. The goal is to obtain the stress history for a critical spot on the flexible structure. The stress history should be recorded for a sufficient amount of time, so that it is representative for the full fatigue life of the flexible structure, and be usable for fatigue life evaluation.

Governing equations
The governing equation for a solid structure is normally stated in the Lagrangian formulation, in which the reference frame moves and deforms with the structural deformations. On the other hand, the governing equations for fluid flow are normally stated in the Eulerian formulation, in which the reference frame remains stationary, independently of the fluid flow. Like the majority of partitioned FSI solvers [15], the FSI toolkit for foam-extend employs the Arbitrary Lagrangian-Eulerian (ALE) formulation for the fluid flow. This means that the reference frame for the fluid domain is independent of the fluid flow (Eulerian), but deforms according to the deformations of the structure (Lagrangian). In practice, this also means that the fluid mesh deforms according to, and conforms to, the deformed shape of the structure.

Governing equations for the fluid flow
The isothermal flow of an incompressible and Newtonian fluid may be fully described by the continuity equation and the Navier-Stokes equations. However, in order to solve these equations accurately with even the smallest turbulent eddies and fastest fluctuations resolved, an extremely fine spatial mesh and very small time steps would be required. Instead, we use large eddy simulation (LES), where only the large-scale turbulent eddies are resolved. This makes it possible to use a coarser mesh and larger time steps, decreasing the computational effort.
In LES, the small-scale turbulent eddies, which are smaller than the cutoff width Δ (normally equal to the cell size), are not resolved [35]. Instead, they are modelled by sub-grid-scale (SGS) stresses, τij, which represent the effects of the unresolved motions. This spatial filtering decomposes the variables φ (pressure and velocity components) into a sum of a filtered component, φ , and a sub-grid-scale component, φ′ , i.e. φ φ φ′ = + [36]. The governing equations for the fluid flow are solved for the filtered components, in which case the incompressible LES continuity equation in the ALE formulation becomes [15]: is the filtered fluid velocity field, um is the velocity of the reference frame (i.e. the mesh), and x = [x1, x2, x3] = [x, y, z] are the Cartesian coordinates. The filtered and incompressible Navier-Stokes equations in the ALE formulation are given as [16,25,28,35,36]: where ρf is the density of the fluid, t is time, p is the filtered fluid pressure field, µf is the dynamic viscosity of the fluid and τij are the SGS stresses. For simplicity, the effects of gravitation are neglected in the simulation.
The SGS stresses are given as [35]: The first term in this equation cannot be calculated from the solution of the filtered governing equations. Instead, it has to be modelled. The simplest SGS model is the Smagorinsky-Lilly model, which is often used, e.g. [6,18,22,36]. However, this model contains a constant which needs to be adjusted on a case-to-case basis [35]. Therefore, the one-equation kSGS eddy viscosity model will be used in the current work. This model includes one additional transport equation, which is solved for the SGS turbulent kinetic energy, kSGS [35]: where the SGS viscosity is given as SGS , the local strain rate of the resolved flow is given as ( ) ∂ ∂ + ∂ ∂ , and the rate of dissipation is given as C ε are constants, while Δ is the filter cutoff width. With the SGS viscosity known, the SGS stresses may be calculated from: where δij is the Kronecker delta. The van Driest damping function [37] has been used in the simulations to reduce the near-wall eddy viscosity.

Governing equation for the structure
The governing equation for the structure is the conservation of linear momentum. For an isothermal structure in the total Lagrangian formulation, taking large deformations into account, this equation may be written as [26,38]: where ρs0 is the initial density of the structure, us is its velocity vector, Ω0 is its initial volumetric domain, Γ0 is its initial surface, and n0 is its outward pointing unit normal vector. S is the applied second Piola-Kirchhoff stress, fb is any body force acting, and F is the deformation gradient, given as: where ds is its displacement vector and 0 ∇ represents the gradient operator at the undeformed geometry.
The second Piola-Kirchhoff stress S may be related to the more commonly used Cauchy stress σ (current force acting on current geometry) by the following equation: For an elastic St. Venant-Kirchhoff material, the constitutive equation is: where E is the Green-Lagrangian strain and µs and λ are the Lamé parameters [38]: E is the elastic modulus and s ν is Poisson's ratio for the material.

Fluid-structure interface
At the fluid-structure interface, the following conditions must be met [15]: Here, uf is the fluid velocity field, us is the structural velocity, n is the unit normal vector of the interface, and σs and σf are the Cauchy stress matrices for the solid and fluid, respectively.

Discretisation of the fluid equations
The FSI toolkit for foam-extend 4.0 is used to solve the governing equations numerically, using the finite volume method. For the fluid flow, the PISO (Pressure Implicit with Splitting of Operators) scheme is used in the present study. The SGS stresses calculated as part of the LES model are similar in magnitude to the numerical truncation errors [35]. It is therefore important to reduce these errors, so that they do not swamp the SGS stresses. This is done by using second-order or higher-order discretisation schemes in time and space. The recommended temporal discretisation scheme for LES is the Crank-Nicolson scheme [39], while the backward scheme may give unacceptable dispersion [22]. The Crank-Nicolson scheme is an implicit scheme, which is second-order accurate in time, by weighting field variables at two succeeding time steps equally [35]. OpenFOAM (and foam-extend) allows the Crank-Nicolson scheme to be mixed with the first-order implicit Euler scheme. It is recommended to use a mix of 90 % Crank-Nicolson and 10 % Euler scheme in OpenFOAM, in order to bound/stabilise the scheme [40]. However, investigations by the authors have shown that using 80 % Crank-Nicolson and 20 % Euler works better in foam-extend, as illustrated by pure CFD simulations in Figure 2. The Courant number is given as Co = |uf|•Δt/Δx, where |uf| is the magnitude of the local velocity vector, Δt is the length of the time step and Δx is the length of the finite volume cell. The maximum Courant number should be kept below 1.0 when using the transient PISO scheme [41]. Here, a constant time step length has been used to keep the maximum Courant number around 0.53.
For the spatial discretisation of the convection and diffusion terms, and the gradient, central differencing with Gaussian integration is used. More specifically, OpenFOAM's scheme "Gauss filteredLinear2 1 0" was used for the div(phi,U) term, in order to remove staggering caused by pressure-velocity decoupling. Further details of the discretisation schemes (in the OpenFOAM branch) may be found in [40].

Discretisation of the structure equations
The finite volume discretisation of equation (6) contains one temporal term, one diffusive term and two source terms [30].
The first-order fully implicit Euler temporal discretisation scheme was first added to the FSI toolkit, and has been used in the current work. The second-order backward scheme has later been added to the FSI toolkit [24], and should be used for increased accuracy [17].
The diffusive term is discretised using central differencing with Gaussian integration, and gradients are computed using the least-squares approach. Further details of the discretisation of the governing equation for the structure may be found in [30].

Discretisation of the fluid-structure interaction
The fluid-structure interaction is solved in a partitioned approach, with strong coupling, as illustrated in Figure 1. The displacement increment and velocity of the structure are transferred to the fluid, by deforming the fluid mesh, while the pressure and viscous forces are transferred from the fluid to the structure. This is the most common procedure, as the reversed exchange of information might lead to instability [15]. Hence, the fluid forces are applied on the structure as natural boundary conditions.

Under-relaxation.
As part of the iterative strong coupling, the fluid mesh is deformed at each iteration step. However, in order to stabilise the iterations and enforce convergence [14], the fluid mesh is not moved the same distance as the predicted structural displacement at the last iteration step. Instead, the iteration is under-relaxed. This means that the displacement of the fluid mesh at iteration step i+1 is calculated as: where di is the displacement of the fluid mesh at iteration step i, 1 i+ d  is the structural displacement predicted in iteration step i, and ωi is called the relaxation factor.
Three relaxation schemes are available in the FSI toolkit for foam-extend [24]: Fixed relaxation (i.e. ωi = constant), Aitken dynamic relaxation and IQN-ILS. The fixed relaxation often leads to slow convergence [14], and dynamic relaxation is therefore employed. The convergence rate of both the dynamic relaxation techniques has briefly been tested, but no significant difference was observed. The Aitken dynamic relaxation method is one of the most popular methods [42], and has been employed in the current work. Details of the method may be found in [43]. In the current work, an initial relaxation factor of 0.4 is used at each time step, and iterations are made until a relative tolerance of 10 −5 is reached. If this tolerance is not reached within 200 iterations, the solver simply moves on to the next time step.

Mesh movement.
The motion of the fluid mesh is computed by a Laplacian approach, in which the points are moved according to the Laplace equation [15,25]: where γ is the diffusion coefficient, which is chosen to be proportional to the squared inverse distance from the moving boundary. The solver "refVelocityLaplacian" is chosen, with 2 non-orthogonal correctors.

Geometry, properties and boundary conditions
The geometry of the problem chosen for this study is shown in Figure 3. It is based on Richter's laminar test case [13], but has been extended to full three-dimensionality, in order to allow turbulence modelling. The geometry of the problem consists of a flexible wall (structure) inside a rectangular duct. The duct has one inlet and one outlet. Its width is 0.8 m and its height is 0.4 m. The flexible wall is 0.4 m wide, 0.1 m thick and 0.2 m high, and is fixed to the floor of the duct, as shown in Figure 3. While Richter's case had an inlet length, IL, of 0.45 m, preliminary studies have shown that this is too short for the turbulent case, and therefore it has been increased to IL = 0.95 m (units are here added to Richter's non-dimensional measures).
As in Richter's case, the duct conveys an artificial fluid of density ρf = 1000 kg/m 3 and kinematic viscosity νf = 10 −3 m 2 /s. In Richter's case, the fluid enters the duct with mean velocity U = 0.15 m/s. This velocity is here increased to U = 5 m/s, in order to get turbulence in the wake of the flexible wall. In this case, the Reynolds number, Re, for the flow entering the duct is 2667, which is within the critical zone between fully laminar (Re < 2000) and fully turbulent (Re > 4000) flow [44]. This value should be sufficient to give turbulence in the wake of the flexible wall, however.
For simplicity, a laminar inlet velocity profile, i.e. a parabolic velocity profile, has been used as an approximation. The peak velocity in the middle of the cross-section is umax = 2U = 10 m/s. The no-slip condition is applied for the bottom, top and sides of the fluid domain, as well as for the interface to the flexible wall. At the outlet, the zero gradient condition is imposed on the velocity, and the pressure has a Dirichlet boundary condition with zero pressure. The properties of the flexible wall have been changed from those in the Richter's case, in order to obtain visible, but yet not extreme deformations. The wall is too thick to obtain significant vibrations with any realistic material. Its density, ρs, has been set equal to 8000 kg/m 3 , Poisson's ratio, νs, has been set to 0.3, while the elastic modulus has been chosen to be E = 4 MPa, which is very soft. The flexible wall is fixed (ds = 0) at the bottom, but is otherwise free to move.

Initial conditions
It is chosen to start the solution from a fully developed turbulent flow, and the structure at rest; a technique which has also been used in [18,42]. The chosen initial velocity field is shown in Figure 4. The structure is allowed to move during the first computational step. This initial condition is preferred, because it requires little computational effort, compared to the alternatives. However, as illustrated in [18,42], the sudden exposure of the flexible structure to a fully developed fluid flow induces an initial transient to its vibration. In the context of fatigue life estimation, we are mostly interested in the stabilised vibrations remaining when this transient effect has been damped out. In order to reduce the transient effect, the inlet velocity could have been slowly increased to the full flow. However, as very many FSI iterations are required even at low velocities, this procedure is computationally demanding. Furthermore, the PISO solver implemented in the FSI toolkit does not contain any ability to automatically adjust the time step to keep the Courant number within a given limit, meaning that the time step would have needed to be manually decreased as the velocity increased.

Computational meshes
While Richter [13] uses finite elements to discretise the computational domain, both Mehl, et al. [14] and Gillebaart, et al. [17] use hexahedral finite volume cells which are approximately 2.5 cm long in each direction, for the flow domain. In a fully three-dimensional model, this corresponds to 40 448 finite volume fluid cells. This is quite a coarse mesh for a large eddy simulation. In the current study, the cell size has therefore been decreased to 1.0 cm in each direction, resulting in a total of 632 000 cells in the fluid domain. In Figure 5, the vibration history for the flexible wall, computed by both meshes, is plotted for 1.5 seconds, and it is clear that the coarser mesh is not able to capture all the small eddies which significantly influence the vibration of the flexible wall. The fine mesh, which is used to obtain the remainder of the results presented, is shown in Figure 6.  The flexible wall has in all cases been discretised by 512 hexahedral finite volume cells, each 2.5 cm long in each direction. Both meshes have been decomposed for 20 processors using the metis algorithm. A computational time step Δt of 250 µs has been used to keep the maximum Courant number around 0.53. Simulation results have been saved every 0.01 s, i.e. at a sampling frequency of 100 Hz.
In order to assess the quality of the simulation, the dependence of these parameters should be further investigated. Such studies could not be conducted within the present work however, due to time and CPU constraints.

Simulation results
Results from the simulations are illustrated in this section, and further discussed in section 6. If otherwise is not mentioned, the figures show the computed fields at the plane z = 0.4 (ref. Figure 3), which cuts through the middle of the flexible wall. Finite volume results are saved in the nodes in the centre of each cell, and the field values are thereby constant throughout each cell. However, volume point interpolation has been used to smooth out the fields shown in the figures. For the remainder of the simulation, the 80-core computer was subjected to an increased workload, and did also randomly crash. Continuous log files are therefore not available between t = 1.5 s and t = 1.88 s. Between t = 1.88 s and t = 2.32 s, the computation time was 96 hours (clock time) using 20 cores, i.e. 130 min 51 s per 0.01 second simulated, on average. As the average number of FSI iterations was more or less the same as for the first 1.5 s, the increased computation time is thought to result from the increased workload on the computer, which had around 120-150 concurrent jobs at this time.

Computation time
For reference, similar simulations without FSI did only require 33 s clock time per 0.01 second simulated, using 20 processor cores.

The fluid flow
The local velocity reaches its maximum at t = 0.079 s. The nearest time step which is saved is t = 0.08 s, and the velocity field at this time step is shown in Figure 8(a). It is clear that the flexible wall has become significantly curved, and that the relatively smooth flow path (continuous red field in Figure 4) has been disturbed.
The velocity field for the last time step is shown in Figure 8(b). The wall is not as curved as at t = 0.08 s, but the amount of turbulence has increased. A lot of mixing is taking place in the wake of the wall, and the velocity profile is no longer smooth upstream the wall either.  Figure 9 for t = 1.80 s and for t = 2.32 s. It is observed that these two fields look relatively similar; indicating that a stabilised flow has been reached. The increase in turbulence due to the vibration of the flexible wall is also evident from Figure 10, where the SGS kinematic viscosity fields are shown for t = 0 and t = 2.32 s. The flow downstream the flexible wall is clearly turbulent, and the upstream turbulence has also been increased. The turbulent field at t = 2.32 s is illustrated by isosurfaces of the second invariant of the velocity gradient tensor, Q, in Figure 11. In addition to the turbulence downstream the flexible wall, it indicates that a spiralling flow pattern has developed near the inlet.

The vibration of the flexible wall
The deflection history for the point (x, y, z) = (0.45, 0.2, 0.4), which is positioned at the middle of the upper front edge of the wall, is shown in Figure 12. The transient initial vibration due to the sudden freedom of the wall at t = 0 is clearly visible. The time history is too short to determine whether these transient effects have disappeared during the simulated interval, but the curves do indicate a possible stabilised vibration starting a bit before t = 2 s. With a sampling frequency of 100 Hz, the deflection magnitude history has been found to be dominated by frequencies around 7.3 Hz, with lesser peaks at 2.2 Hz and 27 Hz.
For the purpose of estimating the fatigue life of a structure, it is essential to know its stress history. The stress state at a point is often represented by the von Mises equivalent stress, which may be calculated as: The von Mises stress distribution in the wall at t = 2.32 s is shown in Figure 13. The asymmetry along the z-axis is assumed to be a result of the turbulence of the flow. As illustrated in Figure 13, the most critically loaded part of the wall is the strip along the front, close to the bottom. It is likely that any fatigue crack will nucleate along this strip. The time histories of the different stress components at the middle of the strip, (x, y, z) = (0.45, 0.0125, 0.4) are shown in Figure 14.  Figures 12 & 14 indicate that the FSI toolkit for foam-extend 4.0 is capable of producing a vibration history, which is appropriate for evaluating the fatigue life of the structure. However, some important aspects need to be taken into account, and will be addressed in this section.

The initial transient
As predicted in Section 4.2, the sudden application of fluid-structure interaction in the fully developed flow induces an initial transient to the vibration history, as observed in Figure 12. As momentum is dissipated from the structure to the fluid, the vibration of the wall, due to the sudden application, is damped, and will diminish over time. This means that a "steady-state" or stabilised vibration will be reached after a certain time. If the wall is subjected to the turbulent flow for most of its life, it is this stabilised vibration which is interesting for fatigue purposes. When the vibration is stabilised, a relatively short vibration response may be extrapolated to represent the full vibration history of the wall. Then, its fatigue life may be estimated by using a S-N curve for the material and applying a cycle counting technique and a cumulative damage law [1], like the well-known linear Palmgren-Miner rule [45] or more recent, nonlinear rules, e.g. [46].
However, as shown in Figures 12 & 14, it is not possible to indisputably assess whether the transient effect has diminished completely within the 2.32 seconds simulated in the current work. There is still a lot of variation. Still, the last three major oscillations of the deflection magnitude shown in Figure 12 seem to have a relatively constant mean value and frequency, with only small deviations in the amplitude. This may indicate that the stabilised vibration either has been reached, or that it is coming. This indication is further strengthened by the observation made of the time averaged velocity fields shown in Figure 9, where the flow seems to have stabilised at about t = 1.80 s.
This indicates that even though the simulated interval of 2.32 s is too short to make any sensible fatigue life estimations, we are not very far away from having sufficient data to make such estimations. It is absolutely possible that having results all the way up to t = 3 or 4 s would have made such estimations possible. The major problem is the time required to obtain these results. Based on the obtained results, somewhere between 44 and 148 hours of computation would be required to reach t = 3 s.

The need for simulation of the fluid-structure interaction
Throughout section 5, it is noted that the amount of turbulence has increased during the time interval from t = 0 to t = 2.32 s, whereas the amount of turbulence stabilised in similar simulations where FSI was not employed. This turbulence increase may be directly related to the fluid-structure interaction. While the turbulence induces vibrations to the flexible wall, the vibrations of the wall do also induce additional disturbances and mixing, i.e. turbulence, to the flow. This illustrates the importance of applying (strong) fluid-structure interaction when simulating the vibration of flexible structures being subjected to a turbulent flow.

Assessment of the simulation results
As no experimental data is available for the case investigated, the validity of the simulation results cannot be confirmed. Furthermore, the authors are not aware of other implementations of LES to this case. The flexible wall is relatively thick and short, and can therefore be analysed neither as a plate nor as a cantilever beam. The next step would therefore be to perform simulations of a more realistic problem, where an analytical solution of the natural frequencies of the structure, and preferably, experimental data, are available.
Several sources of errors exist for the computed results. First and foremost, a full convergence check with respect to mesh refinement and time step length has not been conducted within the present work, due to time and CPU constraints. The results may therefore be sensitive to these parameters. Secondly, while the influence of the inlet length has been investigated, it is possible that the length from the flexible wall to the outlet may be too short. Also, the sampling frequency of 100 Hz may be too low to resolve the highest vibration frequencies. Furthermore, the accuracy of the vibration history may have been increased by using the second-order backward temporal discretisation scheme for the structure.
In addition, as the flow is confined by the duct on four sides, the results rely on the quality of the ν SGS wall function, which may be questionable. Further refinement of the mesh near the wall, so that the first cell centre would lie at a non-dimensional distance, y + , of ~1 from the wall (instead of y + ≈ 12, as used here), would make it possible to omit the use of the wall function [41]. The spiralling flow developing from a unidirectional inlet may indicate that some numerical problem exists.

Computational effort
With the resolution and settings used in the current work, between 64 and 218 hours of computation are required to compute one simulated second, depending on the workload of the computer. The computational cost is therefore high. However, as indicated above, even three to four seconds of data may be sufficient to estimate the full fatigue life of the structure. FSI simulations in foam-extend may be a good alternative to the empirical and analytical tools often used today [7].
Within the current FSI toolkit for foam-extend, a couple of parameters could be tuned which could possibly further decrease the computational effort. In the current work, the computation has been performed in parallel using 20 processor cores, each core handling 31 600 fluid cells and 25-26 structure cells, decomposed by the metis algorithm. It is possible that this is a low amount of cells per core, meaning that excessive time is used to transfer information between the cells. However, computations with even fewer fluid cells per processor core are found in the literature, e.g. Tuković and Jasak's [29] simulation of a droplet using 8 parallel processors, each handling about 7777 cells.
Nevertheless, some computation time could probably be saved by using fewer processor cores to compute the structure. Here, the toolkit's standard setting of using all of the same cores for the computation of both regimes was chosen for simplicity. The fluid and structure are usually computed in serial, i.e. in a staggered pattern, meaning that all the processor cores cannot be fully utilised when running a parallel computation, no matter how the structure is decomposed. Mehl, et al. [14] have proposed two techniques for simultaneous parallel computation of both regimes at each time step, which could be implemented to maximise the utilisation of each processor core.
A relatively large amount of the computation time is spent moving and deforming the fluid mesh [31]. In the current work, the FVM-based "refVelocityLaplacian" solver with 2 non-orthogonal correctors has been used to deform the mesh. The FEM-based solver developed by Jasak and Tuković [25] is also available. It is stated to be more robust, and may influence the computation time.
The number of FSI iterations does also contribute significantly to the computation time. With the Aitken scheme, an initial under-relaxation factor of 0.4 and required relative tolerance of 10 −5 , 17.7 FSI iterations were required at each time step, on average. This is a significant increase from the laminar case. Changing any of these parameters would influence the number of FSI iterations, and thereby the computation time. Limiting the number of FSI iterations will be of major importance for reducing the computational effort for simulating turbulent FSI. van Zuijlen and Bijl [42] have suggested a technique where the flow is solved on alternating coarse and fine meshes during the FSI iteration within each time step, in order to reduce the computational effort for RANS simulations. However, this technique has limited application to LES simulations, where the fine mesh always is required to give a sufficient resolution.
The last solution to reduce the computational effort we will mention is the use of a more sophisticated turbulence model. Implicit large eddy simulation [16], hybrid RANS-LES models [32] and detached eddy simulation techniques [7] have all been used for this purpose to solve turbulent FSI problems.

Conclusions
In this paper, a large-eddy simulation (LES) turbulence model has been implemented to Richter's fluid-structure interaction (FSI) test case, and 2.32 seconds of turbulent fluid-structure interaction has been simulated, using the FSI toolkit for foam-extend 4.0. It has been shown that this toolkit is capable of producing a vibration history, which is appropriate for evaluating the fatigue life of a flexible structure.
The computational effort required to perform such simulations is large, however. Even for the relatively simple case investigated here, 632 000 cells were required to discretise the fluid domain. An average number of 17.7 FSI iterations were required at each time step, due to the relatively slow FSI convergence rate. The FSI convergence rate is significantly lower for turbulent flows modelled by LES, than for laminar flows. Because of this, over 200 computation hours was required for the simulation of 2.32 seconds of turbulent FSI.
This indicates the need of applying more sophisticated turbulence models and/or faster FSI iteration schemes. The FSI toolkit for foam-extend could also require the possibility of automatically adjusting the time step according to the Courant number, and the implementation of additional solvers. Especially a two-phase solver for the simulation of ocean waves and a compressible solver to simulate water hammer induced vibrations of pipelines would be appreciated.
In future work, the factors influencing the computational effort, as mentioned in section 6.4, should be further studied. The studies should also consider a more realistic problem, with experimental and/or analytical data available.