Abstract
Loads due to wave impacts are a limiting factor in the design of liquefied natural gas (LNG) tankers and their insulation. The current methodology considers the load independent from the response of the tank. Better tanks can be designed by knowing the effect of the interaction between the wave loads and the response, however predicting these effects is computationally expensive. In this paper a new application of the non-hydrostatic shallow water equations are presented, namely as a reduced order model (ROM) for fluid structure interaction for wave impacts. Our ROM is compared to a high fidelity model. The proposed ROM is fast and accurately predicts the total impulse and added mass, and therefore the general behaviour of the structure during the free vibration phase. It does however not always accurately predict the maximum force. It is therefore considered an appropriate tool for a first screening of the loads for which fluid-structure interaction is important, after which a more accurate method can be used to evaluate the most interesting cases. A sensitivity study is performed for various impact angles and velocities, showing that the importance of fluid structure interaction depends highly on the specific situation.
Keywords
Wave impacts on inner walls are known phenomena in liquefied natural gas (LNG) sloshing on board of LNG tankers, and is characterized by high variability of the loads due to free surface instabilities before impact [2]. The LNG in the tanks is kept at its boiling point (−162°C), causing phase change in the liquid and cryogenic temperatures for the tank walls which should maintain structural integrity at all times. An insulating cargo containment system (CCS) is applied to the tank walls in order to prevent the steel hull from becoming brittle. For an optimal insulation the tank walls should be as light and thick as possible, for maximum strength they should be as thin and heavy as possible. These contradicting requirements combined with the variability of loads complicate the design of such a CCS.
The design loads of a CCS are typically determined using model tests on a 1:40 scale [11]. A rigid scale model of the LNG tank is fitted with pressure sensors and excited in a similar way as the real tank would. The pressures are statistically post-processed and applied to the structure in patches [21]. In this approach the dynamics and flexibility of the structure is not directly taken into account, only through dynamic amplification factors.
A large test campaign investigating wave impacts for LNG CCS is the ‘Sloshel’ project [2,3,15]. A full scale LNG CCS was placed in a vertical end wall of a flume. Then different shapes of waves loading the wall were obtained by changing positions of the focal point with respect to the wall. Depending on the wave shape, three elementary loading processes (ELPs) were identified [17]:
Direct impact: the first part of the wave hits the wall, and liquid compressibility is important Building jet: the fluid has to change direction and spreads along the wall Pulsating gas pocket: enclosed gas is compressed and causes a varying load
Each of these ELP’s has a different signature, physics, and scaling [8], which makes the aforementioned model tests difficult to translate to real tank applications. However, not all these ELP’s are of equal importance for the response of the structure. The response depends on the loaded area and time signature, as for instance demonstrated in [5]. Furthermore, if it is known which loads (or ELPs) are most important for the structural response, then the hydrodynamicist can focus on accurately predicting these loads.

The structure considered in this paper is the ‘Mark III’ CCS of the company GTT [10,13], which is a sandwich structure with a foam core and plywood face sheets, a cross-section is shown in Fig. 1. On the inside of the tank a stainless steel membrane is fitted to keep the LNG inside, and a secondary membrane is fitted between the foam layers for redundancy. In order to minimize thermal stresses, relaxation grooves are cut just above the secondary membrane, creating ‘subpanels’ of 340 by 340 mm which are to a large extent uncoupled [4]. Note that Fig. 1 represents a slice of the subpanel.
The CCS has two main failure modes: plastic indentation failure of the foam and shear failure of the bottom plywood [6,13]. If either of these two failure modes occur, then there is a high probability that also one of the membranes will fail, which is to be avoided at all times. A reduced order model for the onset on indentation failure is presented in [4], based on [20] and [24]. Prediction of the structural response can be significantly sped up using such a reduced order model, possibly leading to great simplification of the coupled problem.
The problem is however slightly more complex, as shown for instance in [18], where sloshing model tests including FSI are described for a flip-through impact. A 2D tank is used with an elastic flat wall at the end, which represents the CCS. The deformations of the wall are measured and computed (using pressures on the rigid wall). A comparison shows that the wall response requires more than just the pressure on the wall and added mass and damping, especially due to the change of added mass in time. This means for us that the response of the structure to the loads (one way coupling) is different from the fully coupled response (two way coupling). Hence, it can be necessary to take the two way coupling between the load and wall into account when investigating which loads are important.
The two-way coupling can be investigated for the entire tank, as demonstrated in [28]. Here the stiffness of the wall was varied and the tank excited in the same way. As the tank wall becomes increasingly compliant, the pressures decrease and show a more vibratory response. The total time during which the pressure sensor was loaded is longer. This is however numerically expensive and time consuming, especially considering the highly variable loads that require a probabilistic assessment.
For such a probabilistic assessment we look at simplified models, that are fast and sufficiently accurate for the situations in which they are used. The load due to a wave impact and therein encapsulated air can be predicted using two simplified models: the Wagner and Bagnold model [19]. First, the Wagner model can describe the building jet (ELP2) based on an initial shape and potential theory for the flow [25]. The load of the Wagner model starts with a high pressure peak and very locally, then spreading out and then spreads out over the wall while the total force increases. Second, the Bagnold model describes how entrapped gas pulsates (ELP3) based on a gas law and some moving mass representing the surrounding water [1].
In this paper the coupling between the CCS and a wave impact is investigated on the subpanel level. The subpanel is represented by the beam-foundation model proposed in [4], which predicts the response to the first vibration mode of the CCS. ELP 2: ‘building jet’ is chosen as load case, because it is expected to give a response on the length and time scale as the structure. This paper focusses on how the load is modelled and how the load interacts with the structure. It presents a simplified model for ELP2, in the form of non hydrostatic shallow water equations (NHSWE [14]), made suitable for numerical implementation of the coupling. To our knowledge this is the first paper to propose a reduced order model for fluid-structure interaction for wave impacts, that is fast and versatile enough for probabilistic evaluation. Figure 2 shows how the reduced order model from this paper relates to an actual wave impact. The drawn lines approximate the shape of the breaking wave, symmetrical over the horizontal plane. We developed this model to be faster than a high fidelity CFD code, so that it can be used in to perform a sensitivity study of the importance of FSI for various wave impact parameters.

Wave impact from [17]. The orange lines represent the cut-off wedge impact which is simulated in this paper. The cut-off wedge has an initial constant velocity
In the rest of the paper we will derive the equations for the fluid model and couple it to the structure. Then a number of wave impact scenarios are investigated with this code and a more detailed in-house CFD code: EVA. EVA is based on the Navier–Stokes equations for the motion of a mixed fluid with varying properties to model the combination of a liquid and a gas phase. The robust method makes use of a fixed staggered Cartesian grid where for the convective term a second-order upwind scheme is used. The time integration is based on a second-order Adams–Bashforth scheme. The description of the free surface is done with the Volume-of-Fluid (VoF) where the interface is geometrically reconstructed using piecewise-linear line segments. Free surface waves can be generated at boundaries, while reflected waves from the inside of the domain can be absorbed at the same time [26]. By consistently coupling the mass and momentum fluxes, EVA is capable to model high density ratio two-phase flows. The externally added interaction is solved using the Crank Nicolson scheme in the reduced order model and in EVA. A detailed description and validation of the CFD flow solver EVA is given in [23]. Both codes are used in this paper to investigate the effects of fluid-structure interaction for a wave impact.
The starting point is the Wagner impact model, which assumes an initial shape and constant velocity. It then calculates the flow and wall pressures using some perturbation method [19]. However, such a solution may not be feasible for the coupled problem, in which the relative velocity between structure and fluid is not constant. By making the following assumptions we have a similar starting conditions but obtain a model which is easier to solve. The following is assumed:
The wave impacts are perpendicular to the wall, and are symmetric with respect to the axis normal to the wall.
Water height and velocity are prescribed along the normal of the wall, and are only a function of the distance to the wall.
Gravity and air are disregarded in the solution, but kept in the derivation to be in line with literature
With these assumptions we obtain the non-hydrostatic shallow water equations [14,22,27]. In these equations part of the pressure is caused by the vertical acceleration of the fluid. Therefore the assumption of Wagner, disregarding gravity for an impact, can be used, and comparison is possible.
Fluid equations of motion
A similar approach is taken to derive the fluid equations as in [14], only here a flat bottom is assumed representing the symmetry of the impact. The starting point of the Fluid model are the Euler equations:
We can now average the continuity equation over the height:
To summarize, we have the following equations that describe the fluid behavior:
Implementation
The fluid model is discretized with a fractional step method. The first step, from
The second step is to proceed from
Equations (22)–(24) can be combined to yield an equation for
The equations are solved using collocated vertex centered finite differences. Discretization is upwind for convective terms, and central for the terms in Equation (26). The vertices are spaced evenly, dividing the
Structure equations of motion
The CSS is simplified to a mass-spring system following [4], giving an equation of the following form:
Coupled solution
The structure is placed at the edge of the fluid domain
The actual coupling takes place in the pressure correction step of the fluid. The structure is added to the system of equations of the fluid, thus solving both at the same time. This avoids sub-iterations or the choice of under relaxation parameters required for a partitioned solver and gives a robust implementation. The resulting system of equations is directly solved without iterations.
Uncoupled wedge impact
Following the initial investigation of the governing equations, a number of important parameters can be deduced. First, the rise time of the load. At the start the change from wedge to rectangle is at
To better compare the results we decide to present them in a dimensionless form, when applicable. The dimensionless time τ and force ϕ are defined as:
Due to the assumptions made in modelling the fluid flow there is no straightforward validation problem. However, we believe that the comparison between the reduced order model and EVA give enough confidence to trust this method for optimization or probabilistic evaluations. Because EVA is validated for fluid flow it is possible to make the link with reality. Additionally, this model is meant to predict the relevant trends, and the effect of taking hydro-elasticity into account. The next step could be to extend the model in a way that allows direct comparison with experiments.

Dimensionless plot of grid convergence of force on the wall for the two impact angles for EVA.
A convergence study for EVA for the investigated wedge angles is shown in Fig. 3. This figure shows on the horizontal axis dimensionless time and vertical axis the dimensionless total force on the wall. The dimensionless force on the wall is shown to converge to a constant value. The irregularity in the force just after impact, before the force starts to converge, is due to the dynamics of the liquid. Each time a cell originally filled with gas becomes filled with liquid it leads to a label change of that cell and some irregularity in the pressure. Such a label change is also visible near dimensionless time 2.2 for one grid and at time 2.4 for another (

Top row shows the free surface over time, as well as dimensionless force over time. The matrix below gives velocity over height for impact with

Top row shows the free surface over time, as well as dimensionless force over time. The matrix below gives velocity over height for impact with
Figures 4 and 5 show the force on the wall for an impact with
Both codes, NHSWE and EVA, adhere quite well to the expected rate at which the force builds up. At some point the slope decreases and the maximum force is reached and this depends highly on the angle α as mentioned before. Finally, all scaled impact forces converge to unity as expected. If the actual Wagner model was used, the load would keep on rising indefinitely, as in that case there is no end to the wedge. These tests have been repeated using different impact velocities but yield the same curve, meaning our aforementioned scaling factors work. It is clear that EVA shows a much higher peak value, up to 40%. Even disregarding the spikes in the force of EVA, likely caused by the fluid going into the next cell, the maximum force is about 30% higher. However, it should be noted that the difference in total impulse is within 3% percent because the same mass of fluid has to be turned around. In order to explain why there is a difference in maximum force, we need to look at the velocity profiles.
We know that there is one big difference between the models: EVA has multiple cells over the height of the domain, whereas NHSWE has only one. This allows EVA to have more local details, which we know the wave impact problem has [7,25]. The different velocity fields at four places are plotted in Figs 4 and 5, for four times (before max, max NHSWE, max EVA and near steady). These positions and times are indicated in the free surface and force over time plots in these Figures. In the second row of Figs 4 and 5 we can observe why the peak results are so different. The two rightmost graphs show that the horizontal velocity in EVA is higher where it is not obstructed by other parts of the fluid which are in contact with the wall. In EVA we can see the bottom part of the fluid column being stopped by the wall while the top is still progressing. At the same time, the NHSWE model requires the complete height of the column to be stopped at the same time. This requires compared to EVA a higher force from the start, but as total momentum is conserved, this leads to an earlier decay.
When most of the fluid stream is bent off in the bottom row, the initial wedge shape is no longer visible and both codes yield almost the same force. The velocity profiles become similar again, except for some local details.
Summarizing, the upward and downward slopes seem the same, the response is of the same order of magnitude and the impulse is the same because the same amount of fluid is redirected. The difference then comes from the higher resolution in EVA, which models more accurate the velocity profiles due to the sharp bent in the flow, see Fig. 4, which is not as well modelled by NHSWE as by EVA. A downside of EVA are the pressure spikes that increase the absolute maximum. However, as the total impulse is nearly the same, as is the velocity profile at the end of the simulation, we believe it is worthwhile to investigate the coupled response of these wave impact codes with a moving wall.
Overview of the simulations performed with EVA and the simplified NHSWE model, all with
We choose to investigate a subpanel of the Mark III CCS [10,13], with properties derived as described in [4] (
Comparison with EVA
Figure 6 shows the load and response of the CCS to a wave impact with velocity

Dimensionless wall force (a) and response (b) for uncoupled and coupled impacts at 8 and 16 m/s with 60° and 70°.
The response of the wall is shown in Fig. 6b and the response becomes increasingly dynamic with increasing angle α or velocity
We also see that the response of the structure after the impact does not damp out. Over time mass enters the domain with a velocity, therefore adding energy. We expect that the coupled system response will not die out, so a sensible cut-off time should be taken. A simple proposal is the time when the water column in front of the structure exceeds the structure size, or when the uncoupled wall force stabilizes (within a margin) to the expected steady wall force of Equation (32). Here, we took a dimensionless time of
For the 8 m/s impact velocity the impact forces almost coincide for coupled and uncoupled, just the load case with a 70° wedge angle shows a deviation for all time. The 8 m/s, 60° impact is just slightly decreased by the moving structure. With the fast 16 m/s impact a much larger difference is observed, clearly showing the vibration of the structure in the load. We can see that the coupled solution can result in a higher load on the structure than the uncoupled solution. There is a simple explanation: the structure moves along with the fluid but at some point has to come back. When it does, the relative velocity of the fluid with respect to the structure is increased, therefore increasing the load. This all happens at a time scale where the dynamics of the structure is important. In other words: whether the load is increased due to FSI is very dependent on the dynamics of the load and structure.
Finally, it should be noted that Table 1 gives the duration of the final vibration period. We can see a big difference between the coupled and uncoupled structures, where the coupled response is almost the same for both solvers. Physically this means the added mass for both solvers is approximately the same. Furthermore, even though the impact maximum is equal, the NHSWE solver does say something about the effect of FSI and it could be a valuable tool to check the system sensitivity with respect to the coupling.
From Table 1 we can also see the height of the last peak for coupled or uncoupled response is the same. It means that we do predict approximately correctly how much energy is absorbed by the structure, and whether the FSI is amplifying the response or diminishing it. The added benefit of the reduced order model is that it is much faster than the high fidelity method in doing so, in the order of a few minutes versus a day of computation time.

Contour plot of the ratio between coupled and uncoupled force on the wall and displacement of the wall. If the value is larger than unity the force or displacement in the coupled case is higher than for the uncoupled case. Unity is marked with the white contour lines. The black dashed lines are the critical times from equation (30). The initial height is kept constant at
The low fidelity model can be used to check the effect of different angles and impact velocities, as an example of a sensitivity study. Figure 7 shows a study of the difference between coupled and uncoupled impacts for angles between
It is clear from the both graphs that the lines of the ratio have similar alignment as the lines of the critical time. Especially the critical time
It is clear from the left plot in Fig. 7 that the maximum force of the coupled calculation is in almost all cases higher than for the uncoupled calculation. As mentioned before, this is attributed to the structure first moving along with the impact, but after some time coming back, that is when the maximum force is reached. According to this figure this mostly happens when the critical time is between the wet and dry natural period.
The right graph is much harder to generalize. It looks like higher angles or velocities might show more that the FSI coupling decreases the response. But, for the other combinations of angle and velocity there is no clear trend, which means that whether the structure is sensitive to FSI depends on the specific impact parameters.
Wave impacts onto the wall of an LNG tank are limiting for the design. The interaction between the wave and the wall could increase or decrease the response of the structure. High fidelity models for fluid-structure interaction are typically computationally expensive, which is why we propose a reduced order model for the fluid. The reduced order model is coupled to a reduced order model for the structure, being an order of magnitude faster than the high fidelity model. This reduced order model makes it possible to estimate quickly the effect of hydro-elasticity on the dynamics of the system during a wave impact.
We consider as impact case a cut off wedge, as a model for the wave crest hitting the wall. The reduced order model predicts the same total impulse and added mass as the calculations with EVA, which is an in-house developed CFD code. It shows lower force than EVA, yet it is usable to quantify the sensitivity of the structure response to fluid-structure interaction. This is possible because the total impulse and added mass are almost the same, leading to almost the same final excitation. The difference can be explained by the kinematic assumption used to arrive at the reduced order model and its influence on the maximum force.
An example of a sensitivity study is given, where impact velocity and angle are varied for the same structure. From these results it seems like the force in the coupled calculations is generally higher than for the uncoupled calculations. However, there is not such a clear trend for the displacement of the wall. It was not possible to derive a rule of thumb for the importance of FSI in wave impacts, as it depends strongly on the impact shape and velocity.
For future studies it is recommended to generalize the sensitivity analysis in the sense of dimensionless numbers, for instance the ratio of critical time to vibration period. Maybe in such cases it is easier to draw general conclusions, and look for areas where FSI is never important. The reduced order model can be improved adding functions that describe the velocity field in the impacting wave. Therefore the wave impact forces should be modelled more accurately, as the local flow details are better captured. However, this does come at added complexity and effectively brings the reduced order model closer to EVA. Another improvement would be to stretch the domain with the moving wall, improving accuracy for larger wall deformations.
The aforementioned screening method could calculate the response to a wave impact using the simplified models presented in this paper. Then, the most critical load cases can be further investigated using a high fidelity code. Note that in this paper only a cut-off wedge is used for the impact, but the shape can be much more general.
Funding
This work is part of the research programme ‘SLING’ with project number P14-10 which is (partly) financed by the Netherlands Organisation for Scientific Research (NWO). 
Conflict of Interest
M. van der Eijk and P.R. Wellens are Editorial Board Members of this journal, but were not involved in the peer-review process nor had access to any information regarding its peer-review.
