A Simple Model-Based Method for Sloshing Estimation in Liquid Transfer in Automatic Machines

This paper studies the phenomenon of sloshing in the field of automatic machines for packaging liquid products, with specific reference to containers with planar motions. After introducing two equivalent discrete models based on a mass-spring-damper system borrowed from the literature (one linear and one non-linear), a novel method is proposed to evaluate the sloshing height of the liquid, namely the deviation of its free surface at the wall of the container from the equilibrium condition. The merits of this method are that it is easy to use, requiring no experimental evaluation of the system parameters or computationally expensive fluid dynamic simulations, and it gives good results also for highly dynamical motions. Moreover, though this paper focuses on cylindrical containers performing rectilinear movements, the presented technique can be extended to containers of arbitrary shape and generic planar motions. The method is validated by experimental tests using cylindrical containers of different dimensions and many rectilinear motion laws with maximum accelerations up to $12~m/s^{2}$ . The results are compared with those obtained by using other methods of equal complexity available in the literature, showing the effectiveness of the proposed technique.


I. INTRODUCTION
The term sloshing refers to the stirring of a liquid inside a container subject to an external excitation, such as a quick change of velocity. In this paper, we specifically refer to the packaging industry of liquid products, where, on the one hand, the high productivity requires transfer motions of the containers with high accelerations and velocities (roughly up to 10 m/s 2 and 2 m/s), and, on the other hand, it is important to prevent liquids from overflowing and to maintain control on their free surface (e.g., not to dirt the inside of the container, which is a common requirement in the pharmaceutical packaging industry). The simplest way to reduce liquid sloshing is to design a suitable container shape, but this is infeasible in the packaging industry, as the customer chooses the container shape based on market demands.
The associate editor coordinating the review of this manuscript and approving it for publication was Shih-Wei Lin .
For this reason, the problem must be solved by choosing a suitable motion law of the container.
This problem is extensively analyzed in the literature. As for the design of anti-sloshing motion laws, input shaping is applied in [1]- [6] to perform rectilinear or more complex motions of containers with different shapes, Infinite Impulse Response (IIR) filters are exploited in [7] and Finite Impulse Response (FIR) filters are used in [8]. A general study on the last ones is also reported in [9] and [10]. Another common approach to designing the optimal trajectory of a container is the solution of a constrained optimization problem, as in [11] and [12].
An important aspect to consider to design the correct trajectory for a given container filled with a liquid is the evaluation of the maximum sloshing height, namely the highest point reached by the liquid free surface during motion. As a matter of fact, the classical models known in the literature to study the liquid behavior (mass-spring-damper and pendulum models [13]- [17]) do not provide responses with a clear VOLUME 9, 2021 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ practical meaning: they are built with dummy masses whose displacements can be calculated, but these are not immediately related to the actual sloshing height. A simple modelbased computational instrument able to analytically estimate the response of the system in terms of sloshing height for an arbitrary container trajectory would be very useful since it would allow the motion designer to reasonably assess whether the liquid is likely to spill out or surpass a given limit, without experimental testing or computationally-demanding fluid dynamic (CFD) simulations.
In some works, the sloshing height is estimated using a very simple pendulum-based model [1], [18], but no details are provided about the accuracy and the limits of this estimation. In [19] a tuned liquid damper is introduced, and the sloshing height is estimated through the fluid sloshing equations and with a similarity between the tuned liquid damper and a tuned mass damper, but the analysis is developed only for a rectangular container. The same is true for the non-linear model presented in [20], which is described by (complex) equations that are valid only for rectangular vessels. In [21], a two-degree-of-freedom mass-spring-damper is introduced and the sloshing height is computed by solving the dynamic equations of the model: however, some experimental tests or a finite-element model of the container is necessary to obtain all parameters. Experimental or pseudo-experimental results (obtained by fluid dynamic simulations) are also used in [22] and [23], where the liquid behavior is predicted by a model constructed by machine-learning techniques from available data. In [24], an accurate study explicitly related to the evaluation of the sloshing height for intermittent motions in packaging machines is conducted by using the pendulum model for rectangular containers. Here, the limits of applicability of the pendulum model to predict the sloshing height are discussed for different dimensions of the container; in particular, it is shown that the maximum accelerations for the validity of this model are between 1.2 m/s 2 and 2.3 m/s 2 .
To the authors' knowledge, there is no proof in the literature that the mass-spring-damper or the pendulum model may effectively be used to predict the sloshing height for highly dynamical motions. Accordingly, the main objectives and contributions of this work to advance state of the art are the following: • it provides a computational method for the prediction of the maximum sloshing height of a liquid carried by a container subjected to highly dynamical motions: the method is simple yet accurate, namely it does not rely on experimental tests or measurements on the real system, it does not require computationally-expensive CFD simulations, and it is fairly precise for container accelerations up to (roughly) 10 m/s 2 ; • it validates the merits of the aforementioned sloshing-height-prediction method by an extensive experimental campaign, and by comparison with other computational procedures of equal complexity available in the literature.
For our analysis, we will consider low-viscosity liquids (e.g., water) since these are the ones that suffer the most severe sloshing phenomena and, thus, represent a worst-case. Moreover, we will focus on cylindrical containers performing fast rectilinear motions since these are the most common in automatic packaging machines. However, extensions to arbitrary planar motions and generic container shapes will also be discussed.
The paper structure is as follows. Section II introduces two equivalent discrete models for sloshing description. A new model-based method to computationally evaluate the sloshing height is described in Section III. Section IV shows the results of the experimental tests performed to assess the proposed technique, also comparing with other methods available in the literature. Finally, Section V draws conclusions.

II. SLOSHING MODEL
The sloshing phenomenon is well-known in the literature. An analytical development, based on continuum models, is available in [13] and [14], under the following assumptions: the fluid is incompressible, irrotational and non-viscous; capillarity and surface tensions are neglected (which is reasonable for practical applications in the gravitational field); the container is rigid and waterproof, and it undergoes planar translations. For the sake of simplicity, we will focus on cylindrical containers performing rectilinear motions, though the method for sloshing-height estimation that we will present in Section III can be extended to containers of any shape performing arbitrary planar motions. The analysis of the homogeneous response of the continuum model provides the so-called sloshing modes. The natural frequency associated to the generic mode mn is [13]: where R is the radius of the container, h is the filling height, g is the gravity acceleration, and ξ mn (representing a root of the derivative of the Bessel function with respect to the radial coordinate r) is a constant parameter known for every sloshing mode. The values of ξ mn are tabulated in [25]; m denotes the circumferential modes and n the radial ones. Reference [13] shows that circumferential modes greater than one are always negligible. Radial modes have a more relevant contribution to the overall response; however, also radial modes greater than one are usually negligible in many practical applications. Different equivalent mechanical models can be used for a discrete, and thus simpler, description of sloshing. All of them are based on the following assumptions [13]: • they must have the same total mass, moment of inertia, and center of mass of the original system; • they must exhibit the same modes; • the liquid must exert the same wrenches on the container wall.
We can identify three different dynamic regimes for the motion of the liquid free surface, with each one of them corresponding to a different model [13], [16], [17]: • linear motion: valid for limited accelerations of the container and correspondingly limited liquid oscillations, in which the free surface remains approximately planar; • weakly non-linear motion: for oscillations with greater amplitude and closer to the natural frequencies of the system; in this situation, the phenomenon of rotary sloshing appears, which causes the liquid to move out the motion plane of the container; • strongly non-linear motion: valid when the free surface is subjected to very high accelerations and impacts against the container walls. The more complex the motion of the free surface that we want to describe, the more complex the equivalent model. In many practical cases, linear or weakly non-linear models are reasonably accurate so that we will focus on them.
A distinct consideration about damping modeling is in order. In fact, damping is largely due to viscosity, which is neglected when solving the analytical problem of continuum models. For this reason, damping is usually taken into account by using empirical formulas. In our analysis, we will use the formulation found in [13], expressed for a generic radial mode n as: where µ and ρ are, respectively, the dynamic viscosity and density of the liquid.

A. LINEAR DISCRETE MODEL
There are two linear models that are commonly used to describe sloshing: the pendulum model ( Fig. 5 in the Appendix A) and the mass-spring-damper model ( Fig. 1(a)).
In this paper, we use the latter since it gives more accurate results as the container acceleration and the corresponding liquid oscillation amplitude increase (as will be shown in Section IV). In both models, there is a fixed mass m 0 (with moment of inertia I 0 ) in the lower part of the container, and a number of mobile sloshing masses, with the generic one being m n . In the mass-spring-damper model, these masses are connected to the container by springs with stiffness k n and dampers with damping constant c n . Each mass represents a radial mode, while circumferential modes greater than one are neglected. If the container only moves along the horizontal x axis, the motion of m n is described by: where the natural frequency and the damping ratio are given by (1) (setting m = 1 and dropping the corresponding subscript for brevity) and (2), respectively. The termẍ 0 is the acceleration of the container, and x n is the displacement of the nth sloshing mass. The value of the nth sloshing mass m n is obtained by equating, for each mode, the sloshing force acting on the container wall calculated with the continuum model to the sloshing force computed by the equivalent discrete model, for an harmonic excitation, namely [13]: where m f is the total mass of the fluid. In the linear model, the vertical displacement of the liquid center of mass may be neglected, and the free surface remains approximately planar. 1 This means that the shape of the free surface can be described with a plane whose equation in polar coordinates (r, θ) is: where η is the sloshing height at a generic point of the free surface, η n is the sloshing height at the container wall in the motion plane associated with the n-th mode, and the total sloshing height at the container wall ( Fig. 1 The non-linear model used in this Section is borrowed from [15]. In this case, the free surface shape is no longer planar, but it is obtained by the study of the forced response to an harmonic excitation applied to the continuum model: where J 1 is the Bessel function of order 1, and η is: This boundary condition leads the nth sloshing mass, still given by (4), to be constrained to move on a paraboloid 1 Under these conditions, one can compute heights h 0 and h n in Fig. 1(a) [13]: these formulations are not reported, since they have no role in the following study. VOLUME 9, 2021 surface, with a non-linear spring of order w (which gives an elastic force equal to k n x 2w−1 n ) connecting it to the axis of the cylindrical container ( Fig. 2(a)). Reference [15] reports the motion equations of m n along both horizontal axes, x and y. However, experimental tests show that the effect of rotary sloshing (i.e., the movement of the sloshing mass in the direction y perpendicular to the excitation) has minor importance, so that, for the sake of brevity, we focus here only on the motion along the x axis, namely: Here x n = x n /R, ω n and ζ n are still given by (1) and (2), α n is a dimensionless constant of the non-linear spring, and C n = ξ 1n tanh (ξ 1n h/R). Following the instructions given in [15], we choose a cubic spring (with w = 2) and the constant α n in the interval [1/2, 2/3]. These values work well for the dimensions of the containers used in the tests described in Section IV, even if these dimensions differ from the ones of the containers used in [15]. Equation (7) reduces to (3) if linearized.

III. SLOSHING HEIGHT COMPUTATION
Motion equations for the linear or the non-linear model can be solved numerically for every input accelerationẍ 0 of the container, thus giving the evolution in time of the displacements x n of the sloshing masses, but these are not physical parameters. For this reason, it is more interesting to calculate the sloshing height η at the container wall, which is a real, measurable quantity. In the following, we propose a new technique to do so. The method is innovative since: • it is purely model-based; • it is real-time capable; • it can be applied to both the linear and the non-linear discrete model; • it is valid for both linear and weakly non-linear motions (see Section IV); • it can be applied to every container shape, though here we focus on cylindrical ones.
Indeed, the approach proposed in [21] involves some experimental tests or a finite-element model of the container, whereas the formulations used in [19] and [20] are only valid for rectangular container. Also the methods based on machine learning algorithms shown in [22] and [23] require experimental tests or complex fluidodynamic simulations of the liquid behavior. In [18] and [1], the sloshing height is estimated by using the pendulum model (reported in Appendix A), but the latter does not give accurate results for increasing values of the container acceleration (see Section IV).
Our computational procedure is based on the assumption at the basis of the discrete models, namely the equality of mass and center of mass between the continuum model and the discrete ones.
Consider a container with section filled with liquid up to height h. In any one of the discrete models represented in Figs. 1 and 2, the horizontal displacement of the center of mass is given by: When the container is moved along the x axis, the liquid free surface during sloshing is described by a general shape η(x, y, η) in the continuum model, with the position of the center of mass along x being given by: If the shapes of the container and the liquid free surface are known, one can obtain the value of η by equating (8) and (9).
Though in this paper we consider only cylindrical container, (8) and (9) are applicable to containers with an arbitrary shape (as long as a suitable discrete model to describe the system, i.e., to compute x n , is available).

A. LINEAR DISCRETE MODEL
In the linear model the free surface is assumed to be planar and given by (5), so that (9) yields: By substituting this result back in (8) and considering that m f = πR 2 hρ, the sloshing height at the wall of the container can be obtained as a function of the displacements x n of sloshing masses: Appendix A shows that, for small oscillations, the estimation given by the pendulum model is analogous to that offered in (10). However, this is no longer true for increasing values of the container acceleration and liquid oscillation amplitudes, for which the pendulum-model formula is no longer accurate, whereas (10) still is (see Section IV).

B. NON-LINEAR DISCRETE MODEL
Here, we consider the simplified model described by (7). Using the free-surface shape given in (6), (9) gives: Substituting this result back in (8) yields: One of the assumptions of the linear model is the irrelevance of the vertical displacement of the center of gravity, but this assumption is no longer valid for the non-linear model. This means that we can equate not only the displacement of the center of gravity along the x direction but also the one in the vertical direction z. By doing so, we obtain the same result expressed in (11) as further proof of the model coherence. Therefore, only one equation is sufficient to find (11), and it is irrelevant whether the conservation of the center of gravity along x or z is imposed.
By comparing (10) and (11), we notice that the contribution of every mode in the two formulations differs for the scale factor 4/ξ 2 1n . In particular, by only considering the fundamental mode (i.e. n = 1), the relation found for the linear model (10) always overestimates the sloshing height by 4/ξ 2 11 = 1.18 compared to (11), since ξ 11 = 1.84118, [25]. For higher modes, ξ 1n is greater than 2 (ξ 12 = 5.33144, ξ 13 = 8.53632, . . . ). So, for equal displacements of the sloshing masses, higher modes have a larger contribution in the non-linear model than in the linear one; however, higher-mode effect is often negligible, since the values of m n rapidly drops for modes greater than one.
So far, we have considered container motions only along the x axis for the sake of simplicity. The extension of our results to arbitrary planar motions of a cylindrical container subjected to excitation also along the y direction is reported in Appendix B.

IV. EXPERIMENTAL RESULTS AND VALIDATION
The method described in Section III for sloshing-height estimation was tested by using the experimental setup shown in Fig. 3(a), comprising cylindrical containers with different diameters (27mm, 50mm, 92mm) filled with a constant volume of water (15.2 ml, 80.5 ml, 500 ml, respectively), because the liquid height (above a certain minimum value) has a small influence on the system response. The water was mixed with a minimal amount of blue dye to better observe the liquid movement. The quantity of dye was so small that it changed the physical properties of the liquid (density, viscosity, and surface tension) in no appreciable way.
The containers were moved by a linear transport system produced by Beckhoff, called XTS, implementing a 1m-long rectilinear rest-to-rest movement for all trajectories. The time evolution of the liquid sloshing was recorded by a camera mounted on the XTS cart, thus moving with the container. The frames analyzed in the recorded videos are similar to those shown in Fig. 3.
Several trajectories were assigned to the three container types, the liquid response was recorded, and the corresponding sloshing heights were compared with those obtained by the numerical method presented in Section III. In particular, the latter was applied to four models: • linear mass-spring-damper model (L-MSD): the displacement of the sloshing mass was found by solving (3), and the result was introduced into (10) to obtain the sloshing height; • non-linear mass-spring-damper model (NL-MSD): the displacement of the sloshing mass was found by solving (7), and the result was introduced into (11) to obtain the sloshing height; • linear pendulum model (L-PEN): the angle described by the sloshing mass with the axis of the container was found by solving (14), and the result was introduced into (13) to obtain the sloshing height; • non-linear pendulum model (NL-PEN): the angle described by the sloshing mass with the axis of the container was found by solving (12), and the result was introduced into (13) to obtain the sloshing height. In all cases, the numerical problem requires solving a (linear or non-linear) ODE ((3), (7), (14) and (12), respectively), with initial conditions x n (0) =ẋ n (0) = 0 (MSD models) or θ n (0) =θ n (0) = 0 (PEN models). In this paper, we integrated all ODEs by the MATLAB function ode45 (based on an explicit Runge-Kutta (4,5) formula), with specified relative and absolute error tolerances equal to 1e − 8. For all the models (linear and non-linear, mass-spring-damper and pendulum) solving the motion equation to find the evolution in time of the sloshing height only requires a few hundredths of a second. For this reason, for both the linear and the non-linear mass-spring-damper model, 3 sloshing masses are considered to find the responses shown in Fig. 4; as a matter of fact, adding a sloshing mass to the simulation means solving an additional ODE, with a negligible time increment. Fig. 4 shows a few of the obtained numerical and experimental results: many other tests were performed with different durations and motion laws to validate the models. In particular, the graphs in Fig. 4 represent the evolution in time of the sloshing height at the wall of the container in different cylindrical vessels subjected to different movements: • the green line is the trend obtained from the experiments. More information on the trajectories used for the constructions of the graphs shown in Fig.4, are listed in Table 1: • the first column indicates the corresponding trajectory in Fig. 4; • the second column indicates the container to which the motion law is applied: small, medium, or big (vessels characteristics are listed at the beginning of this Section); • the third column indicates which kind of motion law is used: modified trapezoidal acceleration, poly5, or a trajectory derived from these by applying an input shaper or a filter; • the fourth column indicates the motion law duration; • the fifth column indicates the maximum acceleration reached during motion; • the sixth column indicates the maximum velocity reached during motion. Tables 2 and 3 summarize some data extrapolated from the graphs in Fig. 4. Both tables have the same columns, but Table 2 refers to the liquid behavior during the container movement, whereas Table 3 refers to the liquid motion during the container rest phase: • the first column indicates the corresponding trajectory in Fig. 4; • the second column indicates the maximum sloshing height predicted by the linear mass-spring-damper (L-MSD) model; • the third column indicates the maximum sloshing height predicted by the non-linear mass-spring-damper (NL-MSD) model; • the fourth column indicates the maximum sloshing height predicted by the linear pendulum (L-PEN) model; • the fifth column indicates the maximum sloshing height predicted by the non-linear pendulum (NL-PEN) model; • the sixth column indicates the maximum sloshing height reached during the experiment (Exp). The values in brackets are the relative percentage error ( ) between the predicted value (η pred ) and the experimental one (η exp ), i.e.: In general, we may observe that the numerical prediction is always fairly accurate in terms of frequency for all models, but the sloshing height estimated by the pendulum models loses physical meaning for highly dynamical motions. The reason is related to (13); as a matter of fact, this formulation is strictly valid only for small oscillations because, for large ones, θ n tends to reach values in the proximity of 90 o and so η goes to infinity. Better results are obtained with the non-linear pendulum model than with the linear one, but trajectories 5 through 9 and 12 ( Fig. 4(e) -4(i) and 4(l)) show that also the non-linear pendulum model loses adherence to reality for motions with high container accelerations. Instead, the mass-spring-damper models, with the formulation of the sloshing height provided in Section III, give fairly accurate results even in these cases. In particular, by observing the data shown in Table 2, we can see that the error in estimating the maximum sloshing height made by the linear and non-linear mass-spring-damper models is about 30% in the worst case (but usually much smaller) and all predicted values are reasonable. Moreover, as expected, the non-linear model is more precise, thus giving in most cases an error of less than 16%. We notice that trajectories 9 through 12 ( Fig. 4(i) -4(l)) are all anti-sloshing (obtained by applying, respectively, ZVDD input shaping, IIR filters, ZV input shaping, and FIR filters), so even for optimized motion laws (with more complex trends) the mass-spring-damper model gives good results for sloshing-height estimation.
For the sake of completeness, we also report some considerations about the free response of the liquid, i.e., the oscillation of the free surface after the container motion is ended. It is clear from Fig. 4 that all models lose precision at this stage. This is true in particular for the linear models, that cannot catch the oscillations at the end of the optimized trajectories (Traj. 9 through 12), since these are calculated with the objective of damping the residual vibration (for these models, this means that the mathematical response of the first mode goes to zero so that only the oscillations of modes 2 and 3 remain). Both non-linear models perform better, but they make much bigger mistakes in estimating the maximum sloshing height during the rest phase of the container than during its motion, as shown in Table 3. This is perhaps due to the high sensitivity of the free response to initial conditions, namely the exact values that liquid displacement and velocity reach at the end of the container motion, values that our model cannot estimate precisely. In the experiments, rotary sloshing was more prominent when the container was at rest, and the liquid had no longer a preferred excitation direction. Nevertheless, the non-linear mass-spring-damper model still seems to provide the best estimation, though in many cases this is comparable with the one provided by the non-linear pendulum model due to the lower amplitude of the free-surface vibration.
Finally, Table 4 shows the effect of each mode in the responses computed by the mass-spring-damper models: • the first column indicates the corresponding trajectory in Fig. 4; • the second column indicates the maximum sloshing height predicted by the first mode of the linear massspring-damper (L-MSD) model; • the third column indicates the maximum sloshing height predicted by the first mode of the non-linear massspring-damper (NL-MSD) model; Comparison of the maximum sloshing height during the movement predicted with different models and the one measured in the experimental tests. The relative percentage error between the predicted value and the experimental one is shown in brackets.

TABLE 3.
Comparison of the maximum sloshing height during the rest phase (after the container movement) predicted with different models and the one measured in the experimental tests. The relative percentage error between the predicted value and the experimental one is shown in brackets.

TABLE 4.
Comparison of the maximum sloshing height predicted by every mode for both the linear and the non-linear mass-spring-damper model. The percentage value of max(η 2 ) and max(η 3 ) with respect to max(η 1 ) is shown in brackets.
• the fourth column indicates the maximum sloshing height predicted by the second mode of the linear massspring-damper (L-MSD) model; • the fifth column indicates the maximum sloshing height predicted by the second mode of the non-linear massspring-damper (NL-MSD) model; • the sixth column indicates the maximum sloshing height predicted by the third mode of the linear massspring-damper (L-MSD) model; • the seventh column indicates the maximum sloshing height predicted by the third mode of the non-linear mass-spring-damper (NL-MSD) model.
The values in brackets are the percentage values of the maximum sloshing height predicted by the second and the third mode with respect to the one predicted by the first mode, i.e. 100 · max(η 2 )/max(η 1 ) and 100 · max(η 3 )/max(η 1 ), respectively (for both the linear and non-linear mass-spring-damper model). Notice that the maximum sloshing heights in the three modes are not reached in the same instant, so their sum is not the maximum of the overall response expressed in Table 2. The effects of modes 2 and 3 are always almost negligible in the linear model, whereas they are appreciable in the non-linear model.

V. CONCLUSION
This paper analyzed liquid sloshing in containers undergoing fast motions in the field of automatic packaging machines. Based on two equivalent discrete models taken from the literature, one linear and one non-linear, a new method to evaluate the sloshing height of the liquid was proposed for both models. The implementation of the aforementioned models in a numerical solver allows one to easily predict the system's response for every input law assigned to the container without needing experimental tests or computationally expensive fluid dynamic simulations, thus representing a useful tool for the motion-law designer. Experimental results conducted on cylindrical containers of different dimensions proved the effectiveness of the proposed computational method to predict the sloshing height in rectilinear movements even in the case the container performs highly dynamical motions (i.e., accelerations up to 11.8 m/s 2 and velocities up to 2.7 m/s). It must be observed that the prediction of the liquid response is not always accurate during the entire motion for very fast container movements, since the free surface exhibits complex behaviors and it does not preserve a planar or Bessel-function shape, but the estimation of the maximum sloshing height is always rather accurate, with errors always smaller than 30% and in most cases comprised between 1% and 17%.

APPENDIX A PENDULUM MODEL
The main assumption here is that the liquid free surface always remains planar, so that its oscillation can be modeled as a pendulum attached to the line normal to the liquid plane (Fig. 5). In this case, if the container only moves along FIGURE 5. Pendulum model. the x axis, the motion equation is (see [8]): θ n + 2ζ n ω nθn + ω 2 n sin(θ n ) = −ẍ 0 l n cos(θ n ) (12) and the sloshing height is: In principle, (12) can be numerically integrated, for a given excitation, even for high input accelerationsẍ 0 and, thus, large swinging angles θ n . However, the larger the latter, the less reasonable is the initial assumption of a planar liquid free surface, and thus the less accurate are the expected results.
The pendulum dynamics can be linearized for small displacements around the equilibrium point θ n = 0, thus obtaining a formulation equivalent to the one introduced for the mass-spring-damper model in (3): where, for equivalence with the natural frequency in (1): Since, for small displacements, tan(θ n ) ≈ sin(θ n ) = x n /l n , the sloshing height, from (13), becomes: Since the computation of η in (16) is not affected by the sloshing mass m n , it is reasonable to take only the first-mode displacement x 1 into account, otherwise higher-mode displacements would give unjustified large contributions. Therefore, considering only the fundamental mode (n = 1) and comparing (10) and (16), we obtain that (16) underestimates the sloshing height with respect to (10) by the constant term: [25], we have C h = 0.9874 ≈ 1. This proves that (10) and (13) give very similar results, but only under the assumption of small oscillations (and thus small container accelerations).

APPENDIX B SLOSHING HEIGHT FOR PLANAR MOTIONS OF A CYLINDRICAL CONTAINER
The motion equations for the complete 3-dimensional non-linear model introduced in Section III-B, with an excitation along both the x and the y axis (ẍ 0 andÿ 0 ), are [15]: x n + 2ω n ζ n ẋ n + C 2 n x 2 nẋ n + y nẏn x n + C 2 n x nẋ 2 n + x 2 nẍ n + x nẏ 2 n + x nÿn y n VOLUME 9, 2021 + ω 2 n x n 1 + α n x 2 n + y 2 n w−1 +ẍ 0 R = 0 y n + 2ω n ζ n ẏ n + C 2 n y 2 nẏn + x nẋn y n + C 2 n y nẏ 2 n + y 2 nÿn + y nẋ 2 n + y nẍn x n + ω 2 n y n 1 + α n x 2 n + y 2 n w−1 +ÿ 0 R = 0 By solving these two differential equations, we may obtain the displacement of the generic sloshing mass along the x and y axes (x n and y n , respectively). In any instant, the sloshing mass lies in a plane passing through the axis of the cylindrical container and with orientation φ n with respect to the positive x axis (see Fig. 6). The value of φ n may be found by knowing x n and y n as: φ n = atan y n x n FIGURE 6. Top view of the center of gravity and sloshing mass of the 3-dimensional non-linear discrete model. φ n also indicates the plane of the maximum sloshing height for the continuum model. Accordingly, on this plane we can apply the same method described in Section III to calculate η. In (8), we only need to replace x G with r G and x n with r n , thus obtaining: r G m f = n r n m n . Equation (9) is also applicable in the radial dimension to find r G by substituting x with r. Finally, considering that r 2 n = x 2 n + y 2 n , we obtain: η = n ξ 2 1n m n x 2 n + y 2 n πR 3 ρ