3. Calculation of impedance in time

As mentioned earlier, the Laplace-Temps method is a numerical approach that makes it possible to discretize a convolutional integral in time.

3.1. Implementation in code_aster

Equation () shows that the numerical evaluation of weights:math: {Psi} _ {Psi} _ {k} ^ {j} first goes through that of:math: {widehat {Z}}} _ {m}} _ {m} ^ {k}. In turn, this can be obtained from a Cauchy integral applied to a contour:math: `|z|=mathrm {rho} `:

:math:`{Z}_{m}(kmathrm{Delta }t)=frac{1}{2mathrm{pi }i}{int }_{

z

=mathrm{rho }}{widehat{Z}}_{m}left(frac{mathrm{delta }(z)}{mathrm{Delta }t}right){z}^{-k-1}dz`

If we express this integral in polar coordinates \(z=\mathrm{\rho }{e}^{i\mathrm{\theta }}\), then we apply the trapezoidal rule to the phase to discretize it into \(L\) not identical with a value of \(\mathrm{\Delta }\mathrm{\theta }=\frac{2\mathrm{\pi }}{L}\), we obtain:

(3.1)\[ {Z} _ {m} (k∆t)\ approx {Z} _ {Z} _ {m} _ {m} ^ {k} =\ frac {{\ rho}} {L}\ sum _ {l=0} ^ {0} ^ {L-1} ^ {L-1} ^ {L-1} ^ {L-1} ^ {L-1} ^ {L-1} ^ {L-1} ^ {L-1} ^ {L-1} {L-1} {L-1} {\ L-1} {\ L-1} {\ widehat {Z}}} _ {m} ({s} _ {l}) {e} ^ {-i\ frac {2\ frac {2\ pi l} {L} k}\]

where the operator defined in the Laplace domain is sampled on \({s}_{l}=\frac{\mathrm{\delta }(\mathrm{\rho }{e}^{2\mathrm{\pi }il∕L})}{\mathrm{\Delta }t}\) with \(\mathrm{\delta }(z)=\frac{3}{2}-2z+\frac{1}{2}{z}^{2}\) the characteristic polynomial of a second-order multi-step linear method. Note that \(N\) and \(L\) are not necessarily equal, \(N\) being the number of time steps in the time window of interest (for example, as a function of the duration of the seismic signal) and \(L\) being the number of time steps for calculating the impedance in time.

If we consider that \({\widehat{Z}}_{m}({s}_{l})\) is evaluated with variability [2] _ \({ϵ}_{\mathit{CQM}}\), the \({\widehat{Z}}_{m}^{k}\) coefficients will be calculated with precision [3] _

\(O(\sqrt{{ϵ}_{\mathit{CQM}}})\) if \(L=N\) and \({\rho }^{N}=\sqrt{{ϵ}_{\mathit{CQM}}}\). With the values of \(L\) and \(\rho\) fixed, the weights of the series can be obtained from a classical FFT, which reduces the complexity of the calculation algorithm to \(O(L\mathrm{log}L)\) instead of \(O({L}^{2})\). However, the value of \({ϵ}_{\mathit{CQM}}\) is not an easy step to calibrate and depends in part on the nature of the problem to be solved. In fact, values that are too low of \({ϵ}_{\mathit{CQM}}\), such as \({10}^{-20}\) or \({10}^{-30}\), may result in integrative radii \(\mid z\mid =\rho\) that are too small. Since the Cauchy integrand contains a zero singularity, too small values of the radius may result in a poor evaluation of the integral. As such, parametric studies were carried out in to conclude that a good compromise is reached when \({ϵ}_{\mathit{CQM}}={10}^{-10}\) [4] _

.

Likewise, it will be noted that the dependence in \({\rho }^{-k}\) leads, for large \(k\) values, to inaccuracies in the calculated \({\widehat{Z}}_{m}^{k}\) coefficients. Knowing that the use of more precise integration methods, such as Simpson for example, does not correct this problem, the idea is to extend the time range (i.e. \(L=\mathit{mN},m\mathrm{\in }\mathrm{ℝ}\)) over which the impedance is calculated in order to reduce disturbances on the solution of interest. As such, it has been observed that, in general, the inaccuracies of \({\widehat{Z}}_{m}^{k}\) only spread to the response beyond \(t\approx \mathrm{0,7}{T}_{\mathit{FIN}}\), where \({T}_{\mathit{FIN}}\) is the final moment for calculating the impedance [5] _ .

3.2. Key numerical considerations

We saw earlier that, in general, the singular part of the ground impedance can be approximated by a polynomial of order two. If the coefficients of this polynomial are real, it is easy to demonstrate that the ground impedance meets \(\widehat{Z}(\sigma -i\omega )=\mathit{conj}[\widehat{Z}(\sigma +i\omega )]\) in the analytical plane of the Laplace domain [6] _ . It will be noted that when hysteretic damping is used in the ground (case of MISS3D), the imaginary part of its static stiffness is non-zero and the hypothesis of real coefficients is no longer verified.

A second interesting property coming directly from Laplace’s definition of unilateral transform is the following:

(3.2)\[ \ widehat {Z} (s=0) = {\ int} _ {0} _ {0} ^ {0} ^ {0} ^ {N} ^ {N} Z (k\ Delta t)\ approx\ sum _ {k} ^ {N} Z (k\ Delta t)\]

for \(N\) chosen big enough (~ a hundred steps of time). Note that the ground impedance evaluated at \(s=0\) corresponds to the static ground stiffness (real part), which is very different from the instantaneous ground stiffness. This must be taken into account when a dynamic calculation (\(Z(t=0)]\) is assembled at the first limb) is performed following a static calculation (\(\mathrm{\Re }[\widehat{Z}(s=0)]\) is assembled at the first limb) [7] _ .

[1]

Referring to the concept of reproducibility (precision in English).

[2]

Referring to the concept of accuracy (accuracy in English).

[3]

By default, this value initializes the PRECISION keyword for the CALC_MISS operator (cf. [U7.03.12]) when the expected result type is” FICHIER_TEMPS “.

[4]

To do this, CALC_MISS (TYPE_RESU = “FICHIER_TEMPS”) provides for the keyword COEF_SURECH (cf. [U7.03.12]) initialized to 1.35, which is the equivalent of doing \(L=\mathrm{1.35N}\) in equation 23.

[5]

This property is used in code_aster: the impedance is sampled on the upper half-plane, then completed numerically with the conjugate. Therefore, care must be taken with soils with very strong hysteretic damping (greater than 25%).

[6]

Code_aster has a macro-command (PRE_SEISME_NONL, cf. [U4.63.02])) that makes it possible to perform the static-dynamic transition of a non-linear soil-structure interaction calculation.