Method for verifying the swell stability of a buoyancy system associated with a means of transport

Abstract

The invention relates to a computer-implemented method for verifying the swell stability of a buoyancy system (2) comprising at least one float (3) intended to ensure the buoyancy of an associated means of transport (1), comprising steps of: - two-dimensional digital finite-element simulation of a two-dimensional fluid domain comprising at least one two-dimensional liquid domain representative of an irregular swell of given spectrum, - selection of critical swell times, - verification of swell stability via three-dimensional digital finite-element simulation of a three-dimensional model representative of the float (3) and of the associated means (1) in a three-dimensional fluid domain comprising at least one three-dimensional liquid domain representative of a selected critical swell time.

Description

[0001] DESCRIPTION

[0002] TITLE: Method for verifying the wave stability of a buoyancy system associated with a transport device

[0003] Technical field

[0004] The invention relates to the technical field of buoyancy systems enabling the landing and afloat stability of a transport device, more particularly an aircraft, and most particularly a rotary-wing aircraft such as a helicopter.

[0005] Previous techniques

[0006] Any aircraft intended to transport people in maritime areas must be equipped with a buoyancy system which, in the event of an emergency landing, provides sufficient flotation to allow the aircraft to be evacuated.

[0007] Certification regulations specify that an aircraft must be able to land and be stable on water with its buoyancy system. Stability must be proven for free water surface conditions that are defined in these certification regulations through parameters such as wave height and period.

[0008] These free surface states of water, also called "sea conditions", also referred to as "sea state" in English, apply to any liquid surface. The term "landing" also covers the case of an aircraft landing on any free surface of water, whether on the sea or a lake for example.

[0009] Such a buoyancy system may include floats whose deployment is either controlled or triggered automatically by immersion detection. These floats may, for example, be of the textile type and inflatable by an explosive or electrical deployment means.

[0010] Currently, the justification of compliance with normative requirements is carried out mainly by physical tests in basins using reduced-scale models.

[0011] Such physical tests have several disadvantages. For example, the small size of the test tanks existing in the world requires the use of reduced models for which the application of the rules of similarity is not always easy. This is particularly the case for inflatable textile floats, the reduced models of which do not faithfully represent the behavior of the full-size model which is deformable, because the thickness of the fabric cannot be reduced appropriately, which in practice amounts to considering the models of these floats to be non-deformable.

[0012] Moreover, such physical tests are expensive due to the technical and human resources required for their implementation, and the availability of test pools is not always guaranteed, which severely limits optimization by successive iterations.

[0013] To overcome the above-mentioned drawbacks of physical tests, we can use numerical simulations.

[0014] However, numerical simulations have specific drawbacks. On the one hand, numerical simulation software that allows the study of deformable structures of the inflatable textile type does not have native capabilities for generating wave-type loads of a given spectrum. On the other hand, numerical simulation software that allows the generation of wave-type loads of a given spectrum does not have native capabilities for studying deformable structures of the inflatable textile type. There is thus no existing software solution allowing the study of deformable structures of the inflatable textile type subjected to natively generated wave-type load cases of a given spectrum. In addition, the numerical simulation of a large fluid domain and of a long duration can require long computation times that are incompatible with the deadlines accepted for industrial developments.

[0015] Statement of the invention

[0016] In view of the above, the aim of the invention is therefore to propose an improved simulation method, allowing the study of deformable structures, in particular of the inflatable textile type, subjected to load cases of the wave type with a given spectrum generated natively.

[0017] The subject of the invention is a computer-implemented method for verifying the stability in waves of a buoyancy system comprising at least one float intended to ensure the buoyancy of an associated transport device.

[0018] The method comprises steps of: two-dimensional numerical simulation by finite elements of a two-dimensional fluid domain comprising at least one two-dimensional liquid domain representative of an irregular wave of given spectrum, selection of critical wave instants according to the results of the two-dimensional numerical simulation, verification of wave stability by three-dimensional numerical simulation by finite elements of a three-dimensional model representative of the float and the associated apparatus in a three-dimensional fluid domain comprising at least one three-dimensional liquid domain representative of a selected critical wave instant, the initial relative positioning between the three-dimensional model and the three-dimensional fluid domain being carried out according to a master node of the three-dimensional model and at least one reference node of the fluid domain.

[0019] Advantageously, the selected critical wave instant is directly translated into the three-dimensional liquid domain by initial conditions of velocity, free surface and pressure resulting from the two-dimensional numerical simulation step.

[0020] Compared to a conventional three-dimensional numerical simulation implemented by computer, the claimed method requires shorter calculation times and a smaller storage space, thus obtaining an economy of calculation and storage means which exceed the normal physical interaction between a computer program and a computer system.

[0021] Preferably, the initial relative positioning is performed without relative velocity between the master node of the three-dimensional model and the fluid at the reference node of the fluid domain.

[0022] According to one feature, the associated transport apparatus is an aircraft.

[0023] According to another feature, the floats are made of an inflatable textile-type material.

[0024] For example, the given spectrum is a JONSWAP type wave spectrum.

[0025] Advantageously, the swell of given spectrum is obtained from the superposition of elementary waves each having a random phase shift calculated according to an algorithm not introducing a link between the phase shift of a wave and its pulsation.

[0026] Advantageously, the numerical simulation steps use the first order of Stokes theory for irregular swell.

[0027] Preferably, the numerical simulation steps linearly interpolate velocities. Advantageously, the selection of a critical wave instant is carried out as a function of a predetermined threshold value of the maximum slope of a free wave surface or of the maximum slope variation of the free wave surface.

[0028] Preferably, the method comprises an additional step of physical certification testing carried out after the step of verifying wave stability by three-dimensional numerical simulation.

[0029] According to another aspect, the invention relates to a buoyancy system whose stability in waves is verified by a method as described above.

[0030] The invention also relates to a transport apparatus comprising a buoyancy system as described above.

[0031] Brief description of the drawings

[0032] Other aims, characteristics and advantages of the invention will appear on reading the following description, given solely by way of non-limiting example, and made with reference to the appended drawings in which:

[0033] [Fig 1] represents a view of an aircraft equipped with a buoyancy system comprising several floats;

[0034] [Fig 2] illustrates a general flowchart of a method according to the invention;

[0035] [Fig 3] represents the graph of a JONSWAP type wave spectrum;

[0036] [Fig 4] illustrates a fluid domain used in a two-dimensional numerical simulation according to the invention;

[0037] [Fig 5] illustrates the general flowchart of calculation step 200c;

[0038] [Fig 6] illustrates the flowchart of step 218c of calculating the speeds; [Fig 7] schematically represents a fluid domain used in a three-dimensional numerical simulation according to the invention; and

[0039] [Fig 8] illustrates the general flowchart of calculation step 400c.

[0040] Detailed description of at least one embodiment

[0041] Figure 1 shows an aircraft 1 according to the invention. This aircraft 1 is provided with a buoyancy system 2 according to the invention in order to be able to float on the water in the event of a water landing. Such a buoyancy system 2 is provided with at least one float, and for example at least two floats 3, arranged on either side of the fuselage of the aircraft 1.

[0042] The floats may be paired. Therefore, the floats of a pair may be arranged symmetrically on either side of an anteroposterior plane of symmetry of the aircraft.

[0043] In the illustrated example, the floats are inflatable floats. Outside of the landing phases, each float 3 can be folded into a compartment of the aircraft, and in the landing phase each float is inflated by a deployment device (not shown) of the buoyancy system 2.

[0044] We will now describe, with reference to figures 2 to 6, a method 100 for verifying the stability to swell of a buoyancy system 2 according to the invention.

[0045] The method 100 allows the verification of the wave stability of any buoyancy system comprising at least one float subjected to an irregular wave of given spectrum. The method 100 can find applications relating to the evaluation of the wave stability of inflatable boats or aircraft in the event of a water landing.

[0046] Figure 2 illustrates the general flowchart of the process 100. The natural swell appears in the form of irregular undulations in height and period. To describe it, it is convenient to use spectral analysis which decomposes the evolution of the free surface into a superposition of sinusoidal waves, and to deduce the distribution of energy in frequency. This involves performing a transformation from the time domain to the frequency domain, for example by applying a Fourier transform (also called by the acronym "FFT" for Fast Fourier Transform in English) to a time series of free surface level at a point. There are thus several definitions of a swell spectrum such as the JONSWAP spectrum or the Pierson-Moskowitz spectrum.

[0047] For example, the following expression corresponds to the spectrum

[0048] JONSWAP (Figure 3), used in aircraft certification standards to justify wave resistance. (Eq. 1) where: S(w) is the energy spectral density;

[0049] H s is the specific height; w p is the peak pulsation; w is the pulsation; y is the shape coefficient, in particular equal to 3.3 in the standards; a is a parameter defined by the relation where a is a constant.

[0050] Depending on the expression chosen for the spectrum, the method 100 begins with a step 101 of determining the regular waves whose superposition respects the given wave spectrum (figure 2).

[0051] The spectrum is divided into a finite number of intervals of predefined width Aw corresponding to a pulsation increment. Each interval corresponds to a regular wave of pulsation w centered on the interval and of height H according to the following equation: where: Aw is the sampling step of the spectrum.

[0052] In order to avoid a phase superposition of elementary waves, namely regular waves of order 1, a random phase shift cp between -7i and +7i is chosen for each pair H and w defining a regular wave.

[0053] It is mandatory that the random number generation algorithm does not introduce a link between the phase shift and the pulsation of the wave under consideration. An example of an algorithm meeting this criterion is the linear congruential algorithm in which each random number X is calculated according to the following recursive equation:

[0054] X n+1 = a * X n + b~) * mod m (Eq.3) with, in particular Xo=l , a= 1664525, b= 1013904223 and m=2 32 .

[0055] The method 100 continues with a step 200 of numerical simulation by finite elements of a two-dimensional fluid domain. The two-dimensional fluid domain comprises at least one liquid domain representative of an irregular swell of given spectrum.

[0056] Generally speaking, a numerical simulation makes it possible to account for physical phenomena by means of mathematical models. Step 200 of the method 100 consists of numerically simulating a fluid domain comprising a liquid domain subjected to an irregular swell of given spectrum, through three steps, namely: a modeling step 200a, a meshing step 200b and a calculation step 200c. It should be noted that step 200 concerns a simulation of fluid dynamics (also referred to by the acronym "CFD" for Computational Fluid Dynamics in English), consisting of solving the equations governing the flow of a liquid otherwise described by its behavior law and the volumes in which it flows.

[0057] Modeling step 200a thus associates the physical laws with a geometry created by computer-aided design (CAD). This step 200a includes the definition of the domain to be studied, namely: the geometry, the physical characteristics and the conditions applicable to the boundaries of the domain, called boundary conditions.

[0058] In Figure 4, a fluid domain 10 of contour ABCDEF is shown, which in the illustrated example corresponds to a rectangle. The fluid domain 10 comprises a large set of elementary volumes forming a two-dimensional fluid domain of unit thickness and oriented along the XZ plane. The domain 10 is of predetermined length and height, the length being arranged in the direction of flow and oriented along the X axis. Preferably, the length of the domain 10 is greater than the height of the domain 10. The height of the domain 10 measured along the Z axis is constant over the entire length of the domain 10. The length of the domain 10 depends on the wavelength X of the modeled waves.

[0059] The domain 10 comprises a first wave zone 11 corresponding to the contour ABEF and a second zone 12 called damping zone corresponding to the contour BCDE. The wave zone 11 makes it possible to simulate the wave without imposing additional damping. The damping zone 12 makes it possible, by applying viscous forces, to attenuate the energy of the fluid and thus to make the effects of wave reflection at the edges negligible. These viscous forces are proportional to the speed of the fluid and to a damping coefficient.

[0060] The swell zone 11 is in principle long enough to allow the formation of breaking waves and preferably longer than the length of the damping zone 12.

[0061] The fluid domain 10 comprises a liquid domain 14 filled with water. The water height is defined by a surface 15 called the free surface. Any wave is initialized as a variation around a mean depth using a free surface equation applied to the swell zone 11. In the damping zone 12, linear damping is applied. At the beginning of the damping zone, the free surface is not modified and at the end of the zone, the elevation relative to the mean height is zero.

[0062] The fluid domain 10 may also comprise a gaseous domain 13, complementary to the liquid domain 14. The gaseous domain 13 is filled with air at rest or having speeds representative of a given wind speed.

[0063] The boundary conditions are applied to the nodes present on the contours AF and FED on the one hand and on the contour CD on the other hand. On the contours AF and FED, velocities calculated according to Fenton's theory are applied. To limit the number of elements, the liquid domain 14 does not respect the deep sea condition, that is to say that the velocities present at the bottom of the liquid domain 14 are not zero, but calculated so that the velocities at the bottom of the sea are. The bottom of the liquid domain 14 thus corresponds to the lower edge of a fictitious slice of a deep sea. For the nodes of the liquid domain 14 located in the damping zone 12, the velocities calculated according to Fenton's theory are multiplied by a decreasing linear coefficient ranging from 1 to 0 for the nodes corresponding to points E and D respectively.For the nodes located at the exit of domain 10 on the contour CD and located below the free surface and for the nodes on the contour ED of the last two elements near D, a hydrostatic pressure corresponding to the height of water present is applied in order to preserve the mass of water in domain 10. In the case where the simulation takes into account the air present above the free surface, a hydrostatic pressure corresponding to the air pressure is applied to the nodes located on the contour CD and above the free surface.

[0064] In the case where the fluid domain 10 comprises a gaseous domain 13, “non-reflective” type edge conditions applied to the upper edge of the domain 13 are necessary to have stabilized air masses. In the case of a gaseous domain 13 containing air at rest, that is to say where no wind speed is imposed, a “non-reflective” type edge condition is additionally added to the nodes of the contour CD located above the free surface, in order to allow the stabilization of the air without hindering the evacuation of the water.

[0065] Referring again to Figure 2, after the modeling step 200a, the method 100 continues with the step 200b of meshing the domain 10. The meshing corresponds to the discretization of the domain into finite elements. The mesh size is chosen to be finer in the wave zone 11 compared to the damping zone 12. For example, the average mesh size for the wave zone 11 is 300mm x 300mm x 300mm. The mesh size of the elements of the damping zone 12 can vary from a size substantially equal to the mesh size used for the wave zone 11 to a size several times larger in order to simplify the calculations without losing precision.

[0066] After the meshing step 200b, the method 100 continues with the calculation step 200c. Figure 5 represents the general flowchart of the calculation step 200c.

[0067] The calculation step 200c includes a resolution step by a solver 201c which implements algorithms for assembling, storing and matrix calculating all the equations to be solved for the wave simulation and the associated results. The calculated and stored parameters include for each element the fraction of volume filled by the different fluids and make it possible to determine the water level defining the free surface as a function of time. The calculated and stored parameters also include the velocities and pressures at the nodes of the mesh as a function of time.

[0068] The calculation step 200c includes a speed calculation step 210c which is translated computationally into a procedure-type routine which accepts several arguments and which performs several calculation and storage tasks which will be described in more detail below.

[0069] Step 210c of calculating speeds begins with step 212c of declaring local and global variables.

[0070] For example, a global variable of dimension 1 stores all the arguments passed to the routine, such as for example: the wave height Hs, the peak period Tp, the order of the Stokes wave theory equations, the wave definition time, the start time, the start of the damping region, the end of the damping region, the mesh size, the interpolation times (specified below), and the wind speed, if any. A global or persistent variable is a variable whose values ​​are not lost between two solver cycles.

[0071] A persistent variable, of type list containing data structures, called JONSWAP, is also created to contain all the characteristics of the regular waves that superimpose the irregular swell. This variable includes for each wave: the height, the period, the phase shift, the parameters according to Fenton's theory, etc.

[0072] The type of persistent variable used is compatible with parallel resolution.

[0073] It should be noted that it is not necessary to calculate the speeds precisely at all time steps, but it is sufficient to calculate the speeds precisely for two time values ​​called interpolation times and noted tl and t2, and to obtain by linear interpolation the values ​​of the intermediate speeds between tl and t2. This interpolation method greatly lightens the computational load and allows correct results if the time interval between tl and t2 is small enough to admit the approximation by linearization of the speeds.

[0074] To allow interpolation for all time steps between tl and t2, the velocities at each node are stored in a global variable of type structure for tl and t2.

[0075] After the variable creation step 212c, the method 100 continues with an initialization step 214c followed by a reading step 216c, making the variables declared in step 212c available for subsequent calculations.

[0076] The declaration of variables is made at each solver cycle, that is to say at each call to step 210c of velocity calculation. The initialization step 214c is used to initialize, that is to say to calculate and store, the persistent variables. Once the variables are initialized, they are accessible throughout the duration of the calculation without having to be reinitialized at each call to step 210c of velocity calculation, unlike "normal" variables whose values ​​are lost between two calls.

[0077] However, at each call to step 210c of calculating the speeds, it remains necessary to read the persistent variables to be able to access the values ​​(reading step 216c).

[0078] In the next speed calculation step 218c, the speed calculation is done when the current time is greater than the interpolation time t2 of the last node of the calculated block.

[0079] If this is the case, times tl and t2 are updated at step 1001 c: tl takes the value of t2 and t2 is increased by adding an interpolation time step corresponding to t2-tl (figure 6).

[0080] In the next reading step 1002c the JONSWAP variable is read.

[0081] In the next memory allocation step 1003c, the routine allocates memory for a list whose size is the number of regular waves contained in the variable JONSWAP. This list will contain the free surface height of each wave at t1 and t2.

[0082] In step 1004c, the coordinates of the damping zone are updated from the movements of the last node of the block, implicitly assuming that all the nodes have moved in the same way. This update of the coordinates is necessary to allow the domain 10 to follow the drift of a floating object if necessary. The ability of the domain 10 to follow the drift of a floating object makes it possible to reduce the modeled size of the domain 10 and to greatly reduce the calculations.

[0083] The instructions are then executed for all nodes of the block to be calculated (step 1005c).

[0084] In the next update step 1006c, the times tl and t2 of the current node are updated.

[0085] In the next copying step 1007c, the values ​​corresponding to the tl of the velocities, free surface height and immersion state are replaced by the values ​​corresponding to t2. The immersion state can be expressed by numerical values, for example 1 for the immersed state, 2 for a free non-immersed state and 3 for the non-immersed state subject to the wind.

[0086] In step 1008c of calculating the free surface, the free surface of each regular wave is first calculated and then added to the total free surface.

[0087] In the next update step 1009c, the immersion state of the node is updated based on the calculated total free surface.

[0088] The method 100 continues with step 1010c of verifying the immersion state of the node. A node is considered submerged when it is at a level lower than that of the free surface raised by one mesh size. A node is considered located in the wind zone if it is at a level higher than that of the free surface raised by four mesh sizes. A node is considered located in the free zone when it is at a level higher than the free surface raised by one mesh size and lower than the free surface raised by four mesh sizes.

[0089] In the following velocity calculation step 101 1 c, the node's immersion states are taken into account to decide whether wave velocities should be calculated. If, for tl and t2, the node is in the wind zone, then wind velocities are calculated. If, at tl or t2, the node is considered submerged, wave velocities are calculated. If the node is considered to be located in the free zone, the velocity will not be calculated by the routine, but directly by the solver in accordance with the physical modeling used.

[0090] The speed v n of a regular wave n at a given height z is given by the following equation: where: v n is the speed at the node with coordinates x,z of the elementary wave n;

[0091] H n is the height of the elementary wave n; k nis the wave number of the elementary wave n; x is the x coordinate at which the speed is calculated; z is the z coordinate at which the speed is calculated; and (p n is the phase shift of the elementary wave n.

[0092] Equation 4 is a simplified equation of order 1 considering the assumptions of deep water and at t=0.

[0093] However, it makes no sense to calculate the velocity brought by a regular wave above its free surface, because the equation for the velocities of a regular wave is defined only for points below its free surface.

[0094] Therefore, it is necessary to use a method of stretching the speeds in order to avoid regular wave speeds being calculated above the free wave surface.

[0095] For example, one method of recalibrating velocities in a node with coordinates XG, ZG consists of replacing the value z in equation 4 with a recalibrated value z n , given by the following equation: z n = z G - z S LT - z SLN ) (Eq.5) where: ZG is the z coordinate of the node, ZSLT is the z coordinate of a point with coordinate XG and belonging to the free surface of the irregular wave, and, ZSLN is the z coordinate of a point with coordinate G and belonging to the free surface of the regular wave n.

[0096] The z-coordinate recalibration ensures that the velocity on the total free surface corresponds to the sum of the velocities at the free surfaces of each regular wave. The velocity on the total free surface is thus calculated by summing the velocity component at the free surface of each regular wave.

[0097] In the next extrapolation step 1012c, a linear extrapolation is performed to calculate the speed of the nodes located at a level lower than or equal to the free surface raised by three meshes. The speed of these nodes is obtained by linear extrapolation based on the known speed values ​​at the free surface and at the node located just below the free surface and the ratio of the distance measured on the z axis.

[0098] The extrapolation along the z coordinate is calculated using the following equations: where v x is the speed along x; v z is the speed along z; z(...) are the z coordinates used for the extrapolation; zgriiie is the a z coordinate of speed evaluation;

[0099] Coeff a i is the damping coefficient.

[0100] Furthermore, 1 represents the node below the surface, and 2 represents the node on the free surface. The damping coefficient Coeff am depends on the position of the node. The value of the damping coefficient Coeff am varies linearly between

[0101] 1 at the start of the depreciation zone 12 and

[0102] 0 at the end of the depreciation zone 12.

[0103] Outside the damping zone, the coefficient is equal to 1.

[0104] After the extrapolation step 1012c, the interpolation coefficients are calculated and stored in global variables (step 1013c).

[0105] The coefficients are calculated according to the following equations:

[0106] Method 100 continues with step 1014c of interpolating the speeds.

[0107] If the node is considered immersed for tl or t2, the velocities are calculated for each current time te using the interpolation coefficients A and B according to the following equation: v(node, te) = A * tc + B (Eq.10)

[0108] If the node is in a wind zone for tl and t2 and a non-zero wind speed has been given as an argument to the routine, the wind speed along the x axis is calculated from the wind speed given as an argument and the associated amplitude value. Note that in this case no wind speed along the z axis is calculated. In the case where there is no wind speed given as an argument or if this value is zero, no speed is calculated.

[0109] If the node is in a free zone, no speed is calculated.

[0110] Referring again to FIG. 2, the method 100 continues with a step 300 of selecting critical wave instants based on the results of the simulation step 200. A critical wave instant is characterized by a maximum height of a wave, by a time and by a location at which this wave forms, as well as a maximum slope value of the free surface of the wave and / or a maximum variation value of the slope of the free surface of the wave. Such an instant is considered critical with respect to the overturning of a device subjected to this wave. The location at which a critical wave forms is stored in the form of coordinates of a reference node of the liquid domain.

[0111] The criterion for selecting a critical instant may be based on threshold values ​​for the maximum slope of the free surface of the wave or for the variation of the maximum slope of the free surface. Alternatively, the selection step may comprise a sub-step of visualizing the results of the simulation step 200, allowing a user to identify and select critical instants where the waves exhibit breaking-type instabilities.

[0112] It should be noted that the selection of a critical instant is accompanied by the recording of the speed, free surface and pressure values ​​corresponding to this critical instant.

[0113] After step 300 of selecting a critical instant, a second numerical simulation by finite elements makes it possible to verify the stability to wave of a three-dimensional model of at least one float and an associated device in a three-dimensional liquid domain representative of the selected critical instant (step 400).

[0114] The method 100 thus performs a three-dimensional numerical simulation from the results obtained by a prior two-dimensional numerical simulation. This approach has a considerable advantage over a single three-dimensional numerical simulation. Indeed, the prior two-dimensional numerical simulation is carried out over long times which are necessary to obtain the most critical waves of a swell of given spectrum. The three-dimensional numerical simulation which uses the results of the two-dimensional numerical simulation can thus be shortened and concentrated around the identified critical instants, unlike a single three-dimensional numerical simulation which must be carried out over a sufficiently long time for the waves to be able to form.Thus, being able to find the breaking shape of a wave in a three-dimensional numerical simulation without having to simulate a long period before this wave allows considerable computational savings.

[0115] Similarly to the numerical simulation step 200, the step 400 comprises three constituent steps of a numerical simulation, namely: a modeling step 400a, a meshing step 400b and a calculation step 400c.

[0116] The modeling step 400a comprises the definition of a three-dimensional fluid domain 10' shown schematically in FIG. 7 and similar to the fluid domain 10 defined in step 200a. However, unlike the two-dimensional fluid domain 10, the fluid domain 10' has a constant thickness along y comprising a succession of several elementary volumes.

[0117] In a similar manner to the two-dimensional fluid domain 10, the three-dimensional domain 10' comprises a first wave zone 11' and a second damping zone 12'. The three-dimensional fluid domain 10' comprises a liquid domain 14' filled with water and may also comprise a gaseous domain 13', complementary to the liquid domain 14'. The water height is defined by a surface 15' called the free surface.

[0118] It should be noted that the length of the 1 1 ' wave zone can be reduced compared to the length of the 1 1 ' wave zone of the two-dimensional model taking into account the fact that the critical wave instants are introduced directly as initial conditions of the three-dimensional fluid model 10' and that it is not necessary to simulate their formation.

[0119] The boundary conditions of the two-dimensional fluid domain 10 apply analogously to the three-dimensional fluid domain 10'.

[0120] Since waves are unidirectional and have a propagation axis oriented along the X axis, it can be considered that there is no variation in speed along the Y axis upstream of any floating object.

[0121] The modeling step 400a also comprises the definition of a three-dimensional structural model 16' representative of at least one float 3 and an associated transport device. The three-dimensional structural model 16' defines the geometry and the physical characteristics of the elements that compose it. It is thus easy, for example, to define elements with complex behavior by defining, for a given material, behavior laws or mechanical characteristics such as inertia and / or rigidity as a function of the local axes of an element.

[0122] It should be noted that the three-dimensional 16' model uses a Lagrangian description allowing the study of the parameters of the structure at a point which follows its displacement. The 10' fluid domain follows an Eulerian description allowing the study of the flow parameters at a fixed point in space. The coupled Eulerian-Lagrangian analysis takes advantage of the advantages of both descriptions, while limiting their disadvantages.

[0123] The relative positioning between the three-dimensional structure model 16' and the three-dimensional fluid domain 10' is carried out using a master node 17' of the structure model 16' which may correspond, for example, to the center of gravity. The relative positioning of the master node 17' and at least one reference node of the fluid domain 10' is carried out so as to maximize the probability of overturning of the structure model relative to the critical wave identified in the selection step 300. For example, the orientation of the three-dimensional structure model 16' may follow the initial angle of the wave along the Y axis and the anteroposterior plane of symmetry of the model may be oriented perpendicular to the direction of propagation of the waves.

[0124] Alternatively, it is possible to define in advance the three-dimensional structure model 16' and then proceed to position the fluid domain 10' around the structure model 16' by imposing initial conditions of speed, pressure and free surface corresponding to the critical wave instant selected in step 300.

[0125] It is noted that the selected critical wave instant is directly applied to the three-dimensional liquid domain through initial conditions of velocity, free surface and pressure resulting from the two-dimensional numerical simulation step and corresponding to the selected critical wave instant.

[0126] After the modeling step 400a, the method 100 continues with the meshing step 400b. The mesh is refined locally around the areas of interest in order to increase the accuracy of the results in these areas, without increasing the computations in the other areas.

[0127] After the meshing step 400b, the method 100 continues with the calculation step 400c.

[0128] Calculation step 400c will not be detailed because it is analogous to calculation step 200c which has already been described. A general flowchart of step 400c is presented in Figure 8, in which the steps identical to those of step 200c bear the same references. Calculation step 400c includes a speed calculation step 410c.

[0129] Step 410c of calculating speeds begins with step 212c of declaring local and global variables.

[0130] After the variable creation step 212c, the method 100 continues with an initialization step 414c.

[0131] Thanks to the unidirectional wave hypothesis, it is possible to neglect the variation of wave parameters along the Y axis. Therefore, speed calculations (including interpolation for time steps shorter than the interpolation time) are only carried out with respect to the x and z coordinates.

[0132] In the initialization step 414c, a node sorting operation is performed based on the x and z coordinates. Each unique pair of coordinates (x,z) is identified by a unique identification number called ID. When sorting the nodes, a unique ID is assigned to each unique pair (x,z). At the end of the sorting, the nodes having the same coordinates (x,z) are associated with the same ID. The calculations are performed based on the ID and not for each node, which significantly reduces the computational load.

[0133] In the next read step 216c, the variables declared in step 212c are read and stored for subsequent calculations.

[0134] In the next step 416c, the number of IDs is calculated based on the size of the block to be calculated.

[0135] In the next calculation step 418c, the speeds are calculated and interpolated for each ID, identically to the calculations performed in step 218c.

[0136] The numerical simulation carried out in step 400 may have a stopping criterion relating to a critical moment of reversal, for example a rotation of the three-dimensional structural model 16' by 90°, or in the absence of reversal a stopping criterion relating to a predetermined maximum duration of the simulation, for example a value chosen between three and ten seconds.

[0137] If the critical tipping moment is reached, it may be chosen to improve the buoyancy by modifying the geometry of the floats, for example, by increasing their size or by changing their shape or position. A new simulation step 400 is then carried out with a three-dimensional structure model 16' taking these modifications into account.

[0138] If the critical moment of overturning has not been reached, the design of the buoyancy system, characterized in particular by the number, size, shape, and position of the floats, is verified by numerical simulation with regard to stability in the face of an irregular swell of a given spectrum. A physical test in a tank then makes it possible to certify the buoyancy system. The method 100 thus comprises an additional step of physical testing for certification of the buoyancy system and the associated apparatus, carried out after the step of verifying stability in swell by three-dimensional numerical simulation.

[0139] The proposed numerical simulation method significantly reduces the number of physical tests required to design a compliant device, because the design and implicitly the multiple iterative optimization phases are validated by numerical simulation without having to resort to physical model tests. As soon as the wave stability is validated by three-dimensional numerical simulation, a single physical test is sufficient to prove the compliance with the standards of the buoyancy system and the associated device.

Claims

CLAIMS 1. Computer-implemented method (100) for verifying the stability in waves of a buoyancy system (2) comprising at least one float (3) intended to ensure the buoyancy of an associated transport apparatus (1), characterized in that it comprises steps of: - two-dimensional numerical simulation by finite elements of a two-dimensional fluid domain (10) comprising at least one two-dimensional liquid domain (14) representative of an irregular swell of given spectrum, - selection of critical wave moments based on the results of the two-dimensional numerical simulation, - verification of wave stability by three-dimensional numerical simulation using finite elements of a three-dimensional model (16') representative of the float (3) and the associated apparatus (1) in a three-dimensional fluid domain (10') comprising at least one three-dimensional liquid domain (14') representative of a selected critical wave instant, the initial relative positioning between the three-dimensional model (16') and the three-dimensional fluid domain (10') being carried out as a function of a master node (17') of the three-dimensional model (16') and at least one reference node of the fluid domain (10').

2. Method (100) according to claim 1, in which the selected critical wave instant is directly translated into the three-dimensional liquid domain by initial conditions of speed, free surface and pressure resulting from the two-dimensional numerical simulation step.

3. Method (100) according to claim 1 or 2, in which the initial relative positioning is carried out without relative velocity between the master node of the three-dimensional model and the reference node of the fluid domain.

4. Method (100) according to any one of claims 1 to 3, in which the associated transport apparatus is an aircraft.

5. Method (100) according to any one of claims 1 to 4, in which the floats are made of an inflatable textile-type material.

6. Method (100) according to any one of claims 1 to 5, in which the given spectrum is a JONSWAP type wave spectrum.

7. Method (100) according to any one of claims 1 to 6, in which the swell of given spectrum is obtained from the superposition of elementary waves each having a random phase shift calculated according to an algorithm not introducing a link between the phase shift of a wave and its pulsation.

8. Method (100) according to any one of claims 1 to 7, in which the numerical simulation steps use the first order of Stokes theory for irregular swell.

9. Method (100) according to any one of claims 1 to 8, in which the digital simulation steps linearly interpolate speeds.

10. Method (100) according to any one of claims 1 to 9, wherein the selection of a critical wave instant is carried out as a function of a predetermined threshold value of the maximum slope of a free wave surface or of the maximum slope variation of the free wave surface.

11. Method (100) according to any one of claims 1 to 10, comprising an additional step of physical certification testing carried out after the step of verifying stability in waves by three-dimensional digital simulation.

12. Buoyancy system in which wave stability is verified by a method (100) according to any one of claims 1 to 11.

13. A transport apparatus comprising a buoyancy system according to claim 12.

Images