CFD Based Non-Dimensional Characterization of Energy Dissipation Due to Verticle Slosh

Abstract

We present the CFD based non-dimensional characterization of violent slosh induced energy dissipation due a tank under vertical excitation. Experimentally validated CFD is used for this purpose as an ideally suited and versatile tool. It is thus first demonstrated that a weakly compressible VoF based CFD scheme is capable of computing violent slosh induced energy dissipation with high accuracy. The resulting CFD based energy analysis further informs that the main source of energy dissipation during violent slosh is due liquid impact. Next, a functional relationship characterising slosh induced energy dissipation is formulated in terms of fluid physics based non-dimensional numbers. These comprised contact angle and liquid–gas density ratio as well as Reynolds, Weber and Froude numbers. The Froude number is found the most significant in characterising verticle violent slosh induced energy dissipation (in the absence of significant phase change). The validated CFD is consequently employed to develop scaling laws (curve fits) which quantify energy dissipation as a function of the most important fluid physics non-dimensional numbers. These newly developed scaling laws show for the first time that slosh induced energy dissipation may be expressed as a quadratic function of Froude number and as a linear function of liquid–gas density ratio. Based on the aforementioned it is postulated that violent slosh induced energy dissipation may be expressed as a linear function of tank kinetic energy. The article is concluded by demonstrating the practical use of the novel CFD derived non-dimensional scaling laws to infer slosh induced energy dissipation for ideal experiments (with exact fluid physics similarity to the full scale Aircraft) from (non-ideal) slosh experiments.

1. Introduction

The modelling of liquid–gas free-surface interaction (sloshing) within tanks and containers is an area of interest within the aerospace industry, impacting aerodynamic stability and aircraft control characteristics during flight, ground manoeuvres and loads. As such the EU H2020 SLOWD (SLOshing Wing Dynamics) [1] aims to characterise and model the impact of fuel sloshing on the damping characteristics of a wing structure due to vertical excitation where the direction of excitation is normal to the liquid free-surface. The preceding Protospace tests [2,3] proved that liquid sloshing provides significant energy dissipation (damping) to a cantilever based tank undergoing free-vibration due to an initial structural perturbation. In this case, the major frequency of excitation was found to remain constant. Note that the beam tip underwent geometrically linear motion (The verticle displacement is small compared to the beam length.). This article, therefore, focuses on non-dimensionally characterising the violent liquid slosh induced energy dissipation due vertical tank motion at a fixed frequency $\omega$. The tank motion is therefore prescribed as a single-degree-of-freedom (SDOF) system. Further, it is assumed that the liquid density is orders of magnitude larger than that of the gas while the liquid may be assumed as incompressible while the gas is modelled via a weakly compressible formulation. Finally, the flow physics does not involve significant phase change.

Due to the laboratory scale of the cited experiments [2,3], increased understanding of the bias introduced by such geometric scaling is a key objective of this article. To this end we develop suitable similarity parameters to account for the effects introduced by non-ideal (practical) scaling such that the influence of sloshing on the structural dynamics of full scale aircraft wings may be quantified. Once validated, the use of CFD is key in such studies, as it allows for quantification of the influence of a wide range of parameters and their combinations on slosh induced damping. This is not possible experimentally if limited to available fluid properties. In this work, all CFD calculations are conducted with the multi-physics CFD code Elemental^® [4,5,6,7,8].

It is widely accepted that non-dimensional analysis provides a powerful technique to characterise a complex physical system. This enables the study and quantification of the relative importance of various physical parameters/phenomena on such systems [9]. Despite the insights that non-dimensional analyses provide, limited studies of violent slosh through such techniques exist in literature. Thompson and Nein [10] successfully employed dimensional analysis as a method to predict the pressurisation requirements for space vehicles at launch. Summer [11] developed a generalisation of slosh force subject to several excitation frequency parameters for various tank geometries. The above cited work did not deal specifically with slosh induced energy dissipation.

Past work into the application and tuning of liquid damping systems has centred on structural design and civil engineering applications [12,13]. Previous aeroelastic studies into the influence of slosh damping focus on the general trends and changes in flutter boundaries due to various tank fill-level configurations [14]. These campaigns neglected the vertical motion of the tank, employing seismic shaking tables to provide the excitation to the structure.

This article furthers on existing work in a number of respects. First, weakly compressible volume of fluid (VoF) based CFD [15,16] is experimentally validated, predicting violent slosh induced energy dissipation to with-in $2\%$ of the experimental values. The CFD also uncovers that the main source of slosh induced energy dissipation is due liquid impact. Next, the SDOF system is non-dimensionalised such that slosh induced energy dissipation may for the first time be expressed accurately in terms of Froude number and liquid-gas ratio via CFD generated curve fits. These in turn for the first time point to the possibility that slosh induced energy dissipation may be expressed as a linear function of tank excitation kinetic energy. Finally, it is demonstrated how the newly fitted non-dimensional scaling laws may be used with ease to estimate the slosh induced damping of ideally scaled experiments. Note that the development of the aforementioned scaling laws is impossible via laboratory size experimental campaigns due to limitations in available fluids. As such, for this purpose we employ the rigorously validated Elemental^® volume-of-fluid (VoF) CFD software.

This paper is organised as follows:

• We validate the accuracy of the CFD code Elemental^® for modelling violent vertical slosh physics relevant to this article.
• We detail the CFD based energy budget used to quantify effects of scaling non-dimensional properties on the system.
• We present a non-dimensional analysis of a SDOF tank under vertical slosh, isolating the functional relationship characterising slosh induced energy dissipation.
• We define the non-dimensional parameter space of interest for the problem under consideration (to include both experimental and full scales).
• Finally, we develop novel scaling-laws which correlate the slosh induced energy dissipation as a function of the identified non-dimensional parameters. This is done via curve fitting of CFD generated data.
• The developed novel scaling laws are finally applied to quantify ideal (representative of full scale aircraft) experimental slosh induced energy dissipation.

2. Validation of CFD Model

The utilisation of VoF CFD within the study requires accurate, high-fidelity, mesh independent simulations. Specifically, the accuracy with regards to the predicted slosh force is critical as this is the key uncertainty when quantifying slosh induced energy dissipation in this work (see Equation (12)). This is due to this work involving prescribed tank motions (obtained from measurements). The VoF CFD code Elemental^® has been employed extensively to model violent slosh [17,18]. Turbulence is modelled via Large-eddy simulation (LES) employing the Smagorinsky-Lilly model [19] and a weakly compressible gas model [15,16] is employed while the liquid-gas interface initialisation is performed to machine precision accuracy via the AGI tool [20].

The validation of Elemental^® as an accurate tool for this work requires the comparison of the computed results against relevant experimental data. The modelling of the $50\%$ fill-level Protospace experiment [3] has been selected for this purpose as it was designed to be a small scale representation of wing-tank based slosh. The Protospace experiment involved the damped vibration (structural and slosh damping) of a cantilever beam mounted with a seven compartmental fuel tank, as shown in Figure 1. The experiment involved applying an initial deflection and then releasing the cantilever at $t=0s$ followed by free vibration [3].

Using the measured tank accelerations from the experiment as input, CFD simulations were performed on a range of different mesh spacings. For the coarsest mesh considered ($\Delta x={10}^{-3}$), 2D and 3D simulations resulted in near identical computed slosh forces and therefore energy dissipation due to the definition of $E_{Disp}$ (see Equation (12)). In the interest of computational efficiency therefore, 2D was employed for all subsequent simulations. The 2D CFD predicted free-surface interface at various times is depicted in Figure 2. These may be compared to experimental footage (Figure 1). As shown, a range of different slosh regimes are present. These include liquid-solid impact (Figure 1b and Figure 2b), aerated flow Figure 1c and Figure 2c) and weakly linear slosh (Figure 1d and Figure 2d). Further, a reasonable agreement exists between predicted and experimental footage.

The use of LES to model turbulence required estimation of the Taylor micro-scale $\left(\lambda \right)$ [21] as this informs mesh resolution. This was done as part of the SLOWD project [22] which yielded

$$\lambda \approx \sqrt{\frac{10\nu k}{ϵ}}=4\times {10}^{-4}m$$

(1)

where $\nu$ is the fluid kinematic viscosity, $ϵ$ is kinetic energy dissipation rate and k is the turbulent kinetic energy. The above aided in informing suitable CFD mesh spacing as detailed next.

Figure 3 compares the measured and CFD predicted vertical slosh loads (See Equation (11)) for Cartesian meshes, namely $\Delta x=1\times {10}^{-3}$ and $\Delta x=5\times {10}^{-4}$. Completing an error analysis and mesh sensitivity study provides

Table 1. Error analysis.

$\Delta x$	${\left|\right|ϵ\left|\right|}_2$ (%)
$1\times {10}^{-3}$	$1.22$
$6.67\times {10}^{-4}$	$1.212$
$5\times {10}^{-4}$	$1.210$
$4\times {10}^{-4}$	$1.209$

and Figure 4. This highlights engineering precision convergence for $\Delta x\le 6.67\times {10}^{-4}$ as the computed dissipated energy ($E_{Disp}$ in Figure 5 computed via Equation (12)) is found to lie within $2\%$ of the experimental result. This validates the use of Elemental to predict violent slosh induced slosh loads and damping, provided that suitably fine meshes are employed. To ensure direct applicability, the SDOF slosh model used for our non-dimensional characterization will consist of one of the Protospace compartments and be modelled with a 2D Cartesian mesh with spacing $\Delta x=5\times {10}^{-4}$.

3. Energy Analysis

We next derive the relations employed to estimate slosh induced energy dissipation. For the system under consideration the mechanical energy $\left(E_M\right)$ budget is considered for this purpose. This is composed of the work done on the fluid by the domain boundary surfaces or surface work $\left(W_s\right)$, deformation work $\left(E_{Defo}\right)$, kinetic energy $\left(E_{Kine}\right)$, potential energy $\left(E_{Pote}\right)$, viscous dissipation $\left(E_{visc}\right)$, and surface tension energy $\left(E_{ST}\right)$. The full energy budget is defined by

$$\begin{array}{cc}\hfill E_M=W_s=& E_{Defo}+E_{Kine}+E_{Pote}+E_{visc}+E_{ST}\hfill \\ \hfill =& E_{Disp}+E_{Kine}+E_{Pote}+E_{ST}\hfill \end{array}$$

(2)

where $E_{Disp}=E_{Defo}+E_{visc}$ and $E_{Defo}$ includes losses due liquid impact (as described next).

Consider an enclosed fluid domain containing an immiscible liquid–gas mixture (gas is weakly compressible) undergoing a rigid body linear acceleration of ${\mathbf{a}}_d$. The resulting absolute fluid velocity ${\mathbf{u}}_a$ (velocity relative to earth) can be defined as the sum of the domain velocity $\left({\mathbf{u}}_d\right)$ and the fluid velocity $\left(\mathbf{u}\right)$ relative to the non-inertial reference frame (CFD mesh reference frame), therefore ${\mathbf{u}}_a=\mathbf{u}+{\mathbf{u}}_d$. The momentum conservation of the fluid is given by:

$$\begin{array}{cc}\hfill \frac{D}{Dt}\left(\rho {\mathit{u}}_a\right)& =\left[\frac{\partial \rho \mathit{u}}{\partial t}+\mathbf{\nabla }·(\rho \mathit{u}\otimes \mathit{u})+\rho {\mathit{a}}_d\right]\hfill \\ \hfill & =-\mathbf{\nabla }p+\mathbf{\nabla }·\mu \left(\mathbf{\nabla }\mathit{u}+\mathbf{\nabla }{\mathit{u}}^T\right)+\rho \mathit{g}+\gamma \kappa \mathbf{\nabla }\alpha \hfill \end{array}$$

(3)

where p is pressure, t is time, $\alpha$ is VoF volume fraction and the fluid properties are defined in

Table 2. Dimensional Matrix.

Quantity			Units
			M	L	T
$\mu$	Liquid Viscosity	[Pa·s]	1	−1	−1
$\rho$	Liquid Density	[kg·m${}^{-3}$]	1	−3	0
${\rho }_{air}$	Gas Density	[kg·m${}^{-3}$]	1	−3	0
$\gamma$	Liquid Surface Tension	[N·m${}^{-1}$]	1	0	−2
$\theta$	liquid–solid Contact Angle	[rad]	0	0	0
v	Fluid speed	[m·s${}^{-1}$]	0	1	−1
g	Gravity	[m·s${}^{-2}$]	0	1	−2
m	Solid Mass	[kg]	1	0	0
k	Structural Stiffness	[N·m${}^{-1}$]	1	0	−2
c	Structural Damping	[Ns·m${}^{-1}$]	1	0	−1
h	Tank Height	[m]	0	1	0
l	Tank Length	[m]	0	1	0
${\eta }_0$	Height of fluid	[m]	0	1	0
$y_0$	Initial tank offset (spring displacement)	[m]	0	1	0

. Mechanical energy conservation results from taking the dot product of Equation (3) with the volume averaged ${\overline{\mathbf{u}}}_a$ followed by integration over space and time (includes integration by parts of the pressure term):

$$\begin{array}{cc}\hfill & {\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\left[\frac{\partial \rho \mathit{u}}{\partial t}+\mathbf{\nabla }·(\rho \mathit{u}\otimes \mathit{u})+\rho {\mathit{a}}_d\right]dVdt=\hfill \\ \hfill & -{\int }_0^t{\mathit{u}}_d·{\int }_{{\mathcal{A}}_d}p\mathit{n}dAdt\hfill \\ \hfill & +{\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\left[\mathbf{\nabla }·\mu \left(\mathbf{\nabla }\mathit{u}+\mathbf{\nabla }{\mathit{u}}^T\right)\right]dVdt\hfill \\ \hfill & +{\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\rho \mathit{g}dVdt+{\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\gamma \kappa \mathbf{\nabla }\alpha dVdt\hfill \end{array}$$

(4)

where volume averaged quantities are defined as $\overline{\varphi }=\frac{1}{\mathcal{V}}{\int }_{\mathcal{V}}\varphi d\mathcal{V}$, $\mathcal{A}$ is the surface of all mesh finite volumes $\mathcal{V}$, and $\mathit{n}$ is the outwards pointing normal. Further ${\mathcal{A}}_d$ is the fluid domain bounding surface and $\nabla ·{\overline{\mathbf{u}}}_a=0$ due the domain being assumed rigid. The terms in Equation (4) quantify the different forms of energy:

$$\begin{array}{cc}\hfill W_s& =-{\int }_0^t{\mathit{u}}_d·{\int }_{{\mathcal{A}}_d}\left[p-\mu \left(\mathbf{\nabla }\mathit{u}+\mathbf{\nabla }{\mathit{u}}^T\right)\right]\mathbf{n}dAdt={\int }_0^t{\mathit{u}}_d·\mathbf{F}dt\hfill \end{array}$$

(5)

$$\begin{array}{cc}\hfill E_{Pote}& =-{\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\rho \mathit{g}dVdt\hfill \end{array}$$

(6)

$$\begin{array}{cc}\hfill E_{Kine}& ={\int }_{\mathcal{V}}\frac{1}{2}\rho {\overline{\mathit{u}}}_a·{\overline{\mathit{u}}}_{a}dV\hfill \end{array}$$

(7)

$$\begin{array}{cc}\hfill E_{Defo}& ={\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\left[\frac{\partial \rho \mathit{u}}{\partial t}+\mathbf{\nabla }·(\rho \mathit{u}\otimes \mathit{u})+\rho {\mathit{a}}_d\right]dVdt-E_{Kine}\hfill \end{array}$$

(8)

$$\begin{array}{cc}\hfill E_{visc}& =-{\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\mathbf{\nabla }·\mu \left(\mathbf{\nabla }\mathit{u}+\mathbf{\nabla }{\mathit{u}}^T\right)dVdt+W_{s_{\mu }}\hfill \end{array}$$

(9)

$$\begin{array}{cc}\hfill E_{ST}& =-{\int }_0^t{\overline{\mathit{u}}}_a·{\int }_{\mathcal{V}}\gamma \kappa \mathbf{\nabla }\alpha dVdt\hfill \end{array}$$

(10)

where $\mathbf{F}$ is the force from the tank onto the liquid and $W_{s_{\mu }}$ is the domain boundary viscous component included in (5). Therefore, the force is computed via the following expression,

$$F_{Slosh}=-{\int }_{{\mathcal{A}}_d}\left[p-\mu \left(\mathbf{\nabla }\mathit{u}+\mathbf{\nabla }{\mathit{u}}^T\right)\right]\mathbf{n}dA$$

(11)

which is used for the CFD load calculations.

The reason for subtracting the kinetic energy term from the first integral in (8) is to account for the liquid impact related losses (we assume initial fluid velocities of ${\mathit{u}}_a^0=0$). These impact losses (liquid-liquid or liquid-wall) have been shown [23,24] to be accounted for by the liquid incompressible pressure-poisson fractional step [17] method employed in our solver. This is consistent with the accuracy achieved in the preceding section in computing slosh induced energy dissipation. Note that liquid impact was found to account for the dominant source of energy dissipation in this work. From the above, it is instructive to note that $E_{Disp}$ may also be approximated from experimentally measured data as follows:

$$\begin{array}{cc}\hfill E_{Disp}=& {\int }_0^t{\mathit{u}}_d·\mathbf{F}dt-E_{Kine}-E_{Pote}\hfill \\ \hfill \approx & {\int }_0^t{\mathit{u}}_d·\mathbf{F}dt-\frac{1}{2}{m}_f{\mathbf{u}}_d·{\mathit{u}}_d-{\int }_o^t{m}_f{\mathit{u}}_d·\mathit{g}dt\hfill \end{array}$$

(12)

where $m_f$ is the fluid mass in the tank and $\mathbf{F}$ is the force of the tank on the fluid (which may be inferred from load cell measurements). This is provided that accurate tank position (displacement) data are available and noting it was found that $E_{ST}$ accounts for less than $1\%$ of the energy budget.

4. Dimensional Analysis

The remainder of the article will consider a slosh system which is representative dimensionally of the Protospace experimental setup. As such we consider the vertical oscillation (SDOF) of a rectangular 2D tank (single Protospace compartment) supported by a mass-spring-damper system with a spring constant k and damper constant c. All simulations employ air at standard atmospheric conditions as the gas within the liquid-gas slosh system. Therefore the 14 parameters considered relevant to the dimensional analysis are presented in Table 2. Also included are entries of the dimensional matrix for the units of Mass (M), Length (L) and Time (T).

The rank 3 dimensional matrix employs the normalisation of length $\left(h\right)$, mass $\left(\rho h^3\right)$ and time $\left(\frac{h}{v}\right)$ respectively. A straightforward application of the Buckingham Pi Theorem leads to $14-3=11$ dimensionless groups:

$$\begin{array}{cc}\hfill {\pi }_1& =\frac{m}{\rho h^3}=\overline{m}\hfill \\ \hfill {\pi }_4& =\frac{v}{\sqrt{gh}}=F{r}_v\hfill \\ \hfill {\pi }_7& =\frac{c}{\rho h^3}\frac{h}{v}\hfill \\ \hfill {\pi }_{10}& =\frac{{\eta }_0}{h}=\tilde{F}\hfill \end{array}\begin{array}{cc}\hfill {\pi }_2& =\frac{{\rho }_{Air}}{\rho }=\overline{\rho }\hfill \\ \hfill {\pi }_5& =\frac{\rho v^{2}h}{\gamma }=W{e}_v\hfill \\ \hfill {\pi }_8& =\theta \hfill \\ \hfill {\pi }_{11}& =\frac{y_0}{h}={\tilde{Y}}_0\hfill \end{array}\begin{array}{cc}\hfill {\pi }_3& =\frac{\rho vh}{\mu }=R{e}_v\hfill \\ \hfill {\pi }_6& =\frac{k}{\rho h^3}{\left(\frac{h}{v}\right)}^2\hfill \\ \hfill {\pi }_9& =\frac{l}{h}=\tilde{A}\hfill \\ \hfill \end{array}$$

(13)

We define the ${\pi }_1$ group as of the order solid to liquid mass ratio $\overline{m}$ of the one-degree of freedom system. ${\pi }_2$ is the density ratio $\overline{\rho }$ of the gas and liquid phase in the tank. $\pi$ groups 3 to 5 are well known flow physics related dimensionless numbers, namely Reynolds, Froude and Weber. Note that the v subscript is added to the latter to indicate that the standard definition based on the fluid speed was employed.

As we are interested in the energy dissipation of a tank which is supported by a spring-mass-damper system, the fluid non-dimensional numbers are to be written in terms of the resulting excitation. Considering harmonic motion at frequency $\omega$ enables expressing the reference flow velocity in terms of the slosh-mass-spring-damper system response as:

$$v\sim \omega h$$

(14)

This is as in this work the magnitude of tank displacement is of a similar order to the tank height. To consolidate the validity of Equation (14) we note that from experimental observation the first two to three cycles of oscillation are representative of the slosh induced energy dissipation potential for the violent slosh under consideration. Hence only the first 3 excitation cycles are considered for the non-dimensional scaling calculations (to follow) as these are representative of peak loads and slosh induced dissipation. We can therefore modify $\pi$ groups 3 to 5 as follows

$$\begin{array}{cc}\hfill {\pi }_3^{\prime }& =\frac{\rho \omega h^2}{\mu }=Re\hfill \end{array}\begin{array}{cc}\hfill {\pi }_4^{\prime }& =\omega \sqrt{\frac{h}{g}}=Fr\hfill \end{array}\begin{array}{cc}\hfill {\pi }_5^{\prime }& =\frac{\rho {\omega }^2{h}^3}{\gamma }=We\hfill \end{array}$$

(15)

where $Re$, $Fr$ and $We$ are the excitation frequency-based Reynolds, Froude and Weber numbers, respectively. The Froude number as a function of the frequency parameter, therefore, infers a reduced frequency or a Strouhal number which represents the ratio of the characteristic time of the mechanical system and the fluid-characteristic time. Finally, the manipulation of ${\pi }_6$ and ${\pi }_7$ by combining with ${\pi }_1$ provides:

$$\begin{array}{cc}\hfill {\pi }_6^{\prime }& =\sqrt{\frac{{\pi }_6}{{\pi }_1}}=\sqrt{\frac{k}{m}}\frac{h}{v}=\frac{{\omega }_n}{\omega }={\omega }_r\hfill \\ \hfill {\pi }_7^{\prime }& =\frac{{\pi }_7}{2{\pi }_1{\pi }_6^{\prime }}=\frac{c}{2m{\omega }_n}={\xi }_d\hfill \end{array}$$

(16)

Here, we employ the standard definitions of the structural natural frequency of excitation ${\omega }_n$ and the structural damping ratio ${\xi }_d$ for harmonic oscillators, where:

$${\omega }_n=\sqrt{{\displaystyle \frac{k}{m}}}$$

(17)

and

$${\xi }_d={\displaystyle \frac{c}{2m{\omega }_n}}$$

(18)

Further, ${\omega }_r$ denotes the reduced frequency of the system (structure and slosh).

We are justified to choose h (tank height) as the characteristic length, as ${\pi }_9$, ${\pi }_{10}$ and ${\pi }_{11}$ are derived by dividing all other relevant geometrical quantities by the characteristic length. The resulting $\tilde{A}$ is the tank aspect ratio, $\tilde{F}$ the liquid fill ratio and $\widetilde{Y_0}$ the dimensionless initial deflection of the mass-spring system. For the slosh case considered in this article, these quantities are $O\left(1\right)$ and as a result the liquid centre of gravity displacements and wave height. The latter is a consequence of the slosh being within the physical bounds of the tank.

The derivation of a functional relationship to characterise the fluid-structure-system-damping will employ the resulting $\pi$ groups:

$$\xi =f\left(\overline{m},{\omega }_r,\overline{\rho },Re,Fr,We,\theta ,\tilde{A},\tilde{F},{\tilde{Y}}_0\right)$$

(19)

The total system damping $\xi$ can be decomposed into its two parts:

$$\xi ={\xi }_d+\Delta {\xi }_s$$

(20)

where ${\xi }_d$ and $\Delta {\xi }_s$ are the damping due the structure (dry) and slosh damping. To focus the analysis on the slosh induced damping via prescribed tank motion, we shall prescribe $\xi$ (the chosen value will be motivated when used in the relevant sections to follow) and such that the derivation of a functional relationship characterising the slosh damping results in the form:

$$\Delta {\xi }_s=F\left(\overline{m},{\omega }_r,\overline{\rho },Re,Fr,We,\theta ,\tilde{A},\tilde{F},{\tilde{Y}}_0,\xi \right)$$

(21)

From the above, we can identify three distinct groups of dimensionless numbers influencing the slosh damping within the SDOF system. The three groups are flow physics (fluid specific non-dimensional numbers such as $\overline{\rho }$, $Re$, $Fr$, $We$ and $\theta$), fluid-structure interactions (FSI) and geometric non-dimensional factors. This paper employs high fidelity VoF CFD to investigate the influence of a range of flow physics parameters (As defined within Section 6) while the geometric numbers will remain unchanged and the FSI group will be selected such as to be consistent with the flow physics numbers. Specifically, $\omega$ will be dictated by the $Fr$ number.

5. Scaling of Violent Sloshing Systems

To provide a perfect scaled experiment requires complete mechanical similarity with the full scale system. This implies geometric, kinetic and dynamic similarity (identical force ratios for all terms in Equation (3)). While it is possible to achieve geometric and kinetic similarity, ensuring complete dynamic similarity is not practically feasible. This is as it requires identical non-dimensional numbers for all flow physics related ratios ($\overline{\rho }$, $Re$ etc.) which is impossible using available fluids. As a result, experimental similarity of a single non-dimensional number is then targeted in industry, which in the case of slosh is the Froude number. Scaled slosh experiments are therefore designed by ensuring similar Froude numbers between experimental and real systems.

5.1. Froude Scaling Rules

For experimental scale, it is possible to achieve perfect geometric scaling $\lambda =h/h_{fs}$, where $h_{fs}$ is the tank height at full scale, and consequently that this scaling factor applies to all the geometric ratios. As noted above, slosh fluid physics is then commonly scaled using the Froude number:

$$Fr={\omega }_{fs}\sqrt{\frac{h_{fs}}{g}}=\omega \sqrt{\frac{h}{g}}$$

(22)

hence,

$$\frac{\omega }{{\omega }_{fs}}=\sqrt{\frac{h_{fs}}{h}}={\lambda }^{-\frac{1}{2}}$$

(23)

Consequently, the Froude scaling for time t is

$$\frac{t}{t_{fs}}={\lambda }^{\frac{1}{2}}$$

(24)

Ensuring geometrical and Froude similarity constrains two of the three fundamental units (length and time). The rule for mass is generally achieved by simply using a test fluid of the same density as the full scale (in fact when possible the same fluid is used). This assumption leads to

$$\frac{m}{m_{fs}}=\frac{\rho h^3}{{\rho }_{fs}{h}_{fs}^3}={\lambda }^3$$

(25)

All other scaling rules then follow from those derived above for length, time and mass such that:

5.2. Practical Considerations

As a direct consequence of

Table 3. Scaling Rules.

Quantity	Factor
Length	$\lambda$
Time	${\lambda }^{1/2}$
Mass	${\lambda }^3$
Velocity	${\lambda }^{1/2}$
Acceleration	${\lambda }^0=1$
Viscosity	${\lambda }^{1.5}$
Surface Tension	${\lambda }^2$
Force	${\lambda }^3$
Energy	${\lambda }^4$

, using the same fluid for a reduced-scale experiment leads to differences in the Reynolds and Weber numbers to the full scale. In addition, there are no practical constraints in enforcing equality for the other $\pi$ groups.

$$\begin{array}{cc}\hfill \frac{Re}{R{e}_{fs}}& ={\lambda }^{1.5}\ne 1\hfill \\ \hfill \frac{We}{W{e}_{fs}}& ={\lambda }^2\ne 1\hfill \end{array}$$

(26)

By way of example, for 1:5 scale ($\lambda =0.2$), Equation (26) imply viscous effects approximately ten times greater (${\lambda }^{1.5}\approx 0.1$) and surface tension effects approximately twenty-five times greater (${\lambda }^2\approx 0.04$) than those at full scale. An alternative could be to engineer a fluid with the same density as full scale and the appropriate scaling for viscosity and surface tension. The associated time and cost (if at all possible) however make the use of validated numerical models instead ideal.

In practice, even exact Froude scaling is not always feasible due to budget and time constraints. As shown in

Table 4. Protospace Froude Similarity Parameters.

	f [Hz]	$\mathit{F}\mathit{r}$
Required	3.35	1.65
Achieved	7.0	3.44

, the Protospace experiment did not achieve the required $Fr$ number. With the available materials, it was not possible to reach a closer similarity, as this would have required a bespoke cantilever design. The lower natural frequency of $3.35$ Hz would have implied a greater initial displacement of the free-end, possibly leading to a geometrically non-linear oscillation of the beam. Further, the fact that the characteristic velocity achieved is higher than that required by $Fr$ similarity could indicate that the dissipative effects due to liquid impact are overestimated experimentally if compared to an ideally Froude scaled experiment. All of the aforementioned point to the use of CFD to quantify the effect of the various non-dimensional numbers on all key physical characteristics (forces, damping, etc.) of interest. This will result in so-called non-dimensional scaling-laws, from which experimental data may be used to estimate the flow physics of the full-scale tank. This is the objective of the sections to follow and indeed a key contribution of this article.

6. Non-Dimensional Property Parameter Space

The intended aim of the Protospace experiment [3] was to provide a Froude scaled aircraft wing with a geometry scaling of $\lambda =0.2$ while employing water whereas the full-scale model makes use of cold kerosene. The calculated non-dimensional values for the Froude scaled Protospace (employing kerosene), actual Protospace (as achieved in the experimental campaign employing water) and the full-scale model are given in

Table 5. Protospace vs Full Aircraft Non-dimensional Properties.

	Actual	Froude Scaled	Aircraft
	Protospace	Protospace
$\theta$	$1.05$	0	0
$\overline{\rho }$	$1.21\times {10}^{-3}$	$1.39\times {10}^{-3}$	$1.39\times {10}^{-3}$
$Re$	$161.7\times {10}^3$	$4.67\times {10}^3$	$52.2\times {10}^3$
$We$	$5.73\times {10}^3$	$2.78\times {10}^3$	$69.4\times {10}^3$
$Fr$	$3.44$	$1.65$	$1.65$

. Accordingly, the range of the non-dimensional parameter space for which to build the scaling-laws was chosen to contain both the Protospace and full-sized model non-dimensional values, as presented in

Table 6. Non-dimensional property parameter space.

Parameter	Property Range
$\theta$$\in [0.52,1.51]$	−
$\overline{\rho }$$\in [6.49\times {10}^{-4},6.18\times {10}^{-3}]$	$\therefore \rho \in [199,1894]$ kg·m${}^{-3}$
$Re$$\in [42.56\times {10}^3,\infty )$	$\therefore \mu \in [0,3.7\times {10}^{-3}]$ Pa·s
$We$$\in [4.92\times {10}^3,1.37\times {10}^5]$	$\therefore \gamma \in [0.003,0.085]$ N·m${}^{-1}$
$Fr$$\in [0.69,5.16]$	$\therefore f\in [1.4,10.5]$ Hz

. Note that the liquid–solid contact angle lower limit allowed by the software was ${30}^{\circ }$ and thus does not contain the $0^{\circ }$ value for kerosene. It was however found that (see next section) the key importance of contact angle lies in the initial meniscus, i.e., if the contact angle is ${90}^{\circ }$ or not. This is as the $We$ numbers are high while a ${90}^{\circ }$ contact angle results in the absence of Rayleigh–Taylor instability and therefore no violent slosh. Further, turbulence was found to be a minor contributor to energy dissipation (the major being liquid impact) and hence $Re=\infty$ here implies the inviscid limit, i.e., $\mu =0$.

The CFD simulations to follow employ a single compartment tank undergoing a damped sinusoidal vertical acceleration where the vertical displacement is defined by,

$$y\left(t\right)=A{e}^{-\xi \omega t}\left(Bcos\left(\omega t\right)+Csin\left(\omega t\right)\right)$$

(27)

where,

$$\begin{array}{cc}\hfill \omega & =2\pi f\hfill \\ \hfill a_1& ={\left(\frac{f}{7}\right)}^2{a}_0\hfill \\ \hfill A& =\frac{a_1}{{({\xi }^2{\omega }^2-{\omega }^2)}^2+{(-2\xi {\omega }^2)}^2}\hfill \\ \hfill B& ={\xi }^2{\omega }^2-{\omega }^2\hfill \\ \hfill C& =-2\xi {\omega }^2\hfill \end{array}$$

(28)

allowing an approximate acceleration definition such that,

$$a\left(t\right)=a_1{e}^{-\xi \omega t}cos\left(\omega t\right)=\ddot{y}\left(t\right)$$

(29)

Finally, as per the Protospace experiment $a_0=200$ m · s${}^{-2}$ and $\xi =0.06$ for the non-dimensional number sensitivity study (next section). However, in computing the non-dimensional scaling graphs (later section) we set $\xi =0.0$ so as to focus on quantifying slosh induced damping.

7. Non-Dimensional Study Results and Analysis

To enable direct comparison (plotting on the same graph) of computed energy budget and slosh loads for different excitation frequencies and liquid masses, normalisation is required. Hence, the computed energy dissipated $\left(E_{Disp}\right)$ and work in $\left(W_s\right)$ are normalised to the maximum change in potential energy of a simulation $\left(E_u\right)$. Further, computed vertical slosh force $\left(F_y\right)$ is normalised and oscillation count employed as follows:

$$\begin{array}{cc}\hfill \left|\right|E_{Disp}\left|\right|& =\frac{E_{Disp}}{E_u}\hfill \\ \hfill \left|\right|W_s\left|\right|& =\frac{W_s}{E_u}\hfill \\ \hfill \left|\right|F_y\left|\right|& =\frac{F_y}{|F_y(t=0)|}\hfill \\ \hfill n_{oscillation}& =\frac{t}{\omega }\hfill \end{array}$$

(30)

Figure 6 and Figure 7 presents comparative ranges for the influence of fluid properties on the normalised sloshing fluid dissipation energy $\left(\right|\left|E_{Disp}\right|\left|\right)$, surface work $\left(\right|\left|W_s\right|\left|\right)$ and vertical sloshing force $\left(\right|\left|F_y\right|\left|\right)$. The influence range of each non-dimensional number on the aforementioned quantities are presented via the shaded region.

Table 7. RMS Variation due to fluid property scaling relative to the actual Protospace (%).

	$\left|\right|{\mathit{E}}_{\mathbf{Disp}}\left|\right|$	$\left|\right|{\mathit{W}}_{\mathit{s}}\left|\right|$	$\left|\right|{\mathit{F}}_{\mathit{y}}\left|\right|$
$Re$	$7.80$	$4.85$	$13.47$
$\overline{\rho }$	$5.59$	$3.24$	$14.60$
$We$	$9.23$	$4.90$	$12.99$
$\theta$	$12.74$	$6.43$	$14.63$
$\theta ={90}^{\circ }$	$94.06$	$47.82$	$34.27$
$Fr$	$77.46$	$59.08$	$13.74$

lists the percentage root-mean-square difference due to the scaling of each non-dimensional property as compared to the actual Protospace values. The analysis reveals the relatively low influence of almost all fluid properties on $\left|\right|E_{Disp}\left|\right|$ where the RMSE $<15\%$ and negligible on $\left|\right|W_s\left|\right|$ as RMSE $<10\%$. The exceptions are contact angle and $Fr$ number. The former is due a contact angle of ${90}^{\circ }$ resulting in no liquid–gas surface instability (Rayleigh–Taylor) and hence no slosh damping (note that the reported $94.06\%$ variation as opposed to $100\%$ is due numerically induced free-surface instability which occurs late in the analysis). At this point in the scaling analysis, the choice of Froude scaling (as opposed to $Re$ etc.) is vindicated (and consistent with industrial practice) as all other fluid properties are shown to have a relatively small effect. This was further also demonstrated by a recent experimental campaign [25].

We finally consider the effect that Froude scaling (via changing the excitation frequency) has on the system, with the results shown in Figure 7. Interestingly, for the Froude number range tested, the major effect is on slosh damping as shown, i.e., $77.5\%$. This concludes the important insight that slosh induced dissipation is highly sensitive to Froude number, and to a far greater degree as compared to slosh forces.

8. Dissipated Energy Scaling Laws

We next focus on quantifying the change in slosh induced dissipation as a function of the fluid physics non-dimensional numbers. As mentioned before, we now set $\xi =0$ resulting in a vertical displacement and acceleration as defined by:

$$\begin{array}{cc}\hfill y\left(t\right)=& -\frac{a_1}{{\omega }^2}cos\left(\omega t\right)\hfill \end{array}$$

(31)

$$\begin{array}{cc}\hfill a\left(t\right)=& a_{1}cos\left(\omega t\right)\hfill \end{array}$$

(32)

where it is instructive to note that the magnitude of tank displacement $y_0=a_1/{\omega }^2$. It is important that this value remain fixed when varying the fluid physics non-dimensional numbers. Setting $\xi =0$ simplifies the calculation of slosh induced dissipation to:

$$E_{Disp}\approx {\int }_0^t{F}_y{u}_{d_y}dt$$

(33)

The % change in computed slosh induced dissipation energy due a change in a non-dimensional parameter (from the Protospace value $E_{Disp}\left(\mathrm{Protospace}\right)$) is expressed as:

$$\left|\right|\Delta E_{Disp}\left|\right|=\frac{E_{Disp}-E_{Disp}\left(\mathrm{Protospace}\right)}{|E_{Disp}\left(\mathrm{Protospace}\right)|}$$

(34)

Figure 8 depicts the discrete $\left|\right|\Delta E_{Disp}\left|\right|$ values for individually scaled non-dimensional numbers as well as the resulting trend line. Also delineated on the figures are the properties for the Protospace (black line) and full aircraft (green line). While no statistical trend was found for $We$ and $Re$ (as a result not shown), the following cubic trend for $\theta$ was found around ${60}^{\circ }$:

$$\begin{array}{cc}\hfill \left|\right|\Delta E_{Disp}\left(\theta \right)\left|\right|=& 4.14\times {10}^{-4}{\theta }^3-0.091{\theta }^2+6.702\theta -161.7\hfill \\ \hfill r^2=& 0.82\hfill \end{array}$$

(35)

Here, $r^2$ is the weighted sum of the squares of the difference between the fitted curve and the data such that $r^2\in [0,1]$ with $r^2=1$ being a perfect fit. Outliers to this trend occur due to lack of Rayleigh Taylor instability for $\theta ={90}^{\circ }$ which leads to $P_{Disp}\to 0$, while an enlarged meniscus results as $\theta \to 0$ (to be investigated in future work).

A very well-defined trend can however be found for the scaling of $\overline{\rho }$ as:

$$\begin{array}{cc}\hfill \left|\right|\Delta E_{Disp}\left(\overline{\rho }\right)\left|\right|=& -0.54{\overline{\rho }}^{-0.802}+117.39\hfill \\ \hfill r^2=& 0.98\hfill \end{array}$$

(36)

Similarly, a very well defined trend is visible for the effect of the Froude number on slosh induced dissipation. The resulting trend line reads:

$$\begin{array}{cc}\hfill \left|\right|\Delta E_{Diss}\left|\right|\left(Fr\right)=& -6.57F{r}^2-12.42Fr+114.3\hfill \\ \hfill r^2=& 0.99\hfill \end{array}$$

(37)

which is clearly quadratic in nature. This quadratic relationship between slosh induced dissipation and $Fr$ number is indeed a new finding and a key contribution of this article. Also note that the difference in dissipated energy between the Protospace experiment $Fr$ number and that of the actual aircraft (or indeed ideal experiment) are by far the most significant of all quantities. This necessitates the need to use the scaling law to estimate the expected slosh induced damping for the full scale aircraft tank. This is dealt with the in next section.

From the above it is evident that $\Delta E_{Disp}$ scales circa linearly with $\rho$ (due Equation (36)) and with $v^2$ (due Equation (37)). The latter is as the maximum tank velocity $v\propto \omega$ due Equation (14). From this, we postulate that it is possible to express dissipated energy as a linear function of tank kinetic energy viz.

$$E_{Disp}\propto \rho v^2$$

(38)

for the slosh under consideration.

We digress briefly to investigate the difference in fitted scaling curves for each of the three cycles of oscillation. This is to confirm that scaling laws obtained per cycle are consistent with the experimental observation that the slosh damping which occurs over the first two to three cycles of oscillation is representative of the energy dissipation potential. The resulting curves are shown in Figure 9 while similar $r^2$ correlations to previously were found in most cases (see equations to follow). Apart from the $\theta$ curves, there is a reasonable correlation between the curve fits obtained for the three cycles with the $Fr$ number curve fearing the best. The $\theta$ correlation is of less concern due is significantly lower influence on slosh induced damping while $Fr$ is the most important.

Cubic trends for the scaling of $\theta$ are found such that:

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^1\left(\theta \right)\left|\right|=& 2.05\times {10}^{-3}{\theta }^3-0.33{\theta }^{2}16.73\theta -289.1\hfill & \hfill r^2=0.80\end{array}$$

(39)

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^2\left(\theta \right)\left|\right|=& 3.58\times {10}^{-3}{\theta }^3-0.6311{\theta }^2+35.91\theta -647.4\hfill & \hfill r^2=0.66\end{array}$$

(40)

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^3\left(\theta \right)\left|\right|=& -2.32\times {10}^{-3}{\theta }^3+4.28\times {10}^{-3}{\theta }^2+2.56\theta -125.8\hfill & \hfill r^2=0.50\end{array}$$

(41)

Better correlations were found for mass ratio:

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^1\left(\overline{\rho }\right)\left|\right|=& -145.9{\overline{\rho }}^{-0.7322}+118.1\hfill & \hfill r^2=0.94\end{array}$$

(42)

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^2\left(\overline{\rho }\right)\left|\right|=& -96.6{\overline{\rho }}^{-1.12}+98.03\hfill & \hfill r^2=0.84\end{array}$$

(43)

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^3\left(\overline{\rho }\right)\left|\right|=& -169.2{\overline{\rho }}^{-0.7495}+130.9\hfill & \hfill r^2=0.91\end{array}$$

(44)

The best trend was found throughout for $Fr$ number with all trend lines being quadratic in nature:

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^1\left(Fr\right)\left|\right|=& -7.692F{r}^2-6.13Fr+109.9\hfill & \hfill r^2=0.99\end{array}$$

(45)

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^2\left(Fr\right)\left|\right|=& -8.998F{r}^2+2.914Fr+101\hfill & \hfill r^2=0.98\end{array}$$

(46)

$$\begin{array}{ccc}\hfill \left|\right|E_{Disp}^3\left(Fr\right)\left|\right|=& -2.229F{r}^2-34.93Fr+131.4\hfill & \hfill r^2=0.96\end{array}$$

(47)

9. Scaling-Laws Application

The scaling between the full aircraft model and Protospace experiment employs an attempted Froude scaling due being limited by experimental design and fluid availability (as discussed previously). In contrast high-fidelity CFD enables the modelling of both the ideal (employing an ideal fluid which achieves perfect similarity for all physics non-dimensional numbers) and practical (employing kerosene) Froude scaled experiments. These are denoted as Protospace* and Protospace${}^{†}$, respectively with the resulting non-dimensional numbers listed in

Table 8. Fluid Properties due to Froude Scaling. (Nomenclature for scaling methods as defined above).

	Protospace	Protospace${}^{†}$	Protospace *	Aircraft
	(Water)	(Cold Kerosene)	(Ideal)	(Cold Kerosene)
$\overline{\rho }$	$1.21\times {10}^{-3}$	$1.39\times {10}^{-3}$	$1.39\times {10}^{-3}$	$1.39\times {10}^{-3}$
$Re$	$161.7\times {10}^3$	$4.67\times {10}^3$	$52.2\times {10}^3$	$52.2\times {10}^3$
$We$	$5.73\times {10}^3$	$2.78\times {10}^3$	$69.4\times {10}^3$	$69.4\times {10}^3$
$Fr$	$3.44$	$1.65$	$1.65$	$1.65$

. Note that by using kerosene (practical $Fr$ scaling or Protospace${}^{†}$) both density ratio and $Fr$ numbers equal that of the full scale Aircraft tank. Further, by manipulating fluid properties in the CFD model it is possible to achieve perfect correlations for all shown non-dimensional numbers as shown for ideal scaling (Protospace*).

We now seek to quantify the difference in computed dissipated energy for the different cases and have two options in which to achieve this viz. CFD or using our scaling laws. Starting with CFD first, both practical and ideal scaled cases were simulated and the resulting computed dissipated energy $E_{Disp}$ compared in

Table 9. CFD simulated dissipated energy. (Nomenclature for scaling methods as defined above).

	Protospace	Protospace${}^{†}$	Protospace *
$E_{Disp}\left(W\right)$	$-37.21$	$-5.32$	$-5.59$

. As expected, the actual Protospace experiment overestimates the slosh induced dissipation while practical $Fr$ scaling (Protospace${}^{†}$) yields a slosh induced dissipation which is within $5\%$ of the ideal experiment.

Considering finally the use of the scaling laws, we aim to scale the actual Protospace experiment to the ideal experiment. As both $\overline{\rho }$ and $Fr$ are to be accounted for, a response surface of $\left|\right|\Delta E_{Disp}\left|\right|$ is created using Kriging interpolation which is shown in Figure 10. The resulting scaling computed dissipated power is $4.83\mathrm{W}$ which is within $14\%$ of the CFD computed value.

Please note that the above derived scaling laws allow for scaling of a SDOF system with a geometric scaling $\lambda =1$. This is as the excitation velocity is affected directly by geometric scaling as per Table 3, i.e., velocity is expected to scale by ${\lambda }^{\frac{1}{2}}$. For a larger aircraft size tank this would increase the impact velocity of the slosh and hence also the energy dissipation.

10. Conclusions

This article presented the first ever volume of fluid (VoF) computational fluid dynamics (CFD) based non-dimensional analysis which characterises violent slosh induced energy dissipation of a tank under vertical excitation. This in the absence of significant phase change. The use of rigorously validated CFD for this purpose was essential as doing so via experimental means is not possible due limitations in available liquids. The work resulted in several new insights which are deemed of significant value to the field of liquid slosh modelling and liquid impact dynamics. Each of these key contributions are summarized next.

The VoF CFD software was first validated via modelling a representative slosh experiment. A sophisticated two phase model was employed which accounted for surface tension, turbulence and gas compressibility (via a weakly compressible formulation). The validation study concluded that the slosh induced energy dissipation is predicted to within $2\%$ of the experimental values. The CFD based mechanical energy budget also identified liquid impact as the major contributor to slosh induced energy dissipation. These constitute the first novel findings of this work.

Next a non-dimensional analysis was detailed. The specific emphasis was on developing a functional relationship between energy dissipation and fluid physics. The latter was characterised by Reynolds number, Weber number, contact angle, gas-liquid density ratio and Froude number. The validated CFD software was then used to estimate the influence of each of these non-dimensional numbers on slosh induced energy dissipation. As would be expected, Rayleigh–Taylor instability was found key to inducing liquid impact due purely verticle excitation. Hence, an appropriate contact angle is key to inducing violent verticle slosh induced energy dissipation. Apart from Froude number, all other non-dimensional numbers were for the first time found to have a small (circa $10\%$) influence on slosh induced energy dissipation. Hence, the Froude number was concluded to be the most important non-dimensional number in characterising slosh induced damping.

The validated CFD was then employed to characterise the change in slosh induced dissipation as a function of the most important non-dimensional numbers for the case of constant amplitude oscillatory tank excitation. Clear correlations (scaling laws) were found for both liquid–gas density ratio and Froude number. The latter scaling laws are linear and quadratic, respectively, which is deemed a significant new finding. From this it was postulated for the first time that violent slosh induced energy dissipation may be expressed as a linear function of tank kinetic energy. The final novel contribution of the article entailed demonstrating the application of the developed scaling laws for Froude number and density ratio to estimate slosh induced energy dissipation for an ideally scaled experiment. Thus, using the developed scaling laws, an accuracy of $86\%$ was achieved when compared to the CFD computed value.