Mini-AUV Hydrodynamic Parameters Identification via CFD Simulations and Their Application on Control Performance Evaluation

Abstract

This manuscript presents a fully detailed methodology in order to identify the hydrodynamic parameters of a mini autonomous underwater vehicle (mini-AUV) and evaluate its performance using different controllers. The methodology consists of close-to-reality simulation using a Computed Fluid Dynamics (CFD) module of the ANSYS™ Workbench software, the processing of the data, obtained by simulation, with a set of Savistky–Golay filters; and, the application of the Least Square Method in order to estimate the hydrodynamic parameters of the mini-AUV. Finally, these parameters are considered to design the three different controllers that are based on the robot manipulators theory. Numerical simulations are carried out to evaluate the performance of the controllers.

1. Introduction

The increasing implementation of underwater vehicles in ocean exploration for scientific or commercial benefits, rescue maneuvers, and transporting operations implies an efficiency improvement in matters of design, as well as the development of new control techniques and methodologies, which, in conjunction, evoke significant advances in terms of technology for underwater applications [1,2,3].

The identification of the hydrodynamic parameters of autonomous underwater vehicles (AUVs) at the early stages of conception stands as the main motivation of the present work, since the hydrodynamic parameters could potentially serve as a measure of the design efficiency [4,5,6,7].

In order to estimate the hydrodynamic parameters of a prescribed underwater vehicle, different experimental methodologies have been reported in the literature [4,5,6,7,8,9,10,11,12,13,14,15,16,17], yet, the majority of them require the implementation of top-technological expensive gadgets and tools, which makes the design task difficult for those researchers whose budget and facilities are an issue to overcome. In this regard, the implementation of Computed Fluid Dynamics (CFD) software and other computational resources has emerged as an alternative to conceive the vehicles at the first stages of research.

Recent publications provide evidence of the efforts to improve hydrodynamic parameters estimation techniques, emphasizing the increasing importance of the implementation of computational resources. In [18], a maximum likelihood identification algorithm in conjunction with CFD calculations is proposed and validated via experimentation. On the other hand, the hydrodynamic parameters are treated as a function of the angular position in the performance of sharp maneuvers, as exposed in [19,20]. A control-oriented modeling approach is proposed for a low-speed complex-shaped semi-AUV at [21], and a low-cost, efficient CFD procedure is introduced to compute the parameters.

The methodologies that were adopted by a scientist to accurately estimate the vehicles hydrodynamic characteristics are mainly those that have been used in experimental studies; thus, the issue to overcome, besides the computational cost, is to imitate the real scenario conditions in the software environment; issue that has been widely discussed in the literature and that has gathered the attention of researchers in early decades [6,7,8,9,10,17].

The authors consider that the efforts for modeling the hydrodynamic effects in a virtual environment have produced a faster conception and evaluation processes of the vehicles, as such, an action facilitates the recursive improvement process, similar to that of the control algorithms validation [1,2,3,22].

In the control vein, the knowledge of the hydrodynamics coefficients and the physical properties of the system could lead to a potential misconception of the control laws, whereas they are used for position stabilization [23] or trajectory tracking and path planning [22], among other complex tasks [24,25,26]. Nonetheless, the estimation techniques and observers theory have been adopted in maritime engineering, since they allow for the vehicle to know the external disturbances and the hydrodynamic effects that deteriorate the system itself, extending the concept of “Autonomous” vehicles.

At the early stages of the design of underwater vehicles, the stability of the system is a concern that draws the attention of researchers, since the success of the majority of other maneuvers highly depends on the capability of the robot to stay at a given desired point, and to remain at a specific orientation, as could be the study of sunken ships, species observation, and the monitoring of underwater structures, to mention a few [1,2,3].

A brief survey on this topic gives evidence of the diverse techniques that were implemented to ensure the stability of underwater vehicles [23,24,25,26], being the PD and PID controllers, with different variations, the most commonly adopted. The aforementioned controllers have been implemented in order to accomplish tasks as depth control under continuous waves conditions [27,28,29,30] or the path following adversarial environmental conditions [31]. Nonetheless, robust control techniques have been adopted in order to improve the behavior of the system, as in [32,33,34], where the corresponding authors show the effectiveness of the Sliding Mode Control (SMC) applied over these vehicles, extending the capabilities of the theory to the observation domain [35]. Moreover, some techniques [36] consider an iterative identification method to estimate the parameters of the system and the design of the controller at the same time.

In all of the cited works, the control design considers that the physical properties and the hydrodynamic coefficients of the system are known, even when some observers are implemented to compute such coefficients on-board; these are needed before the real operation takes place, as the estimations must be validated. Thus, characterizing the vehicle in terms of the hydrodynamics has a considerable impact on the overall performance of the AUVs. Problem Statement and Paper Contributions The early identification of the hydrodynamic effects that a mini-AUV experiences during motion is the issue to overcome in this manuscript by:

• The adaptation of a yet-virtually-translated experimental test set [11,12,13] to the environment of the CFX on ANSYS™ Workbench software (ANSYS CFX) may be referred to as the industry-leading computational fluid dynamics software for turbo-machinery applications.

• The simulation of the mini-AUV as a submerged body in the ANSYS CFX environment for determining the forces, accelerations, and velocities during motion.
• The conception of the control laws and the stability study, including the estimated hydrodynamic parameters of the mini-AUV, through the robot manipulator control theory, in order to validate the behavior of the system.

It must be highlighted that, even when Finite Volume Method (FVM) simulations are most commonly adopted for hydrodynamics and fluid applications, this work considers, at the same time, the FEM simulation, as it is more appropriate for structural applications, which include forces definition. The ANSYS CFX module manages both methods, which allows for an accurate simulation of both phenomena [37,38]. Nonetheless, recent works on the literature report, for specific cases of study and simulation conditions, slight differences in computation efficiency between the methods mentioned above, and the results do not vary significantly [39,40].

In this regard, the current manuscript is structured by six sections and an appendix; the description of the mini-AUV, including its dynamics, is presented in Section 2. Section 3 establishes the procedure in order to estimate the corresponding hydrodynamic parameters and the conditions programmed within the software environment. In Section 4, the conception of the controllers and corresponding Lyapunov stability analysis are presented, while Section 5 presents the results of the hydrodynamic parameters identification process and the position control task. The concluding remarks and the proposal of future work are comprised in Section 6. Finally, Appendix A contains all of the data generated during the computing process of the hydrodynamic parameters that are exposed in Section 3.

2. Mini-AUV Mathematical Modelling

The current work considers the mini-AUV prototype that was conceived at the Institute and is depicted in Figure 1, in which the four propellers propulsion system is easily identified. Furthermore, the inertial reference frame^I _{O^{I xI yI}z}and the body–fixed reference frame^b _{O^{b xb yb}z}are defined. It is worth highlighting that the origin of the body–fixed reference frame coincides with the center of gravity_Cgof the mini-AUV._Cbstands for the buoyancy center of the vehicle, which is supposed to be aligned with_Cg.

The maximal dimensions of the vehicle (established as0.6×0.42×0.24m in^bx,^by, and^bz , respectively) place it in the “mini vehicle” category [23,41]. According to the CAD software that was used for its conception (Solidworks™), the volume of the mini-AUV is 0.006 m³, which in conjunction with its mass,m=7.5kg, defines the principal moments of inertia about the axis^bx,^by, and^bzas_Ix=0.11^kg·m2,_Iy=0.199^kg·m2, and_Iz=0.251^kg·m2, respectively.

The additional information that is provided by Figure 1b,c regarding the position of the motors, with respect to the body–fixed reference frame, is used in the following paragraphs to define the dynamics of the mini-AUV.

According to the Newton–Euler method [42,43], the dynamics of a rigid body of massm>0∈IRwith six-Degrees of Freedom (DOFs), and whose principal axis of inertia (with their respective moments of inertia_Ix,_Iy,_Iz>0∈IR) coincide with the axis of the body–fixed reference frame, and the corresponding origin matches with the center of gravity, as in the case of the mini-AUV, is described by:

_{Mrb q˙rb}+_{Crb qrb}+_Grb=_τrb

such that:

_qrb=ξ˙ω,_Grb=mg₀₃,_τrb=_Ftτt∈^IR6;_Mrb=m_I30303I˜,_Crb=₀₃₀₃₀₃SωI˜∈^IR6×6

whereξ=^xyzT∈^IR3, andω=^pqrT∈^IR3denote the position of the vehicle in the space and the rotational velocities in the body–fixed reference frame, respectively,g=^00−gT∈^IR3stands for the gravity vector,₀₃∈^IR3does for the zero vector,_I3∈^IR3×3and₀₃∈^IR3×3stand for the identity and the zero matrices,I˜=diag_Ix,_Iy,_Iz∈^IR3×3does for the rotational inertia tensor, andS★is the skew-symmetric operator over the vector★∈^IR3. In this regard,_Ftand_τtinclude, correspondingly, the forces and torques that are produced by the propulsion system_fp,_τp∈^IR3, and the hydrodynamic effects due to the motion of the vehicle_fh,_τh∈^IR3, not to mention the influence of the external disturbances, which is neglected in this paper. In this sense:

_τrb=_Ftτt=_fpτp+_fhτh

Taking the geometry of the mini-AUV depicted in Figure 1, and the forces produced by each motor_fi∈IR(withi=1,2,3,4), into consideration, it leads to:

_fp=_fxfyfz=_{Cθ CψSϕ Sθ Cψ}−_{Cϕ SψCϕ Sθ Cψ}+_{Sϕ SψCθ SψSϕ Sθ Sψ}+_{Cϕ CψCϕ Sθ Sψ}−_{Sϕ Cψ}−_{SθCθ SϕCθ Cϕf2}+_f40_f1+_f3=R_f2+_f40_f1+_f3

_τp=_τpτqτr=_lyf3−_f1+d_f2+_f4lzdf2+_f4lyf2−_f4+d_f1+_f3

whereR∈^IR3×3provides a vector mapping from to the body–reference frame to the inertial frame and it corresponds to the rotation matrix that is described by the Euler angles contained in the vectorη=^ϕθψT∈^IR3whose rate of change,η˙∈^IR3, is related toω, such that:

ω=Wη˙withW=10−_Sθ0_{CϕSϕ Cθ}0−_{SϕCϕ Cθ}∈^IR3×3

Additionally,d∈IRstands for the proportional factor that relates the force that is produced by the propeller to its corresponding exerted free moment_τi=d_fi∈IR.

On the other hand, the hydrodynamic forces and torques,_fhand_τh∈^IR3, respectively, are produced by complex phenomena, which include the added mass effects, radiation-induced potential damping, and restoring forces. The magnitude of these terms depends on the velocities of the body expressed in the body–reference frameυ=^RTξ˙andω ; thus, they are provided in the reference frame that is mentioned above and modeled by an equation of the form of Equation (1) with an added mass inertia matrix_Ma>0∈^IR6×6, a Coriolis terms matrix_Ca∈^IR6×6, a damping matrix_Dd>0∈^IR6×6, and a vector of restoring forces_Gr∈^IR6, which comprises the buoyant forceB∈IR . A fully extended treatment of these terms can be consulted at [11,12,13,44,45,46,47]. Furthermore, the SNAME notation, regarding the hydrodynamic parameters, is adopted now-on, and the hydrodynamic effects are included within the dynamics of the mini-AUV, such that Equation (1) can be rewritten as:

_{Mt q˙rb}+_{Ct qrb}+_{Ddm qrb}+_Gt=_fpτp

with_Mt=_Mrb+_Mam ,_Ct=_Crb+_Cam ,_Gt=_Grb+_Grm . It must be highlighted that the matrices_Mam ,_Cam ,_Ddm and the vector_Grm correspond, respectively, to_Ma,_Ca,_Dd, and_Grwith the proper conversion of the translational terms to the inertial frame.

Because the mini-AUV dynamics is highly nonlinear and coupled, several assumptions are suggested to be considered throughout the literature [11,12,13,46,47]. In this regard, the off-diagonal elements of the matrix_Maare neglected, as they are smaller than those of the diagonal; the mini-AUV is assumed to be naturally stable in roll and pitch (ϕ=θ=^0∘andp=q=0deg/s∀t≥0s), andv=v˙≈0, as the actuators only produce forces along the body axis^bxand^bz. Additionally, the mini-AUV operates at low speeds, such that the Coriolis and centripetal effects, and the quadratic damping terms, are small enough to be neglected yet, the factor_Xuuis still considered, sinceu>>wandu>>r; thus, the six-DOFs model of the mini-AUV can be reduced to a four-DOFs dynamical model of the form:

_Mqq¨+_Cqq˙+_Dqq˙+_Gq=^_fpTτrT

_Mγγ˙+_Dγγ+_Gγ=_τγ

where Equation (8) holds for the motion that is described in the inertial frame, and Equation (9) does the proper for the motion described in the body reference frame. Additionally,q=^xyzψT∈^IR4,γ=^uvwrTand

_Mq=J_Mγ ^JT;_Mγ=diagm+_Xu˙,m+_Yv˙,m+_Zw˙,_Iz+_Nr˙∈^IR4×4

_Cq=−_{Sψ CψYv˙}−_Xu˙−_{Sψ2 Yv˙}+_{Cψ2 Xu˙}00_{Cψ2 Yv˙}+_{Sψ2 Xu˙Sψ CψYv˙}−_Xu˙0000000000ψ˙∈^IR4×4

_Dq=J_Dγ ^JT;_Dγ=diag_Xu+_Xuuu,_Yv,_Zw,_Nr,∈^IR4×4

_Gq=^00B−mg0T;_Gγ=^J−1 _Gq∈^IR4

_τγ=^{J−1 _fpT_τrT}=^XYZNT∈^IR4

J=_Cψ−_Sψ00_SψCψ0000100001∈^IR4×4

Notice that the v dynamics has been considered to keep the matrices within^IR4×4.

Equations (8) and (9) are respectively used in the following sections to identify the hydrodynamic parameters of the mini-AUV, and to conceive the controllers and their stability analysis; meanwhile, Equation (7) was used to carry out the corresponding numerical simulations to validate the controllers.

3. Hydrodynamic Parameters Identification

The experimental methodology that was presented in [11,12,13,14,15] was considered to identify the hydrodynamic parameters of the vehicle. This methodology implies the study of the performance of the vehicle in only one DOF at a time. The vehicle is constrained on the other 5-DOFs.

The process to identify the parameters of the mini-AUV consists of finding the numeric values of the unknown constants of Equation (9). To this end, ANSYS CFX was used to simulate the real environmental conditions in order to obtain the forces acting over the mini-AUV and the corresponding velocities. The accelerations were computed using a Savitzky–Golay filter and the Least Squares Theory was implemented to treat the overall data [44,45].

3.1. Basis of the Least Squares Method

Through the basis of the Least Squares Theory, an estimation of unknown parameters can be obtained by following the set equations:

_θ^LS=^{HTH−1 HT}ϑ;_σ^θ=diagcov_θ^LS;

cov_θ^LS=^{σ^2 HTH−1};^{σ^2}=^{ϑ−H_θ^LST}ϑ−H_θ^LSdimϑ−dim_θ^LS

where_θ^LSis the vector of estimated parameters, H is the deterministic matrix model,ϑthe measurement vector (forces vector),_σ^θis the variance, andcov_θ^LSis the covariance; meanwhile,^{σ^2} stands for the standard deviation of the studied data set [48,49].

The Equation (9) model can be extended in order to make the equations of motion satisfy the required format to apply the Least Squares Methodology. Thus, it follows, for the surge (^bx) motion representative example, that:

_mxu˙+_Xuu+_Xuuuu=X

with the constant_mx=m+_Xu˙.Furthermore, the mini-AUV is observed to have no symmetry with respect to the^by^bzand^bx^by planes (Figure 1), which implies that the hydrodynamic parameters for the positive and negative motion must be estimated [11,44,45]. Notice that, in yaw motion, the parameters remain the same for turning right or left due to the symmetry of the vehicle in the^bx^bzplane.

Equation (18) can be expressed in a matrix shape, for n data points, as:

_X1X2X3⋮_Xn=H_mxXu˙XuuwithH=_{u˙1u1u1u1u˙2u2u2u2u˙3u3u3u3}⋮⋮⋮_u˙nununun∈^IRn×3

that fits to Equations (16) and (17) as_θ^LS=^_mxXu˙XuuTandϑ=^{_X1X2X3⋯_XnT}.

Next, the conditions that are considered in the software ANSYS CFX to perform the simulations are described. 3.2. ANSYS CFX Simulation Environment

As a first step into the simulation software, the domains were defined, one for the fluid domain and the other for the rigid body, i.e., the mini-AUV. The control volume or fluid domain was established to be a cube with dimensions15m×15m×15 m. Afterward, the CAD model of the mini-AUV was imported to the software and then placed at the center of the cube to ensure the fully submerged condition, as shown in Figure 2.

The sweep method and two different sizings were programmed to generate the fluid mesh. The sizings were applied to the edges of the cube to obtain more elements near the mini-AUV [8,10]. Both hlsizings were set to give 70 divisions on the edges in^Iyand^Iz . The change rate, which defines the ratio size between the inner and the outer elements of the edge, was set to 5, obtaining 343,000 elements and 357,911 nodes with a mean element quality of 0.752476. Figure 2b depicts the result of this configuration.

Regarding the mesh of the mini-AUV, the Hex dominant method with all quad option was used. A sizing of 10 mm was also programmed with a change rate of 1.1 and a curve sensitivity of^25∘ , resulting in 32,708 elements and 25,428 nodes whose mean element quality was 0.751376. Figure 2c offers a visual reference of the mini-AUV mesh.

In total, 375,708 elements and 383,339 nodes were generated with a global mean quality of 0.751476, which coincides with the literature reports [8,9,10,37], and that allows to consider such a mesh to be valid to perform the simulation.

The simulation time is defined as 30 s with two different time steps of0.2and0.25s in order to obtain 150 and 120 data points, respectively, to have a rough approach to the time step influence.

The default properties of the water available on ANSYS CFX were used to simulate the behavior of the mini-AUV in^bx [37]. The buoyant effect was activated, while the mini-AUV was simulated as an immersed solid with the properties provided in Section 2. The scaling function was set to 100, as suggested in [37]. The submerged body option was also selected, allowing to apply the force, Figure 3, in order to produce the forward and backward motion of the mini-AUV.

Lastly, at the bottom and the sides of the fluid domain, the wall conditions were established; meanwhile, an opening was assigned at the top surface. The translational motion of the mini-AUV in^Iyand^Iz was constrained, as experimentally suggested by [11]. The results of the simulation and the data treatment are introduced and discussed in Section 5.1.

4. Controllers Conception

The controllers were designed with the mini-AUV hydrodynamic parameters when considering the control theory of robot manipulators exposed in [43]. For these ends, Equations (8) and (9) are cited, highlighting that the current analysis stands exclusively for the position control case, which implies that the desired point_qd=^_xdydzdψdT∈^IR4does not vary with time, such that_q˙d=_q¨d=0∈^IR4∀t≥0.

4.1. PD Controller

The PD control law with gravity compensation [43] can be adapted to the mini-AUV dynamics in Equation (8), such that:

_τq=_KPq˜−_KVq˙+_Gq

whereq˜=_qd−q∈^IR4corresponds to the position error vector,_KP>0∈^IR4×4stands for the proportional gains matrix, and_KV>0∈^IR4×4does for the derivative gains matrix. The expression above can be related to_τγby the matrix J, such that_τγ=^JT _τq ; thus, Equation (20) becomes:

_τγ=^JT _KPq˜−^JT _KVJγ+_Gγ

then, in conjunction with Equation (9), defines the closed-loop equation of the system as:

_Mγγ˙+_Dγγ+_Gγ=^JT _KPq˜−^JT _KVJγ+_Gγ

with the state of the system being defined asx=^q˜TγTT∈^IR8, and whose equilibrium point_xe∈^IR8is located at^0T0TT∈^IR8 that can be computed by considering the closed-loop Equation (22) and solving:

x˙=−q˙γ˙=−Jγγ˙=₀₄₀₄

The stability of the system can be further addressed considering a Lyapunov candidate function of the form:

Vq˜,γ=12^γT _Mγγ+12^q˜T _KPq˜≥0

whose time derivative is given as:

V˙q˜,γ=−^{γT JT} _KVJγ−^γT _Dγγ

Accordingly,V˙q˜,γ≤0∀t≥0, since_KVand_Dγ are positive definite matrices by definition, thus the equilibrium is said to be stable [43].

Moreover, the LaSalle theorem leads to conclude that such stability is global asymptotic, since the setΩ is strictly defined, from Equations (22) and (25), as [42,43]:

Ω=x:V˙q˜,γ=0=q˜=0∈^IR4,γ=0∈^IR4

Following this vein, a PID controller conception, following the robot manipulators control theory, is introduced next. 4.2. PID Controller

The PID control law is proposed with a base on the theory available at [43]; in this regard, and compared with the PD controller, it includes a diagonal integral gain matrix_KI>0∈^IR4×4such that:

_τq=_KPq˜−_KVq˙+_{KI ∫0t}q˜dt+_Gq

The integral action of the PID control law in Equation (26) introduces an additional state variable that is denoted herein byφ∈^IR4, and whose time derivative isφ˙=q˜∈^IR4, leading to:

_τq=_KPq˜−_KVq˙+_KIφ+_Gq

Equation (27) is related to the model in Equation (9), as explained in Section 4.1, in the sense that:

_Mγγ˙+_Dγγ+_Gγ=^JT _KPq˜−^JT _KVJγ+^JT _KIφ+_Gγ

To proceed with the stability analysis, the state vector, for this case of study, is defined, as follows:

χ=^φTq˜TγTT∈^IR12

The only equilibrium point_χe∈^IR12 of the system matches the origin. Nevertheless, a Lyapunov stability analysis is conducted to prove it is asymptotically stable. In this regard, a global change of coordinates is executed to obtain a new state vector and its corresponding time derivative, respectively, as [42,43]:

^χ*=^φ*q˜γ=ιφ+q˜q˜γ∈^IR12;^χ˙*=ιq˜−Jγ−Jγγ˙∈^IR12

from where the closed–loop equation can be written as:

_Mγγ˙=^JT _KPq˜−^JT _KVJ+_Dγγ+^JT _KIι^φ*−q˜

withι∈IRbeing a small positive constant.

Equation (31) results easy to be programmed due to the coordinates change as the integral term has been substituted, but the Lyapunov stability analysis is extensive. In this regard, Equation (28) can be translated to the inertial reference frame in order to develop the stability analysis; in this manner, the control law closes the loop, expressed as:

_Mqq¨+_Cqq˙+_Dqq˙+_Gq=_KPq˜−_KVq˙+_KIφ+_Gq

With this representation of the system, a new vector of states is defined by^x*∈^IR12, such that:

^x*=^κTqTφTT∈^IR12

whereκ=_Mqq˙∈^IR4 is the generalized momentum of the mini-AUV [46]. Hence, it has been shown by the time differentiation of the Lyapunov candidate function of the form [43]:

V^x*=^x*T_Mq−1ι_I404ι_I4KPKI04KIι_KI^x*

with_I4∈^IR4×4and₀₄∈^IR4×4standing for the identity and the zero matrices, respectively; thatV˙^x*≤0, and thatqconverges to_qdif:

_λminKV>_λmaxMq;_λminKI>0;_λminKP>_λmaxKV+2ι_KI

andιis small enough, such that the following expression holds:

121−ι_KP−ι_Mq+ι2∑i=14_q˜i∂_Mq∂_qi>0

Notice that this solution only guarantees local stability in a limited region about the origin of Equation (33), [42,46,47]. A robust Sliding Mode Control approach is also proposed in order to carry out the mini-AUV behavior comparison.

4.3. Sliding Mode Controller

Based on the theory exposed in [43,50], and the results exposed in [42], the sliding surface vectorS∈^IR4is defined as:

S=−Jγ+_λmq˜

with_λm∈^IR4×4being a diagonal matrix of gains_λqi >0.

The attractive sliding surfaces vectorS˙∈^IR4is established as:

S˙=−_KZsignS

with_KZ∈^IR4×4a matrix of gains_Kzi >0.

Equation (35) corresponds, at the same time, to the time derivative of Equation (34), which could also be defined as:

S˙=−J˙γ−Jγ˙−_λmJγ

Thus:

_KZsignS=J˙γ+Jγ˙+_λmJγ

The stability analysis of this controller is developed according to the theory that was described in [43], which corresponds to robot manipulators control theory.

Let the state of the system be described by the vectorS , from Equation (37), so that it can be concluded that the origin of the system is the only equilibrium point [42,43,50].

From the definitions ofSandS˙ in Equations (34) and (35), respectively, it holds that:

^STS˙≤0

The Lyapunov candidate function is defined as:

VS=12^STS≥0

The time derivative corresponds to Equation (38), thus Equation (39) is a Lyapunov function, and the origin of the system is said to be stable [43].

The system shown in Equation (37) is autonomous; thus, the setΩ might be defined to apply the LaSalle theorem [43]:

Ω=S∈^IR4:V˙S=0

It is straightforward to conclude, from Equation (38), thatS=0 is the only initial condition in (40); thus, the stability of the origin is asymptotic [42,43,50].

The three controllers that are studied here are validated in the following section, in which the CFX simulation results are also included. 5. Simulation Results In this section, the ANSYS CFX simulation results and the treatment of the data obtained are exposed. In the second subsection, the hydrodynamic parameters identification results are used in order to design the controllers and show their accuracy by numerical simulations. 5.1. ANSYS CFX Simulation Results and Data Treatment

The information given by the fluid dynamics simulations was processed while using Matlab™. Figure 4 shows that the position, velocity, and acceleration of the mini-AUV, as provided by the ANSYS CFX module, are depicted, highlighting that such results correspond to its motion in^bx . Recalling [11,12,13,14,15,44,45], the plots can be considered to describe a close-to-reality behavior whose evidence is shown in Figure 5 [11,12,13] (Figures taken from [11,12,13] for comparison purposes only, no copyright infringement is intended.).

From Figure 4, it can be concluded that the time step variation of 0.05 s does not have a notorious impact on the displacement and velocity simulation results. On the other hand, notice that the position and velocities of the system follow the curves that were obtained by real experiments, with the differences due to the design and the dimensions of the vehicles being evident. Nevertheless, the accelerations seem not to have the same behavior, as the one provided by ANSYS CFX is noisy, in such a manner that it needs to be filtered, as suggested in the methodologies provided throughout [11,12,13,14,15]. In this vein, the acceleration results that are depicted in Figure 4 were used to estimate the parameters as an isolated possible scenario, as such behavior would potentially lead to a larger estimation error, as explained in the following paragraphs.

The Savitzky–Golay filter was used in order to filter the velocity and acceleration signals since, according to [11,51,52], the main filter characteristic relays on an accurate curve tracking and following, with an attenuation of abrupt signal changes.

In Matlab™, the function of the Savitzky–Golay filter is given by the instructionsgolayfilt(x,k,f). The input parameters of the function are: the vector x, which is the signal to be filtered, the polynomial grade k, and the size of the windowf>k with this being an odd number [51,52].

The data were analyzed with a set of different filter parameters to evaluate the performance of the filter. Moreover, according to the documentation [51,52], the acceleration can also be computed from the velocity by applying the Savitzky–Golay method, which filters the input signal and numerically computes its first derivative. In Matlab™, this method is comprised of the functionsgolay(N,F) , where N is the polynomial grade, and F the size of the window. In Table 1, the parameters of the Savitzky–Golay filter and the corresponding case identifier number (label) are established.

Case 1, as exposed in Table 1, corresponds to the data without any filtering process, i.e., as it was previously obtained from the ANSYS CFX simulation. Cases 2–5 stand for the cases in which a Savitzky–Golay filter was independently applied to the velocity and the acceleration. Figure 6 and Figure 7 depict the results of this filtering stage with the values that are provided in Table 1.

Figure 7 shows an appreciable attenuation of the acceleration; meanwhile, it can be appreciated that the velocity curves have been smothered. Figure 6 suggests that, mainly, the filter at case 4 has a better performance in this matter. In Figure 8, the velocities and the accelerations concerning case 3 (randomly selected) of Table 1 are isolated in order to present a comparison to appreciate the differences between the filtering process with 120 and 150 data points.

Smoother signals are appreciated in the cases where 150 data points were available for the analysis, which, alongside the results presented in [13,14,15], led to the conclusion that the more data points available, the better the estimation process. Nevertheless, conclusions cannot be firmly established at this point of the data treatment process, as the Least Squares Method has not been used yet.

A second filter was applied to the already filtered accelerations obtained at cases 2–5, as suggested by [11,12,13]. For such a task, the second filter parameters, as shown in Table 2, were established, and new labels for these cases were assigned.

Figure 9 depicts the double-filtered accelerations, where the attenuation of the signals is notorious. Figure 10 shows case 8 of the data treatment strategy.

Two supplementary cases were considered. Case 10 corresponds to the implementation of thesgolayfunction to filter the velocity and compute the acceleration from the velocity data provided. Finally, case 11 corresponds to case 10 yet, with an additional Savitzky–Golay filter applied to the acceleration signal that was previously computed by thesgolay function. Table 3 introduces the parameters for these two cases. Figure 11 and Figure 12 depict the corresponding results.

With the data acquired and arranged in the form of Equation (19), the expressions in Equations (16) and (17) were applied in order to obtain the estimated values of the hydrodynamic parameters.

In Appendix A, the tables containing the statistical data for each analysis are presented. Based on the data set that was obtained with 120 data points concerning the positive surge motion of the mini-AUV, the added mass value is observed to adopt a magnitude near9kg, while the linear and quadratic damping coefficients are close to0.1kg/sand6kg/m, respectively.

By definition, the added mass and the damping matrices must be strictly positive definite, implying that_mx>m,_Xu>0, and_Xuu>0 [46,47]. Most of the parameters estimated with the 120 iterations simulation results do not satisfy the conditions mentioned above. On the other hand, on matters of the hydrodynamic coefficients that are computed with the data of the ANSYS CFX simulation with 150 iterations, it was found that more results fit the specifications, which suggests a close relationship between the accuracy of the estimation and the number of data points in such a manner that the more information is obtained, the more accurate the estimation may be. Several works in the literature have also given evidence of this result [8,9,10,37]. The fact that the standard deviation of the linear damping coefficient is close to the value of the parameter could be attributed to the data points used to compute the parameters [4,5,6,16,37].

Based on the lowest variance and standard deviation of the parameters in the Appendix A, and in comparison with the results exposed by [11,12,13,14], the estimation that is selected to be valid is that of the case 11; thus, the hydrodynamics coefficients that characterize the positive^bxmotion of the vehicle correspond to:

_Xu˙=1.63082585kg,_Xu=0.142306865kg/s,and_Xuu=5.153776644kg/m

The same identification procedure was followed in order to estimate each of the hydrodynamic coefficients. The overall set of hydrodynamic parameters of the mini-AUV is available at Table 4.

Upon conclusion of the parameters identification procedure, the information obtained was used to define and simulate the controllers proposed in Section 4.

5.2. Control Simulations Results

Because the vertical motion of the mini-AUV is decoupled from the horizontal plane, the control algorithm is divided into two parts: one for the control depth and the other for the motion in the horizontal plane; as depicted by the block diagram in Figure 13. The motion of the mini-AUV in the^Ixand^Iydirections is coupled to theψdynamics. Subsequently, given a desired point in the space^[_xdydzd]T∈^IR3, the desired angle_ψdis computed when considering the initial positions_x0and_y0, such that:

_ψd=^tan−1_yd−_y0xd−_x0

First, the mini-AUV is driven to the desired angle_ψdand, once it has reached the commanded set-point, it starts to move forwards to reach the desired point_xd,_yd. In the meantime, the depth (z) control drives the mini-AUV to the desired vertical position.

For the simulations, the geometric parameter_lywas set to be0.16 m, and the estimates of the hydrodynamic parameters, as introduced in Table 4, were used. Additionally, the gains of the PD, PID, and SM controllers correspond to those in Table 5. Such gains ensure a slow-motion of the vehicle, as it had been assumed in Section 2 and Section 4.

The simulation was run for a total timetsim=240s with a constant sampling timedt=0.001s in Matlab/Simulink™2018b while using a computer that was embedded with an 8 GB RAM and an Inter Core™i5-8250 CPU @ 1.60 GHz & 1.80 GHz processor. The initial conditions of the system were set to 0 and the reference point to reach was selected to be at_xd=−1m,_yd=1m,_zd=1m . The results of such control studies are exposed through Figure 14, Figure 15, Figure 16 and Figure 17.

In Figure 14 and Figure 15, it is possible to appreciate the performance of the controllers regarding the translational horizontal motion of the mini-AUV. In this matter, the difference in the settling time is notorious.

The mini-AUV under the PD and SM control techniques has an oscillatory behavior, while the system under the influence of the PID controller does not, due to the effect of the integral term. Nevertheless, such oscillatory behavior can be attenuated by a detailed tuning procedure, according to the desired performance or task requirements.

Figure 15 gives evidence of the steady-state response of the mini-AUV, which has a better performance when implementing the PID and SM controllers, since the error tends to 0 as times goes to infinity, which is produced due to the integral term in the case of the PID controller, and to the robust approach and conception of the Sliding Mode controller.

Based on the gains presented in Table 5, it can be observed that the PID and PD control drive the mini-AUV to_ψdat similar times; on the other hand, the SM controller drives the system to_ψd more slowly, yet, in the three cases, the desired orientation seemed to be reached, which is reinforced by the result that is depicted in Figure 16.

Concerning theψdynamics, no oscillation behavior was desired as an arbitrary requirement, thus the set of control gains previously introduced ensures such a requirement. Not a notorious difference is appreciated when the PD and PID controllers are applied, yet the response of the mini-AUV, when controlled by the SMC, is slower, which cannot be considered to be an issue, since it can be fixed by selecting a new set of gains.

In Figure 17, for the motion in^Iz, a difference between the transient response of the three controllers is appreciable, yet, similar results can be obtained by selecting the proper set of parameters, depending on the controller.

Some considerations have to be taken, even when the three controllers accomplish the goal of controlling the position of the vehicle. The results suggest that these controllers can successfully drive the system to the desired point, yet external disturbances influences have not been considered; thus, this scenario must be validated in simulations and reality. Besides, it can be commented that a better performance is expected from the SM controller, since it is a robust control technique. 6. Concluding Remarks and Future Work In this paper, the issue of identifying the hydrodynamic parameters of a mini-AUV in the early stages of the design was treated. In this regard, the dynamic model has been described and simplified in order to study the position stability problem under three different control approaches. The robot manipulator control theory was used to design the three controllers and study the stability of the closed-loop system. The adaptation mentioned above may be considered to be a new control design methodology for this kind of vehicle.

The hydrodynamic parameters of the mini-AUV were estimated by close-to-reality hydrodynamics simulations. A real experimentation procedure and the environmental conditions were adapted; additionally, the Immersed Solid option in ANSYS CFX to analyze the behavior of these vehicles may establish new options and capabilities to be explored and exploited in the upcoming works. Such options, alongside the recent founds and improving methods on the Least Squares estimation technique [20,49,53,54,55], could improve the estimation of such parameters and lead to more efficient simulation algorithms. Nevertheless, some issues that are described in [6] can emerge, thus a careful use of the CFD software is suggested, yet it is a useful tool in early design stages, the estimation of the parameters should take place in real and virtual scenarios.

A detailed and refined 3D model of the mini-AUV can be improved in the upcoming projects to increase the efficiency of the vehicle in operation, based on the results that were provided by the ANSYS CFX software, since the hydrodynamic parameters can be considered to be a tool for measuring the design efficiency [4,5,6,16]. A continuous improvement process also implies exploring more options and several parameters configuration to develop the simulations and obtain more accurate estimations, not to mention the real experimentation and problems that come with the data acquisition task.

The use of a Work Station Computer is strongly encouraged, since this would allow for performing faster simulations; thus, more data points can be obtained, and the precision of the estimation can be improved as well as the quality of the mesh. Performinf more simulations considering different real approaches and new methodologies is part of the upcoming work, as suggested by [4,5,6].

Numerical simulations validated the controllers; however, to control the mini-AUV, when considering the whole phenomena in its dynamics and uncertainties, is opened to new projects at the Institute. This manuscript pretends to serve as an introductory work to the conception procedure of underwater vehicles and their hydrodynamic parameters identification.

Case	u Filter		u˙Filter
	k	f	k	f
1	-	-	-	-
2	1	3	1	3
3	1	5	1	5
4	1	7	1	7
5	3	7	3	7

Case	u Filter		u˙1st Filter		u˙2nd Filter
	k	f	k	f	k	f
6	1	3	1	3	1	7
7	1	5	1	5	1	7
8	1	7	1	7	1	7
9	3	7	3	7	1	7

Case	uSgolay		u˙Savitzky–Golay Filter
	N	F		k	f
10	4	11		-	-
11	4	11		1	7

Parameter	For Positive	For Negative
	Velocity	Velocity
Xu˙kg	1.630825	1.932838
Xukg/s	0.142306	0.091036
Xuukg/m	5.153776	5.868182
Zw˙kg	2.588308	3.908386
Zwkg/s	7.157339	7.632676
ZgN	−14.309872	−14.309872
Nr˙kg·m2	0.0240417	0.0240417
Nrkg·m2/s	0.013791	0.013791

PD	Value	PID	Value	SMC	Value
kpx (kg/s2)	0.2	kpx (kg/s2)	0.012	Kzx (m/s2)	0.05
kpy (kg/s2)	0.2	kpy (kg/s2)	0.012	Kzy (m/s2)	0.05
kpz (kg/s2)	1	kpz (kg/s2)	1	Kzz (m/s2)	1
kpψ (Nm)	0.1	kpψ (Nm)	0.01	Kψ(s−2)	1
kvx (kg/s)	2	kvx (kg/s)	0.005	λx(s−1)	0.7
kvy (kg/s)	2	kvy (kg/s)	0.005	λy(s−1)	0.7
kvz (kg/s)	5	kvz (kg/s)	3	λz(s−1)	0.07
kvψ	0.5	kvψ	0.5	λψ	0.15
		kix (kg/s2)	0.0003
		kiy (kg/s2)	0.0003
		kiz (kg/s2)	0.0005
		kiψ (Nm)	0.00005

Author Contributions

Conceptualization, J.J.C.-Z., and K.A.C.-G.; methodology, J.J.C.-Z., and G.I.P.-S.; software, J.J.C.-Z., and G.I.P.-S.; validation, J.J.C.-Z., G.I.P.-S., and K.A.C.-G.; formal analysis, J.J.C.-Z., G.I.P.-S. and K.A.C.-G.; investigation, J.J.C.-Z., G.I.P.-S., J.R.-R., L.A.M.-H., and K.A.C.-G.; resources, G.I.P.-S., J.R.-R., and L.A.M.-H.; writing-original draft preparation, J.J.C.-Z., G.I.P.-S., and K.A.C.-G.; writing-review, and editing, J.J.C.-Z., G.I.P.-S., J.R.-R., L.A.M.-H., and K.A.C.-G.; supervision, G.I.P.-S.; project administration, K.A.C.-G.; funding acquisition, J.J.C.-Z. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

Simulation data are available on request to José J. Castillo-Zamora.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A. Additional Data Obtained at the Hydrodynamics Parameters Estimation Process

Appendix A provides the statistical data obtained during the estimation of the hydrodynamic parameters exposed throughout the paper. Table A1 and Table A2 correspond to the analysis of the positive motion of the vehicle in^bx , and Table A3 stands for the negative motion in^bx. The results concerning the treatment of the data for the positive and negative motion of the vehicle in^bz are exposed in Table A5 and Table A4, respectively. Lastly, Table A6 shows the parameters estimation process results regarding the rotation motion about^bz(ψ ). The data highlighted with bold stands for the estimated hydrodynamic parameters that were considered valid as they fulfill the requirements discussed in Section 5.

Table

Table A1.Estimation of hydrodynamic parameters for the positive motion in^bxwith 120 iterations.

Table A1.Estimation of hydrodynamic parameters for the positive motion in^bxwith 120 iterations.

Case	1	2	3	4
mxkg	5.765440244	7.722155889	8.703062008	9.068444505
Xukg/s	0.394198185	0.087986129	-0.053646189	-0.119567686
Xuukg/m	4.102304926	5.591797761	6.326064633	6.70061939
Variance	0.018398838	0.009261478	0.00571343	0.006470697
Stnd Devmx	0.373086205	0.318102822	0.271556001	0.303454831
Stnd DevXu	0.346590739	0.248908439	0.19783262	0.212992615
Stnd DevXuu	1.782154508	1.282058946	1.022656101	1.106084553
Case	5	6	7	8
mxkg	7.64703005	9.21518067	9.385690628	9.566054714
Xukg/s	0.135204424	-0.052398199	-0.082941606	-0.123623731
Xuukg/m	5.361670536	6.225837822	6.430655651	6.684074805
Variance	0.008921403	0.006293474	0.00595281	0.006132446
Stnd Devmx	0.308064692	0.303419582	0.29964355	0.31056563
Stnd DevXu	0.243119075	0.205816942	0.202121283	0.207350628
Stnd DevXuu	1.250463857	1.059446073	1.044496664	1.076580065
Case	9	10	11
mxkg	7.64703005	8.391070197	9.206834393
Xukg/s	0.135204424	0.073448775	0.004971469
Xuukg/m	5.361670536	5.665845761	5.944058517
Variance	0.008921403	0.006211526	0.006128335
Stnd Devmx	0.308064692	0.277198866	0.301857089
Stnd DevXu	0.243119075	0.210389465	0.209432835
Stnd DevXuu	1.250463857	1.071664569	1.06616972

Table

Table A2.Results of the estimation of hydrodynamic parameters for the positive motion in^bxwith 150 iterations.

Table A2.Results of the estimation of hydrodynamic parameters for the positive motion in^bxwith 150 iterations.

Case	1	2	3	4
mxkg	5.117152947	7.637939366	8.401739147	8.810009611
Xukg/s	-0.688362153	-0.369721766	-0.266771427	-0.190308172
Xuukg/m	-2.615927203	-4.150283399	-4.677750899	-5.082018777
Variance	0.023813223	0.009599942	0.006820511	0.006529014
Stnd Devmx	0.363228922	0.286740198	0.258334352	0.264353731
Stnd DevXu	0.34730935	0.222358949	0.188752273	0.186091016
Stnd DevXuu	1.776793931	1.138640192	0.968753797	0.957926171
Case	5	6	7	8
mxkg	7.569930729	8.966171898	9.072344669	9.222521444
Xukg/s	-0.393675317	-0.256729862	-0.235866828	-0.186813331
Xuukg/m	-4.031800872	-4.672459037	-4.812534801	-5.084410331
Variance	0.009589836	0.006068915	0.005989863	0.006152794
Stnd Devmx	0.283968253	0.257984418	0.259231581	0.267621685
Stnd DevXu	0.221672771	0.177160111	0.176990795	0.180647372
Stnd DevXuu	1.134017318	0.906868497	0.908260782	0.929839948
Case	9	10	11
mxkg	7.569930729	8.305017728	9.13082585
Xukg/s	-0.393675317	-0.231001365	-0.142306865
Xuukg/m	-4.031800872	-4.7674745	-5.153776644
Variance	0.009589836	0.006973362	0.005825496
Stnd Devmx	0.283968253	0.260816606	0.259120036
Stnd DevXu	0.221672771	0.195171189	0.178821518
Stnd DevXuu	1.134017318	0.989958155	0.906535646

Table

Table A3.Results of the estimation of hydrodynamic parameters for the negative motion in^bxwith 150 iterations.

Table A3.Results of the estimation of hydrodynamic parameters for the negative motion in^bxwith 150 iterations.

Case	1	2	3	4
mxkg	4.490766485	7.476342027	8.518860195	9.237964526
Xukg/s	0.693665481	0.292592974	0.111440887	0.035208908
Xuukg/m	2.916991222	4.905349207	5.801569548	6.211827813
Variance	0.028267535	0.012738937	0.008886547	0.006633462
Stnd Devmx	0.37350401	0.3343355	0.304900289	0.279102261
Stnd DevXu	0.398094845	0.270160539	0.228000587	0.198653983
Stnd DevXuu	2.121582311	1.440839276	1.218244546	1.064598718
Case	5	6	7	8
mxkg	7.396818781	9.435143478	9.602678894	9.806286025
Xukg/s	0.290536196	0.049737835	-0.028148947	-0.013924468
Xuukg/m	4.901303762	6.069494533	6.449949306	6.422202462
Variance	0.01305013	0.006110739	0.006417925	0.006205586
Stnd Devmx	0.336010698	0.272154653	0.284726263	0.285326144
Stnd DevXu	0.272728409	0.187983118	0.194393045	0.192392557
Stnd DevXuu	1.452850112	1.001725205	1.037955359	1.030669723
Case	9	10	11
mxkg	7.396818781	8.447696759	9.432838931
Xukg/s	0.290536196	0.190128043	0.091036771
Xuukg/m	4.901303762	5.400394316	5.868182533
Variance	0.01305013	0.007716584	0.006124236
Stnd Devmx	0.336010698	0.280816445	0.274928891
Stnd DevXu	0.272728409	0.215873225	0.192786902
Stnd DevXuu	1.452850112	1.140846127	1.018376008

Table

Table A4.Results of the estimation of hydrodynamic parameters for the negative motion in^bzwith 150 iterations.

Table A4.Results of the estimation of hydrodynamic parameters for the negative motion in^bzwith 150 iterations.

Case	1	2	3	4
mzkg	2.292646188	5.696928677	8.155996872	9.469217358
Zwkg/s	6.823495905	7.115760391	7.339611974	7.476632539
ZgN	-13.8553377	-14.0199111	-14.145331	-14.2203024
Variance	0.449265112	0.29264543	0.193421042	0.168973397
Stnd Devmz	0.377118115	0.495312033	0.493076209	0.516926018
Stnd DevZw	0.279542905	0.228175279	0.18747121	0.176925072
Stnd DevZg	0.138675692	0.113488356	0.093421496	0.0882791
Case	5	6	7	8
mzkg	5.938688654	10.10754985	10.76495451	11.4083869
Zwkg/s	7.130143256	7.454997333	7.543491172	7.632676342
ZgN	-14.0309875	-14.2150142	-14.2626364	-14.3098727
Variance	0.274305565	0.149793176	0.131641757	0.12178331
Stnd Devmz	0.48346965	0.506240918	0.496024952	0.499068899
Stnd DevZw	0.22076476	0.165605526	0.156247509	0.151410383
Stnd DevZg	0.109874512	0.082761713	0.078125774	0.075744148
Case	9	10	11
mzkg	5.938688654	7.039846207	9.477346809
Zwkg/s	7.130143256	6.671642701	6.75833619
ZgN	-14.0309875	-13.7838234	-13.8376442
Variance	0.274305565	0.326174877	0.280832183
Stnd Devmz	0.48346965	0.631375399	0.732671098
Stnd DevZw	0.22076476	0.23520344	0.218829433
Stnd DevZg	0.109874512	0.117923137	0.109872371

Table

Table A5.Results of the estimation of hydrodynamic parameters for the positive motion in^bzwith 150 iterations.

Table A5.Results of the estimation of hydrodynamic parameters for the positive motion in^bzwith 150 iterations.

Case	1	2	3	4
mzkg	1.051614679	4.752953968	6.587469056	8.306230779
Zwkg/s	5.293985579	6.066173258	6.397884543	6.798312912
ZgN	-13.3074485	-12.9363122	-12.7760468	-12.5819159
Variance	0.324015369	0.220089	0.180641832	0.142403005
Stnd Devmz	0.314127547	0.530858752	0.578877071	0.584207504
Stnd DevZw	0.362108568	0.312915507	0.289395855	0.264295339
Stnd DevZg	0.158922394	0.138398903	0.128488759	0.117945516
Case	5	6	7	8
mzkg	4.990263133	8.906760485	9.787928045	10.08830879
Zwkg/s	6.094918316	6.850862097	6.994152502	7.157339123
ZgN	-12.9242288	-12.5459698	-12.480379	-14.309872
Variance	0.213614138	0.137502744	0.117958461	0.110153575
Stnd Devmz	0.534406829	0.605189466	0.588705344	0.578252497
Stnd DevZw	0.308440059	0.260545523	0.243063523	0.23836489
Stnd DevZg	0.136373993	0.116724808	0.108895784	0.107017636
Case	9	10	11
mzkg	4.990263133	5.767758094	8.199786253
Zwkg/s	6.094918316	4.877488253	4.995036925
ZgN	-12.9242288	-13.4829903	-13.4206186
Variance	0.213614138	0.291733449	0.262690659
Stnd Devmz	0.534406829	0.821382471	0.976860702
Stnd DevZw	0.308440059	0.322705147	0.307549992
Stnd DevZg	0.136373993	0.142323817	0.13591481

Table

Table A6.Results of the estimation of hydrodynamic parameters for the motion inψwith 150 iterations.

Table A6.Results of the estimation of hydrodynamic parameters for the motion inψwith 150 iterations.

Case	1	2	3	4
Izzkgm2	0.09150314	0.154938563	0.205543484	0.233317348
Nrkg/s	0.01451316	0.014305644	0.014146321	0.014072976
Variance	0.000415892	0.000331019	0.000266467	0.000238701
Stnd DevIzz	0.014166958	0.016301838	0.016890015	0.017264417
Stnd DevNr	0.007894602	0.007063666	0.006365152	0.006056496
Case	5	6	7	8
Izzkgm2	0.152784249	0.24495803	0.258931023	0.274659776
Nrkg/s	0.014292885	0.01387251	0.013782836	0.01379101
Variance	0.000331061	0.000224127	0.000208268	0.000192268
Stnd DevIzz	0.016080323	0.017150475	0.017042359	0.016952814
Stnd DevNr	0.007051294	0.005810967	0.00562722	0.005435597
Case	9	10	11
Izzkgm2	0.152784249	0.181382573	0.243559836
Nrkg/s	0.014292885	0.01718263	0.017010375
Variance	0.000331061	0.000299337	0.000228551
Stnd DevIzz	0.016080323	0.016875525	0.017346636
Stnd DevNr	0.007051294	0.006720781	0.005872666