Method of simulation of unsteady aerodynamic loads on an external aircraft structure

09767229 · 2017-09-19

Assignee

Inventors

Cpc classification

International classification

Abstract

The invention relates to a method for simulating the unsteady aerodynamic loads being exerted on the external structure of an aircraft, notably in the context of a simulation of the buffeting of a wing surface in an airflow. The method includes a step of measurement of pressure (410) at different points of a grid, a step of calculation of the spectral density at these points (420) followed by extrapolation/interpolation operations to calculate the missing measurements (430), a step of estimation of the pressure coherence (440) for each pair of points of the grid, a step of estimation of the interspectral pressure density (450) for these same pairs of points, from the coherence thus estimated, and a step of calculation of the aerodynamic loads (460) by summing the real part of the interspectral density for the area of the wing surface having a separation of the boundary layer.

Claims

1. A method of simulation by computer of unsteady aerodynamic loads being exerted on an external structure of an aircraft, positioned in an airflow, the method comprising: (a) a step of measurement of pressure by multiple sensors located on the surface of the external structure, and positioned at first plurality of points of a grid having a first direction in a direction of the airflow, and a second direction, separate from the first direction, where the pressure signals of the sensors are stored in a memory; (b) a step of calculation of a spectral density of the pressure signals at the first plurality of points of the grid as a spectral pressure density at the first plurality of points; (c) a step of extrapolation and/or interpolation of the spectral pressure density at the first plurality of points in order to obtain a spectral pressure density at second plurality of points of the grid, the second plurality of points having no sensor; (d) a step of estimation of coherence of pressure for multiple pairs of points of the grid, from a predetermined model of a coherence function; (e) a step of calculation of an interspectral pressure density for the multiple pairs of points of the grid, from the coherence of the pressure estimated for these pairs of points, and a spectral pressure density being calculated for the points of these pairs; and (f) a step of calculation of the unsteady aerodynamic loads by summing the real part of the interspectral density in the area of the external structure having a separation of the boundary layer of the airflow.

2. The method of simulation according to claim 1, wherein in step (c) the spectral pressure density of the first plurality of points is extrapolated in linear fashion in the second direction and the spectral pressure density of the first plurality of points is interpolated using an interpolation polynomial of degree greater than or equal to three in the first direction.

3. The method of simulation according to claim 1, wherein in step (d) the model of the coherence function γ.sub.nm(f), for a pair (n,m) of points of the grid, is given by: γ nm ( f ) = exp ( - 2 π f ( .Math. x m - x n .Math. α x ( f ) + .Math. y m - y n .Math. α y ( f ) ) ) .Math. exp ( j φ ( f , x m - x n , y m - y n ) ) φ ( f , x m - x n , y m - y n ) = 2 π f ( x m - x n V x + y m - y n V y ) where x.sub.n, y.sub.n are the respective coordinates of the first point of the pair of points in the first and second directions; x.sub.n, y.sub.n are the respective coordinates of the second point of the pair of points in the first and second directions; V.sub.x is a speed of propagation of instability in the first direction; V.sub.y is a speed of propagation of the instability in the second direction, f is a frequency; α.sub.x (f) is a decoherence factor in the first direction and α.sub.y (f) is a decoherence factor in the second direction.

4. The method of simulation according to claim 3, wherein an estimation {circumflex over (α)}.sub.x(f) of the decoherence factor in the first direction is calculated by the following average: α ^ x ( f ) = - 2 π f Δ x K .Math. k = 1 K 1 ln ( .Math. γ k Δ x ( f ) .Math. ) where γ.sub.k.sup.Δx(f) is coherence at frequency f between pressure signals measured respectively by two sensors separated by a distance Δx in the first direction, and in that an estimation {circumflex over (α)}y(f) of the decoherence factor in the second direction is calculated by the following average: α ^ y ( f ) = - 2 π f Δ y K .Math. k = 1 K 1 ln ( .Math. γ k Δ y ( f ) .Math. ) where γ.sub.k.sup.Δy(f) is coherence, at frequency f, of the pressure signals measured respectively by two sensors separated by a distance Δy in the second direction, and where K is the number of measurements considered for the average.

5. The method of simulation according to claim 3, wherein an estimation {circumflex over (V)}, of the speed of propagation of the instability in the first direction is calculated from the following average: V ^ x = 2 π Card ( Ω x ) .Math. n , m Ω x Δ x nm λ nm where λ.sub.nm is the gradient of a coherence phase, as a function of frequency, of pressure signals measured respectively by two sensors n,m separated by a relative distance Δx.sub.nm=x.sub.m−x.sub.n in the first direction, Ωx is a first set of pairs of sensors, and Card (Ωx) is a cardinal of Ωx; and in that an estimation {circumflex over (V)}.sub.y of the speed of propagation of the instability in the second direction is calculated from an average of the following gradient: V ^ y = 2 π Card ( Ω y ) .Math. n , m Ω y Δ y nm μ nm Where μ.sub.nm is the gradient of a coherence phase, as a function of frequency, of the pressure signals measured respectively by two sensors n,m separated by a relative distance Δy.sub.nm=y.sub.m−y.sub.n in the second direction, and Ω.sub.y is a second set of pairs of sensors.

6. The method of simulation according to claim 5, wherein the speed of propagation of the instability in the first direction and/or in the second direction is obtained for pairs of sensors the measured pressure signals of these pairs of sensors having a coherence module greater than a predetermined threshold (χ.sub.Th).

7. The method of simulation according to claim 1, wherein in step (e) an interspectral pressure density S.sub.nm(f) for a pair of points of the grid is obtained by:
S.sub.nm(f)=γ.sub.nm(f)√{square root over (S.sub.nm(f))}√{square root over (S.sub.nm(f))} where γ.sub.nm(f) is a pressure coherence between the first and second plurality of points of the pair of points; S.sub.nm(f) is the spectral pressure density at the first point of the pair of points and S.sub.nm(f) is the spectral pressure density at the second point of the pair of points.

8. The method of simulation according to claim 1, wherein in step (f) the aerodynamic loads on the external structure are obtained from: F excitation ( f ) = .Math. n , m Ω Re ( S nm ( f ) ) σ .fwdarw. n σ .fwdarw. m where S.sub.nm(f) is an interspectral pressure density between a first point n and a second point m of the grid, f is a frequency, Re(.) designates the real part, where the surface of the external structure is divided into elementary surfaces associated with the different points of the grid, {right arrow over (σ)}.sub.n is the vector with norm equal to the elementary surface associated with a first point n and orthogonal to the external structure at this point, {right arrow over (σ)}.sub.m is the vector with a norm equal to the elementary surface associated with second point m and orthogonal to the external structure at this point, where the summation of the interspectral density is accomplished for a set Ω of pairs of points of the grid belonging to the area of the external structure having a separation of the boundary layer of said flow.

9. A method of simulation by computer of the buffeting of an external structure of an aircraft in an airflow, wherein a simulation is made of unsteady aerodynamic loads, excitation F.sub.excitation(f), being exerted on this external structure according to the simulation method of claim 1, and in that the equation is resolved to the following vibrations:
M{umlaut over (z)}+C′ż+K′z=F.sub.excitation where M is the mass matrix of the external structure, C′ and K′ are respectively the damping and rigidity matrices of the external structure in the wind, and z is the coordinate in the direction orthogonal to the plane of the external structure.

Description

BRIEF DESCRIPTION OF THE ILLUSTRATIONS

(1) Other characteristics and advantages of the invention will appear on reading a preferential embodiment of the invention, made in reference to the attached figures, among which:

(2) FIG. 1 represents diagrammatically a grid of pressure measurement sensors for an aircraft wing surface;

(3) FIG. 2A represents the phase of the coherence function between pressure measurements of two sensors separated in the direction of the wing surface's chord;

(4) FIG. 2B represents the phase of the coherence function between pressure measurements of two sensors separated in the direction of the wing surface's span;

(5) FIG. 3A represents the real part of the interspectral pressure density for two points located close to the leading edge of a horizontal empennage;

(6) FIG. 3B represents the real part of the interspectral pressure density for two points located close to the trailing edge of a horizontal empennage;

(7) FIG. 4 represents in the form of a block diagram a method of simulation by computer of the unsteady aerodynamic loads on an external structure of the aircraft, according to one embodiment of the invention.

DETAILED DESCRIPTION OF PARTICULAR EMBODIMENTS

(8) The method of simulation by computer described below relates to the estimation of the aerodynamic loads on an external structure of an aircraft, also designated more simply as the “structure” in the remainder of the description. The expression “external structure of an aircraft” is understood in this case to mean an element of the aircraft subject in the course of the flight to an airflow at its surface, such as a wing surface or an empennage, for example. This method uses, as above, a semi-empirical model.

(9) More specifically, a model of the structure, equipped with pressure sensors, is firstly subjected to a measuring campaign in a wind tunnel. This campaign enables the unsteady pressures at different points of a measuring grid on the surface of the said structure to be determined. The grid has a first direction Ox in the direction of the flow, and a second direction Oy separate from the first. For example, if the structure is a wing surface the first direction could be the chord and the second direction the wing span, as illustrated in FIG. 1. The lines of the grid are aligned in the first direction and the columns of the grid are aligned in the second direction.

(10) The coordinates of a point at the intersection of a column i and of a line j of this grid will be noted (i,j), δx is the grid interval in the first direction and δy is the grid interval in the second direction. In what follows the points of the grid are, in addition, indexed by n=1, . . . N. The sensors will be indexed by the indices of the points of the grid where they are located.

(11) For each sensor n present at a point (i,j) of the grid the spectral density S.sub.nn(f) of the pressure signal is firstly determined by:

(12) S nn ( f ) = lim T .fwdarw. ( 1 T E ( P ~ n ( f ) .Math. P ~ n * ( f ) ) ) ( 4 )
where E(.) designates the expectation, {tilde over (P)}.sub.n*(f) is the conjugate of {tilde over (P)}.sub.n(f) and {tilde over (P)}.sub.n(f)=TF(P.sub.n(t)) is the Fourier transform of the pressure signal, P.sub.n(t) measured by sensor n, over a time interval T. It is recalled that the spectral density of a signal gives the frequency distribution of the energy of this signal. In practice it is supposed that the hypothesis of ergodicity is checked, and the expectation is calculated in (4) by averaging {tilde over (P)}.sub.n(f).Math.{tilde over (P)}.sub.n*(f), over several time intervals T.

(13) For each pair of sensors n,m, the interspectral density of the pressure signals measured respectively by these sensors is then determined. The interspectral density of the pressure signals of sensors n,m, S.sub.nm(f) is given by:

(14) S nm ( f ) = lim T .fwdarw. ( 1 T E ( P ~ n ( f ) .Math. P ~ m * ( f ) ) ) ( 5 )
where {tilde over (P)}.sub.n(f)=TF(P.sub.n(t)) is the Fourier transform of the pressure signal, P.sub.n(t), measured by the first sensor and {tilde over (P)}.sub.m(f)=TF(P.sub.m(t)) is the Fourier transform of the pressure signal, P.sub.m(t), measured by the second sensor. It is recalled that the interspectral density of two signals gives the distribution of the energy common to these two signals according to the frequency. Due to the combination in expression (5) it also gives, for each frequency, information concerning the phase shift between the two pressure signals. In the same manner as for the spectral density, the expectation in (5) is calculated, in practice, by averaging {tilde over (P)}.sub.n(f).Math.{tilde over (P)}.sub.m*(f), over several time intervals T.

(15) The coherence function between the pressure signals of two sensors is defined by:

(16) γ nm ( f ) = S nm ( f ) S nn ( f ) S mm ( f ) ( 6 )

(17) Coherence module |γ.sub.nm(f)| between the two signals is between 0 and 1. When the coherence module is equal to 1, the spectral distributions of the two signals coincide perfectly. Conversely, when the coherence module is equal to 0, the two signals are entirely decorrelated.

(18) The man skilled in the art will understand that the calculation of the spectral pressure of the measuring signals and of the interspectral power of the pairs of measuring signals is equivalent to calculating the interspectral matrix of these signals. Indeed, the interspectral matrix of these signals is expressed in the following form:

(19) 0 .Math. ( f ) = ( S 11 ( f ) S 12 ( f ) .Math. S 1 N ( f ) S 21 ( f ) S 22 ( f ) .Math. S 2 N ( f ) .Math. .Math. .Math. S N 1 ( f ) S N 2 ( f ) .Math. S NN ( f ) ) ( 7 )
Matrix Σ(f) is hermitian and positive.

(20) In the same way, the coherence matrix of the signals may be defined by:

(21) Γ ( f ) = ( 1 γ 12 ( f ) .Math. γ 1 N ( f ) γ 21 ( f ) 1 .Math. γ 2 N ( f ) .Math. .Math. .Math. γ N 1 ( f ) γ N 2 ( f ) .Math. 1 ) ( 8 )
which, like matrix Σ(f), is hermitian and positive.

(22) When the spectral density of the signals of the different pressure sensors has been calculated, this density is estimated for the missing measurements, in other words for the points of the grid with no sensor.

(23) More specifically, the spectral densities for the points of the grid with no sensor are obtained using an operation to extrapolate/interpolate the spectral densities available at points which have a sensor. The spectral density of the missing measurement at 120 may thus be obtained by extrapolation of density of the measurement at 121, using the ratio of the spectral densities obtained in an adjacent column. In other words, spectral density S.sub.120(f) at point 120 is obtained, for example, by

(24) S 120 ( f ) = S 121 ( f ) S 120 ( f ) S 121 ( f ) ,
where S′.sub.120′(f), S.sub.121(f), S.sub.121′(f) are the spectral densities of the pressure signal at points 120′, 121, 121′.

(25) Once each column of the grid has been completed in this manner any missing columns are then added. Unlike the method of interpolation in relation to the spectral amplitude which is known from the prior art, in this case a polynomial of interpolation of a degree at least equal to 3 is used, for example a cubic spline which gives satisfactory regularity of the spectral density distribution in the first direction of the grid. In other words, to estimate a spectral density at a point 120 (FIG. 1), an interpolation function (of a degree higher than or equal to 3) is applied, passing through the spectral densities of points 131 and 132.

(26) An estimate is then made of the interspectral pressure density for each pair of points of the grid.

(27) This estimate is based on an original modelling of the coherence function. Indeed, it has been possible to show, advantageously, that the coherence function could be modelled in the following form:

(28) γ nm ( f ) = exp ( - 2 π f ( .Math. x m - x n .Math. α x ( f ) + .Math. y m - y n .Math. α y ( f ) ) ) .Math. exp ( j φ ( f , x m - x n , y m - y n ) ) ( 9 )
where α.sub.x(f) and α.sub.y(f) are decoherence factors in the first direction and the second direction respectively. α.sub.x(f) and α.sub.y(f) are also called longitudinal and transverse Corcos coefficients.

(29) The multiplicand of expression (11) reflects the decoherence between the pressure signals as a function of the distance between the sensors (or more generally between the points of the grid when the spectral and interspectral densities are obtained by extrapolation/interpolation). It was possible to show that, in quite a surprising manner, the decoherence factors α.sub.x(f) and α.sub.y(f) depended only slightly on the geometry of the structure and on the external conditions of the flow.

(30) The multiplier is a phase term representing the propagation time of the points of instability (i.e. of the vortexes) between the two points of the grid. The φ(f, x.sub.m−x.sub.n, y.sub.m−y.sub.n) phase depends on the frequency and on the differences of the coordinates of the points of the pair considered:

(31) φ ( f , x m - x n , y m - y n ) = 2 π f ( x m - x n V x + y m - y n V y ) ( 10 )
where V.sub.x and V.sub.y are respectively the speeds of propagation of the instability (of the vortexes) in the first direction and the second direction.

(32) According to a first variant the decoherence factors are calculated once only after a first measuring campaign.

(33) According to a second variant the decoherence factors are calculated in each campaign, i.e. for each new structure, as follows:

(34) α ^ x ( f ) = - 2 π f Δ x K .Math. k = 1 K 1 ln ( .Math. γ k Δ x ( f ) .Math. ) ( 11 - 1 ) α ^ y ( f ) = - 2 π f Δ y K .Math. k = 1 K 1 ln ( .Math. γ k Δ y ( f ) .Math. ) ( 11 - 2 )
where γ.sub.k.sup.Δx(f) is the coherence function calculated experimentally from expression (6) for two sensors separated by a distance Δx (therefore Δx>0) in the first direction of the grid, γ.sub.k.sup.Δy(f) is the coherence function calculated experimentally from expression (6) for two sensors separated by a distance Δy (therefore Δy>0) in the second direction, k is the index of the measurement and K is the number of measurements considered.

(35) It is important to note that the calculations of the decoherence factors according to (11-1) and (11-2) can be made validly only if the measurements of both sensors n,m have a high coherence module, in other words |γ.sub.nm(f)|≧γ.sub.th, where γ.sub.th is a threshold value dependent on the measuring time and the processing accomplished. It was able to be shown that a value of γ.sub.th=0.4 was suitable in practice.

(36) The propagation speeds of the instability, V.sub.x and V.sub.y, are also determined experimentally from the pressure measurements.

(37) More specifically, for a pair n,m of sensors separated by a relative distance Δx.sub.nm (positive or negative) in the first direction, phase φ.sub.nm.sup.Δx(f) of the coherence function is determined in relation with the frequency, for different measurements (where each measurement gives a pressure signal which changes over time). The speed of propagation of the instability in the first direction is estimated from an average of gradient α.sub.nm of phase φ.sub.nm.sup.Δx(f) as a function of the frequency:

(38) V ^ x = 2 π Card ( Ω x ) .Math. n , m Ω x Δ x nm λ nm ( 12 - 1 )
where Ωx is the set of the pairs of sensors considered for the average in question, and Card(Ωx) is the cardinal of Ωx, i.e. the number of pairs of sensors in considered.

(39) Similarly, for a pair of sensors n,m separated by a relative distance Δy.sub.nm (positive or negative) in the second direction, phase φ.sub.nm.sup.Δy(f) of the coherence function is determined, together with its gradient μ.sub.nm, as a function of the frequency, and the speed of propagation of the instability in the second direction is deduced by means of:

(40) V ^ y = 2 π Card ( Ω y ) .Math. n , m Ω y Δ y nm μ nm ( 12 - 2 )
where Ωy is the set of the pairs of sensors considered for the average in question, and Card(Ωy) is the cardinal of Ωy, i.e. the number of pairs of sensors considered.

(41) FIG. 2A represents phase φ.sub.nm.sup.Δx(f) of the pressure coherence function between two points of the grid separated by a distance Δx.sub.nm in the first direction. FIG. 2B represents phase φ.sub.nm.sup.Δy(f) of the pressure coherence function between two points of the grid separated by Δy.sub.nm in the second direction. In both cases frequency f has been represented in the form of the Strouhal number, V. It is recalled that the Strouhal number is a normalised frequency

(42) V = f f
with

(43) f = V L ,
where V.sub.∞ is the speed of the undisturbed flow and L is a length characteristic of the structure.

(44) It is remarked in the example of FIG. 2A that curve φ.sub.nm.sup.Δx(f) has a linear portion up to a relatively high Strouhal number. Conversely, in the example of FIG. 2B, curve φ.sub.nm.sup.Δy(f) has a smaller linear portion. The loss of linearity of the phase is due essentially to the loss of coherence as a function of frequency, since the loss of coherence is more pronounced in the direction of the wing span than in the direction of the chord. It is also noted that the gradient of curve φ.sub.n.sup.Δx(f) is higher than that of φ.sub.n.sup.Δy(f).

(45) Coherence function γ.sub.nm(f) of the pressure for two pairs of points of the grid is determined from equations (9) and (10), after the decoherence factors have been estimated by means of (11-1) and (11-2) and the speeds of propagation of the points of instability have been estimated by means of (12-1) and (12-2).

(46) It is then possible to deduce from this the interspectral pressure density between two points n,m of the grid according to (6), or in other words:
S.sub.nm(f)=γ.sub.nm(f)√{square root over (S.sub.nn(f))}√{square root over (S.sub.mm(f))}  (13)

(47) Lastly, the unsteady aerodynamic loads on the wing surface may be estimated by:

(48) 0 F excitation ( f ) = .Math. n , m Ω Re ( S nm ( f ) ) σ .fwdarw. n σ .fwdarw. m ( 14 )
where the summation is made for the set Ω of pairs of points of the grid belonging to the separation area of the boundary layer, Re(.) means the real part and {right arrow over (σ)}.sub.n is a vector normal to the surface of the structure at point n having a norm equal to an elementary surface σ.sub.n around this point. In other words, the separation surface is divided into elementary surfaces associated respectively with the different points of the grid, and the summation of the real part of the interspectral pressure density is made for these elementary surfaces.

(49) It is understood according to expression (14) that the accuracy of the estimate of the unsteady aerodynamic loads is directly related to that of the real part of the interspectral pressure density.

(50) In FIGS. 3A and 3B the real part of the interspectral pressure density has been represented for a pair of points located on the horizontal tail plane or HTP of a T-tail tail plane configuration. In each of these figures the interspectral pressure density estimated by the simulation method of the prior art has been designated by 320, the interspectral pressure density estimated by the method of simulation according to one embodiment of the present invention by 310, and the real interspectral pressure density obtained from the pressure measurements of the sensors located at these points by 330. The interspectral densities are given according to the Strouhal number.

(51) In FIG. 3A the points considered are located close to the leading edge of the tail plane, whereas in FIG. 3B the points considered are located close to the trailing edge.

(52) It is observed in FIG. 3A that the simulation by the method of simulation according to the invention, 310, gives a very satisfactory approximation of the real spectral density, 330, and in any event one which is much better than that obtained with the simulation method of the prior art, 320.

(53) It is observed in FIG. 3B that the interspectral density obtained by the simulation method according to the invention, 310, leads to an overestimate of the real interspectral density, 330. However, the general appearance of the curve is followed, and the approximation is, once again, much better than the one obtained with the simulation method of the prior art, 320.

(54) FIG. 4 represents, in the form of a block diagram, a method of simulation by computer of the unsteady aerodynamic loads on an external structure of the aircraft, according to one embodiment of the invention.

(55) This simulation method begins with a step 410, during which a campaign of measurement of the unsteady pressure at the surface of the structure in the course of wind tunnel testing is conducted. The pressure measurements are made using multiple sensors located on the structure and positioned at the first points of a grid. This grid has a first direction in the direction of the flow, and a second direction distinct from the first. The grid interval in the first direction is δx and the grid interval in the second direction δy. Each pressure signal obtained by a sensor is stored in a memory after being sampled.

(56) In step 420, for each of the sensors, or in other words for each of points n of the grid equipped with a sensor, the calculation is made of the spectral density of the pressure signal at this point. The calculation of the spectral density at a point is made from the expression (4), i.e., in practice, by taking segments of duration T of the pressure signal at this point, by applying the Fourier transform of this signal, and by calculating the average of {tilde over (P)}.sub.n(f){tilde over (P)}.sub.n*(f) for these segments.

(57) In step 430 operations are made to extrapolate/interpolate the spectral density of the pressure in order to determine the spectral pressure density at the points of the grid which are not equipped with sensors.

(58) More specifically, extrapolation is first applied in the second direction. When a sensor is missing at a point of a column of the grid, the spectral density at this point is obtained by extrapolating the spectral density at an adjacent point in the same column. The extrapolation coefficient is determined from the ratio between spectral densities of the corresponding points in the column closest to the one considered.

(59) An interpolation is then made of the spectral density in the first direction, using an interpolation polynomial of degree higher than or equal to 3, for example a cubic spline. The coefficients of the interpolation polynomial may differ according to the line considered.

(60) In step 440 the coherence of the pressure is estimated for each pair of points of the grid or, at the least, for each pair of points belonging to the area of the structure having a separation of the boundary layer. This estimate is based on the original model of the coherence function given by expressions (9) and (10).

(61) To accomplish this, decoherence factors are used in the first and second directions respectively, as calculated in a previous campaign. Alternatively, the decoherence factors are estimated using expressions (11-1) and (11-2), by calculating the coherence of the pressure signals measured by two sensors. The respective propagation speeds of the instability in the first and second directions are also estimated, using expressions (12-1) and (12-2), by calculating the phase gradient of the coherence of the pressure signals measured for pairs of sensors.

(62) In step 450 the interspectral pressure density is estimated for each pair of points of the grid. For a pair of points n,m of the grid the interspectral pressure density is determined, by means of expression (13), from the pressure coherence between these two points, estimated in step 440, and from the spectral pressure density at each of these two points, as calculated in step 420 (if the point is equipped with a pressure sensor), or obtained in step 430 by interpolation (if the point is not so equipped).

(63) Lastly, in step 460, the aerodynamic loads being exerted on the structure are calculated by means of expression (14), or in other words by summing the real part of the interspectral density for the area of the structure having a separation of the boundary layer of the flow.

(64) The determination of the aerodynamic loads F.sub.excitation(f) may form part of a method of simulation of the buffeting of the structure in an airflow. In this case, after having estimated the aerodynamic loads as described above, the reduced elasticity equation (2) is resolved. This equation is advantageously resolved in the frequencies space (by applying a Laplace transform), such that it is not necessary to apply an inverse Fourier transform of F.sub.excitation(f). To accomplish this, the DLM method or the generalised Euler method may be used, in a known manner.