# Ellipsometry of surface acoustic waves using 3D vibrometry for viscoelastic material characterization by the estimation of complex Lamé coefficients versus the frequency

Aziz Bouzzit, Loïc Martinez, Andrés Arciniegas, Stephane Serfaty, Nicolas Wilkie-Chancellier

► To cite this version:

Aziz Bouzzit, Loïc Martinez, Andrés Arciniegas, Stephane Serfaty, Nicolas Wilkie-Chancellier. Ellipsometry of surface acoustic waves using 3D vibrometry for viscoelastic material characterization by the estimation of complex Lamé coefficients versus the frequency. Applied Acoustics, 2025, 228, pp.110312. 10.1016/j.apacoust.2024.110312. hal-04724297

HAL Id: hal-04724297

https://hal.science/hal-04724297v1

Submitted on 8 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Ellipsometry of surface acoustic waves using 3D vibrometry for viscoelastic material characterization by the estimation of complex Lamé coefficients versus the frequency

Aziz Bouzzit*, Loïc Martinez, Andres Arciniegas, Stéphane Serfaty, Nicolas Wilkie-Chancellier

Laboratoire SATIE (UMR CNRS 8029), CY Cergy-Paris Université, 5 mail Gay Lussac, 95031 Neuville Sur Oise, France

## ARTICLE INFO

**Keywords:**

Viscoelastic materials

Ellipsometry

3D vibrometer

Complexes Lamé coefficients

Rayleigh wave

## ABSTRACT

Surface acoustic waves (SAW) are adequate regarding material characterization because they have low geometric attenuation compared to bulk waves. SAW can be generated easily by normal excitation using contact transducers or power lasers and have also a unique elliptic polarization, characterized by two parameters: the ellipticity (H/V) ratio between the horizontal and the vertical components of the elliptic motions and the orientation angle (θ) between the horizontal axis of the ellipse and the surface. In the case of a viscoelastic isotropic material, a complete characterization is achieved by the association of the quantitative measurement of the polarization and the propagative characteristics, the complex wavenumber, of the SAW. In practice, this operation is performed using 3D lased vibrometry for propagation monitoring in space and time. The post-processing is carried out by Quaternion Fourier Transform, the Prony algorithm and the complex Lamé coefficients identification for the theoretical model of propagation on the material. Good agreement is observed between the obtained results and the ones of the pulse-echo method.

## 1. Introduction

The study and characterization of viscoelastic materials hold great significance since they are encountered when one investigates the propagation of acoustic waves through various media and in different application domains. Materials used for structural applications may manifest viscoelastic behaviors, which have a significant impact on their performance. In engineering applications, the viscoelastic behavior of a material may exhibit as an unintentional side effect, or for other applications the viscoelasticity of a material may be deliberately exploited in the design process to achieve comportmental goals. Within branches of the community of applied mathematics, the mathematical foundations of viscoelasticity theory spark interest. More, in materials science, metallurgy, and solid-state physics, there is a strong interest in visco-elasticity because it is causally related to a variety of microphysical processes and can be used as an experimental probe of these processes. The relationship between viscoelasticity and microstructure is utilized in viscoelastic testing as an inspection method and material design. Many applications of the viscoelastic materials characterization can be found when dealing with natural materials such as stone, earth, and wood in the case of building construction and the monitoring of their

integrity. One can see the diversity in application domains and the usefulness of the characterization of the viscoelastic aspect of materials, such as bone or soft human tissue monitoring in the case of biomedical engineering and medical diagnosis, the study of the behavior and the integrity of concrete structures, in the fabrication of tires and the filament of light bulbs, or in the study of the relaxation of the musical instrument strings [1].

The viscoelastic materials behavior is demonstrated through its response to dynamic excitation. Various techniques are employed to estimate these properties, based on measuring the different material deformations resulting from a controlled excitation. Methods like tensile testing [2], dynamic mechanical analysis [3] and parallel-plate rheometry [4] are of this type, these lasts are limited to low frequencies. Methods based on wave propagation in and on the material can be used to overcome this limitation. Such as shear wave dispersion ultrasonic vibrometer whose method uses the phase difference of the shear wave between two locations along the propagation path to qualitatively estimate the dispersion of the velocity. The viscoelastic properties of the medium can then be estimated using an inverse model [5]. This method can be extended to plate-like structures: in this case, shear waves can be replaced by Lamb or Rayleigh waves [6]. This approach uses only the

<footnote>* Corresponding author.
E-mail address: aziz.bouzzit@cyu.fr (A. Bouzzit).</footnote>
measurement of the wave velocity and does not make use of the decrease in the amplitude of the wave (attenuation). Monitoring the propagation of the normal component of the surface Rayleigh wave by a 1D laser vibrometer, the complex wavenumber can be estimated [7], and via the invert problem, the complex shear/elastic moduli can be estimated also.

Compared to the shear waves, which are volume waves, Rayleigh waves are more accessible in terms of generation and monitoring as both aspects of the wave can be done on the surface, and remotely using a power laser for the excitation [8] and Laser vibrometers for the vibration detection and monitoring [9]. Also, the theoretical modelization of the propagation of Rayleigh waves on the surface of an isotropic linear and homogenous viscoelastic medium is well studied over the last century, starting with the work of Scholte in 1947 [10] up to the latest advancements presented by Sharma [11]. In this work, Sharma presented a complex analysis technique to simplify and then solve the secular equation that governs the propagation of the Rayleigh wave, also the expressions for the particle motion are presented in the sagittal plane of the propagation. These expressions are of great interest as they define the polarization of the surface wave.

In this study the characterization by monitoring both components is completed using a 3D laser vibrometer, extracting the complete polarization in the sagittal plane. This operation have been done using a 1D laser with an angle of incidence, as demonstrated by [9], but the experiment needs to be done twice for each extraction position with different angles of incidence, this may affect the repeatability of experiment as the small variation can induce large deviation on the results. Also, the positions of the extraction need to be on the sagittal plane of the wave, out of that plane the three components of the deformation are mandatory in order to reproject the signals on the sagittal plane. Therefore these limitations are surmounted by the use of a 3D laser vibrometer. In addition to the propagative properties of the wave (complex wave number), and by utilizing the theoretical modeling of the propagation and the polarization, an inverse problem allows the complete characterization of the viscoelastic medium by estimating the complex Lamé coefficients as a function of frequency.

This paper is organized in three parts. Firstly, the theoretical modeling of the direct problem allows the estimation of the complex wave number and the polarization from the mechanical properties, i.e. by the complex Lamé coefficients. Secondly, the inverse problem and the algorithms allowing the measurement of the different propagative parameters of the surface wave are presented. Finally, the experimental validation is carried out with the study of a Rayleigh wave propagating on a block of Epoxy.

## 2. Theoretical background

It is well known that the elastic deformation of an isotropic material is the superposition of two elementary components, the volume-conservative shear deformation that can be estimated by the shear modulus G and volume-non conservative elastic deformation that can be characterized by bulk modulus K. Lamé coefficients (λ, μ) are a combination of G and K and can also be used to model elastic materials, as they can simplify the notations. Indeed, when working with Hooke's law in 3D, this last formalism is used for the study of the propagation of Rayleigh waves. When subjected to dynamic loading, a viscoelastic material experience a delayed deformation from its initial position, attributed to internal viscous friction within the material. When exposed to harmonic forcing, this delay is evident as a phase shift between the applied load and the resulting deformation. The magnitude of this shift is directly related to the viscous losses within the material. As a result of this phase lag, Lamé coefficients can be treated as frequency-dependent complex functions [12].

$$ \begin{cases} \lambda^*(\omega) = \lambda'(\omega) + i\lambda''(\omega) \\ \mu^*(\omega) = \mu'(\omega) + i\mu''(\omega) \end{cases} \qquad (1) $$

The real parts {λ'(ω), μ'(ω)} are related to the elastic behavior of the material and the imaginary part {λ''(ω), μ''(ω)} to the viscous one. The theoretical expression of the variation of the complex Lamé parameters as function of the frequency are given by different models, for instance by Voigt model, standard Linear model also known as Kelvin or Zener model, Maxwell model and fractional order model [13,14].

These models express the real and imaginary parts of Lamé coefficients as a function of the frequency of the harmonic excitation and material-related properties. Being able to measure the variation of the complex Lamé parameters as a function of frequency allows us to retrieve the mechanical properties by choosing the adequate model to be fitted to the measurements.

Let's consider a viscoelastic material characterized by its complex Lamé coefficients and its density (μ*, λ*, ρ), occupying the half-space region $x_3 \le 0$ (as illustrated in Fig. 1a). The propagation of a Rayleigh wave with an angular frequency (ω) on the free surface $x_3 = 0$, following the $x_1$ direction, is governed by the secular equation [15]:

$$ R^{*3} - 8R^{*2} + \left( 24 - 16 \frac{\mu^*}{\lambda^* + 2\mu^*} \right) R^* - 16 \left( 1 - \frac{\mu^*}{\lambda^* + 2\mu^*} \right) = 0 \qquad (2) $$

where:

$$ R^* = \frac{\rho}{\mu^*} \left( \frac{\omega}{k_R^*} \right)^2 \qquad (3) $$

The roots of this equation determine the complex value of R, and hence the wavenumber $k_R^*$ of the Rayleigh wave propagating on the surface of the half-space. $k_R^*$ contains the propagative information of the Rayleigh wave, the phase velocity ($\omega/k_R^*$) and the attenuation ($k_R^*$). This information can be estimated by monitoring the propagation of the wave on the propagation surface.

Unfortunately, the inverse problem cannot be only constructed based on these two values, as the complete centralization of the medium is achieved by estimating the four coefficients (λ*, μ*, λ'' and μ'') of Eq. (1). One is in a situation of an undefined inverse problem that can accept an infinity of solutions. Other measurable parameters need to be added to overcome this situation. These parameters are related to the polarization of the surface wave. As it is well known that surface waves have an elliptical particle motion trajectory or polarization [16,17], this polarization can be characterized by two parameters: (i) the ellipticity ($\chi = \arctan(H/V)$) calculated from the ratio between the in-plane or horizontal (H) and out-of-plane or vertical (V) components of the particle motion and (ii) the orientation of the ellipse ($\theta$) which is the angle between the major axis of the ellipse and the horizontal axis $x_1$, as illustrated in Fig. 1b).

The theoretical expressions of the polarization parameters ($\chi, \theta$) as function of the complex Lamé parameters can be extracted from the components ($u_1, u_3$) of the particle displacement $\vec{u}$ in the ($x_1, x_3$) plane. For a harmonic excitation with an angular frequency (ω) and on the free surface ($x_3 = 0$), these components are given by the following expressions [15]:

<figure><img src="image_3.png" /><figcaption>Fig. 1. a) Representation of the wave propagation and b) Rayleigh wave polarization.</figcaption></figure>
$$u_n^* = A_n^* e^{i(k_R^* x_1 - \omega t)}, n \in \{1, 3\} \qquad (4)$$

where:

$$\left\{ \begin{array}{l} A_1^* = \frac{R^*}{R^* - 2} \\ A_3^* = \frac{R^*}{2(1 - R^*)^{1/2}} \end{array} \right. \qquad (5)$$

$A_1^*$ and $A_3^*$ are the complex normalized amplitudes of the two components. Eq. (5) gives the expressions of these amplitudes as function of $R^*$, the solution of Eq. (2), and complex Lamé coefficients. The $(u_1, u_3)$ components of the particle deformation correspond to the real part of the complex $(u_1^*, u_3^*)$ amplitudes.

The theoretical expressions of the polarization parameters $(\chi_{theo}, \theta_{theo})$ as a function of $A_1^*$ and $A_3^*$ are given by the following formulas [15]:

$$\left\{ \begin{array}{l} \chi_{theo} = \left| \frac{u_1^*}{u_3^*} \right| = \left| \frac{A_1^*}{A_3^*} \right| \\ \theta_{theo} = \frac{1}{2} \arctan \left( \frac{A_1^* A_3^* + A_3^* A_1^*}{A_1^* A_1^* - A_3^* A_3^*} \right) \end{array} \right. \qquad (6)$$

For a specific angular frequency $(\omega)$, these equations allow us to estimate the elliptic parameters of the polarization using the complex material properties.

$\chi$ and $\theta$ are measurable parameters that can be extracted from the monitoring of the propagation of the two components on the sagittal plane of the surface Rayleigh wave. Adding to these two parameters $(\chi, \theta)$, the extraction of the propagative properties of the wave $(k_R^* = k_R^* + i k_R^*)$, we will have four parameters that will allow us to construct the inverse problem to estimate the four complex Lamé coefficients $(\lambda', \lambda'', \mu', \mu'')$.

In the present case, the inverse problem is constructed using the Eqs. (2), (3) and (6). For a specific value of the angular frequency $(\omega)$ and the density $(\rho)$, one finds the Lamé coefficients $\{\lambda^*, \mu^*\}$ that corresponds to the measured values of $\{\chi, \theta, k_R^*\}$. The optimization process is implemented using a derivative-free method [18] to find the minimum of the multivariable error function Err. This function is constructed from the absolute value of the relative difference between the theoretically estimated and the experimentally measured values of the polarization and the propagation properties of the surface Rayleigh wave:

$$Err = \left\{ \left| \frac{\chi_{exp} - \chi_{theo}}{\chi_{theo}} ; \left| \frac{\theta_{exp} - \theta_{theo}}{\theta_{theo}} ; \left| \frac{k_R^* \exp - k_R^* theo}{k_R^* theo} ; \left| \frac{k_R^* \exp - k_R^* theo}{k_R^* theo} \right| \right| \right| \right\} \qquad (7)$$

The next section proposes a signal processing method that allows us to estimate the polarization $(\chi_{exp}, \theta_{exp})$ and the propagative $(k_R^* \exp)$ properties of a propagating wave as a function of frequency and wavelength, starting from a bivariant space-time signal extracted from the monitoring of the propagation of the wave using a 3D lase vibrometer.

### 3. Polarization ($\chi_{exp}, \theta_{exp}$) estimation of bivariant space-time signals using 2D quaternion Fourier transform (2DQFT)

The monitoring of the propagation of the surface Rayleigh wave is performed by the measurement of the two components $(u_1, u_3)$ of the particle deformation, in different time instants $(t)$ and space positions $(x_1)$ in the direction of propagation along the sagittal plane of the wave. The retrieved signal is a bicomponents or bivariant space-time 2D matrix, and can be represented either as a vector-valued signal $S(x_1, t) = [u_1, u_3]$, or as a complex-valued signal:

$$S(x_1, t) = u_1(x_1, t) + i \cdot u_3(x_1, t) \qquad (8)$$

In this complex representation (8), the bivariant signal can be

considered as a special type of quaternion-valued signal. This allows us to process the signal on its polar form using the 2D quaternion Fourier transform (QFT), giving a direct access to the polarization parameters of the bivariant signal as a function of the frequency and wavenumber [19,20].

Let's present here the quaternion space and the 2D QFT method which are first introduced by Sir William Rowan Hamilton in 1843 [21]. They are the generalization of complex numbers for which two components are defined: the real and the imaginary part. The quaternion is then defined by four components, one real part and three imaginary parts, enabling any quaternion to be represented in a hypercomplex Cartesian form as:

$$q = a + b \cdot i + c \cdot j + d \cdot k \qquad (9)$$

where $a, b, c, d \in \mathbb{R}$ are its components. Imaginary units $i, j, k$ complex operators (generalization of complex operator $i$, also denoted $j$) satisfy the fundamental formula for quaternions multiplication:

$$i^2 = j^2 = k^2 = ijk = -1 \qquad (10)$$

and

$$\begin{aligned} ij &= -ji = k \\ jk &= -kj = i \\ ki &= -ik = j \end{aligned} \qquad (11)$$

For this property, it should be noted that the multiplication of quaternions is not commutative. Care is therefore necessary in using quaternions and in the coding of quaternion-based algorithms to ensure that the ordering of operands is correct. However, the usual operations such as addition, scalar multiplication, and equality behave similarly to complex cases [19]. Similar to complex numbers, any quaternion $q$ can be written in the Euler polar form as:

$$q = |q|e^{i\theta}e^{k\varphi}e^{j\varphi} \qquad (12)$$

where $|q|$ is the amplitude and $\varphi \in [-\pi, \pi]$ is the phase of the quaternion, these parameters are equivalent to those of complex numbers. $\theta \in [-\pi/2, \pi/2]$, $\chi \in [-\pi/4, \pi/4]$ are the already defined polarization parameters. The transition from the cartesian to the polar form and vice versa is possible and the theoretical expressions can be found in literature [19].

Considering two harmonic signals $u$ and $v$ with the same angular frequency $(\omega)$ and different amplitude and phase, they can define a bivariant signal with specific polarization parameters. This signal can be embedded by a quaternion in the same form as in Eq. (12), which gives direct access to the polarization parameters that we are looking for in this study.

This operation can be carried out by the quaternion Fourier transform (QFT), that takes the signal in its complex form as $S = u + v \cdot i + 0 \cdot j + 0 \cdot k$, and transform it into the Euler polar form for different frequencies. The QFT was first studied by Jamson in 1970 [22]. In this case, the signal is 2D space-time matrixes, which means that one has to perform the 2D version of the QFT. This was introduced for the first time by Ell for the study and the analysis of linear time-invariant partial differential systems [23]. For an arbitrary bivariant space-time signal, the 2D quaternion Fourier transform (2DQFT) is defined by:

$$\begin{aligned} \hat{S}(k_1, \omega) &\triangleq \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} S(x_1 t) \cdot e^{-\Lambda(k_1 x_1 + \omega t)} dx_1 dt \\ &= a(k_1, \omega) e^{i\theta(k_1, \omega)} e^{-k\varphi(k_1, \omega)} e^{j\varphi(k_1, \omega)} \end{aligned} \qquad (13)$$

This definition is similar to the usual Fourier transform, at the difference of two fundamental elements. The first is related to the position of the Fourier atom $e^{-\Lambda(k_1 x_1 + \omega t)}$ with respect to the quaternion values signal $S(x_1, t)$ which is crucial due to the non-commutative nature of the product in the quaternion space. It is chosen to place the Fourier atom on the right side of the signal for convenience and to agree with the usual
Fourier transform. The second element of difference is the so-called the axis of the QFT, Λ, this last is a free parameter. It is only restricted to be a pure unit quaternion. The choice of this axis has been argued in literature [19], it is chosen to be Λ = j in our case as it was recommended in the reference.

The implementation of the 2DQFT can be performed using the classical 2D Fourier transform (2DFT), and can be executed efficiently by using the 2DFFT algorithm. This is possible because the quaternion-valued signal $S(x_1, t)$ can be decomposed into a pair of complex valued signals as follow:

$$S(x_1, t) = S_1(x_1, t) + \Lambda_{\perp} S_3(x_1, t) \quad (14)$$

where $\Lambda_{\perp}$ is a pure unite quaternion orthogonal to the axis of the 2DQFT ($\mu$).

By linearity of the 2DQFT, one gets:

$$\hat{S}(k_1, \omega) = \hat{S}_1(k_1, \omega) + \Lambda_{\perp} \hat{S}_3(k_1, \omega) \quad (15)$$

where $\hat{S}_1$ and $\hat{S}_3$ are the 2DFT complex-valued of $S_1$ and $S_3$.

This confirms the possibility of obtaining the 2DQFT by combining two standards 2DFT [19,24].

The use of the 2DQFT on the measured signal $S_{exp}(x_1, t)$ allows to estimate the polarization properties of the Rayleigh wave ($\chi_{exp}, \theta_{exp}$) as a function of the frequency ($\omega$) and wavenumber ($k_1$), as previously shown in the theoretical expressions (6). ($\chi_{exp}(\omega), \theta_{exp}(\omega)$) are the first two elements of the invers problem: the following part proposes an algorithm for the estimation of the remaining two elements.

## 4. Estimation of the propagative properties ($k_R^* \exp$) of Rayleigh wave using the Prony algorithm

Estimation of the propagative properties of Rayleigh wave consists in the measurement of the two parts of the complex wavenumber for different values of frequency:

$$k_R^* \exp(\omega) = k_R^* \exp(\omega) + i k_R'' \exp(\omega) \quad (16)$$

These parameters need to be extracted from the space-frequency signal ($S_{exp}(x_1, \omega)$), which is a simple 1D time Fourier transform of the measured space-time signal ($S_{exp}(x_1, t)$). From Eq. (4), it can be seen that the two components of the wave undergo the same attenuation and propagate with the same wavenumber during propagation, which means that the choice of one component over the other is indifferent in terms of the expected results.

The principle of the Prony method is to identify the unknown parameters ($k_R^*$) by minimizing the quadratic difference between the model $S_{mod}(x_1, \omega)$ and the experimental data $S_{exp}(x_1, \omega)$. Considering that only one surface wave is propagating, the theoretical expression of $S_{mod}(x_1, \omega)$ is as follow [25–27]:

$$S_{mod}(x_1, \omega) = \delta^{-1}(x_1) \cdot H(0, \omega) \cdot e^{-k_R^*(\omega)x_1} \cdot e^{i k_R^*(\omega)x_1} \quad (17)$$

where $H(0, \omega)$ is the spectral amplitude at the position $x_1 = 0$ the origin of observation, and $\delta^{-1}(x_1)$ is the Heaviside step function. This minimization problem is not linear in space domain, it can be linearized by using the Z transform. On Z domain $S_{mod}(z, \omega)$ is an infinite impulse response filter that can be written as:

$$S_{mod}(z, \omega) = \frac{A(\omega)}{1 - p(\omega)z^{-1}} \quad (18)$$

Form this expression, one can see that the optimization problem comes to finding the complex values of $A(\omega)$ and $p(\omega)$ that minimized the quadratic error between $S_{mod}(x_1, \omega)$ and $S_{exp}(x_1, \omega)$. The propagative properties can be extracted from $p(\omega)$ by using the following equations:

$$\left\{ \begin{aligned} k'_R(\omega) &= \frac{\arg(p(\omega))}{2\pi \cdot \Delta x_1} \\ k''_R(\omega) &= -\frac{\log(|p(\omega)|)}{\Delta x_1} \end{aligned} \right. \quad (19)$$

where $\Delta x_1$ is the distance between two successive measurement positions.

The use of the Prony algorithm on the experimental space-frequency signal allows to estimate the propagative properties of the Rayleigh wave, with these measurements and the polarization properties estimated using the 2DQFT, the algorithm of the inverse problem is complete and can be executed to estimate the complex Lamé coefficients of a viscoelastic material as a function of frequency. An experimental work has now to be performed for the validation of the proposed method.

## 5. Experimental validation

In this section, one details the experimental setup and materials used to monitor the propagation in space and time of the Rayleigh wave to extract the signal $S_{exp}(x_1, t)$. Afterword, the estimation of the Lamé coefficients from the application of the proposed method is compared with the results of other characterization methods.

A block material and a transducer are mounted in the experimental setup (Fig. 2a): the side on which the wave is propagating, is facing the vibrometer. The excitation is done by a synchronized JSR Ultrasonics DPR3 pulse generator. The time signals of the three ($u_1, u_2, u_3$) components of the particle displacement at different positions are extracted using Polytec® Laser Vibrometer PSV 500-3D-V with a sampling frequency of $f_s = 63MHz$. The measurement positions are located on a line in the sagittal plane of the wave, with a distance of $\Delta x_1 = 1.3mm$ between two successive points, with a total number of points $N_x = 62$.

A $10 \times 10 \times 3cm^3$ Epoxy block, with a density of $\rho = 1170kg/m^3$ is used as a support for the propagation of the waves (Fig. 2b). The longitudinal and transversal waves velocities are respectively $V_L = 2503m/s$ and $V_T = 1119m/s$. $k_L'' = 19.1Np/m$ and $k_T'' = 46.4Np/m$ are respectively the attenuation of longitudinal and transversal waves. These values are measured using the time signals extracted from pulse-echo method [28], the waves are generated using two different types of transducers with the same central frequency $f_c = 500kHz$, The spectral analysis of the measured temporal signal revealed a predominant frequency content centered around $f_c = 180kHz$. This frequency shift is attributed to the viscous attenuation of the epoxy, a phenomenon previously studied in the measurement of shear wave attenuation in soft tissues [29,30] and seismic attenuation characterization [31]. The measured values of $V_L, V_T, k_L''$ and $k_T''$ correspond to the central frequency of the acquired time signals. Since the pulse-echo method provides measurements only at a single frequency, these values are assumed to be constant across frequencies for simplicity. This assumption holds true for epoxy within a narrow frequency range around $f_c$. Outside this range, epoxy exhibits attenuation that increases linearly with frequency. The slope of this linear dependency is one of the characteristics of the used Epoxy [32].

The values of the complex wavenumbers $k_T^*$ and $k_L^*$ are calculated from $\omega, V_L, V_T, k_L''$ and $k_T''$, and the coefficients ($\lambda_{theo}^*, \mu_{theo}^*$) are deduced from:

$$\left\{ \begin{aligned} \lambda_{theo}^* &= \frac{\omega^2 \rho}{k_L^{*2}} - 2\mu_{theo}^* \\ \mu_{theo}^* &= \frac{\omega^2 \rho}{k_T^{*2}} \end{aligned} \right. \quad (20)$$

The surface wave is generated using a shear wave transducer mounted on the side of the block, with a central frequency of $f_c = 500kHz$ z (Fig. 2b). The transducer applies a force tangent to the surface, this will generate a Rayleigh wave in the adjacent face [33], the signals are
<figure><img src="image_4.png" /><figcaption>Fig. 2. a) Experimental setup b) epoxy block used in the experiment.</figcaption></figure>

extracted following a line in the sagittal plane of the wave, Fig. 2b. The two components of the measured space-time signal $S_{exp}(x_1, t)$ are shown in the Fig. 3.

One can see a surface wave propagating but the elliptic polarization of this last is note obvious. In order to visualize this polarization, the time signals measured at position $x_{1i} = 21.7mm$ are plotted versus the time, and then $u_3$ as a function of $u_1$ (Fig. 4).

In the time signal, the phase delay between the two components can be noted, given its unique elliptic trajectory of the particle motion. This last can be visualized in the second graph ($u_3$ versus $u_1$).

The geometric properties of this ellipse ($\chi_{exp}$, $\theta_{exp}$) parameters are of interest: their estimation as a function of frequency can be performed using the 2DQFT applied on $S_{exp}(x_1, t)$.

As expressed in Eq. (14), the Euler polar former of the 2DQFT give direct access to the polarization parameter as function of the frequency and wavenumber. Fig. 5 shows $\chi_{exp}(k_1, \omega)$, $\theta_{exp}(k_1, \omega)$ and $a_{exp}(k_1, \omega)$ extracted from the 2DQFT of $S_{exp}(x_1, t)$.

In the $a_{exp}$ representation, the maximum energy of the propagating modes is observed at $f_{R_{exp}} = 111kHz$. For the maximum energy the propagating wavenumber is $k_{R_{exp}} = 104m^{-1}$. From the characteristics of the Epoxy block, the theoretical value of the wavenumber of Rayleigh wave can be computed for the frequency $f_{R_{exp}} \cdot k_{R_{theo}} = 106m^{-1}$. The difference between $k_{R_{exp}}$ and $k_{R_{theo}}$ can be considered low relative to the values (<2%).

Artifact surrounding the area of interest in the images of $\chi_{exp}$ and $\theta_{exp}$ are often present, and they are related to the equations used to calculate these last [19], these equations are using the division operation, which is known to be numerically sensitive and unstable when used with small values. As a solution regarding the extraction of the values of

polarization, $a_{exp}$ can be used as a mask on $\chi_{exp}$ and $\theta_{exp}$, with this method the polarization parameters are directly retrieved for frequencies and wavenumbers that have energy, as shown in Fig. 5. For the frequencies and wave numbers with energy above 10 % of the maximum of $a_{exp}$, the $H/V$ ration are around 0.6. The values of orientation $\theta_{exp}$ are varying around 1.6rad. The extracted parameters are plotted versus the frequency in Fig. 6.

The retrieved frequencies are in the interval $f_{exp} = [57-182]kHz$. For these frequencies, the theoretical values of $\chi_{theo}$ and $\theta_{theo}$ are computed using Eqs. (2) to (7) and the mechanical properties of the Epoxy block ($\lambda_{theo}^*$, $\mu_{theo}^*$). The theoretical values are plotted and compared to the experimental ones in Fig. 6. One can see that, for the ellipticity the theoretical and experimental results varies around the value of $\chi = -0.628$, and for the orientation the values vary around 1.57rad. The relative difference between the experimental and theoretical values is plotted in the third graph, the results are in a good agreement, as the value of the error is laisse then 5%.

The estimated polarization parameters need to be completed with the propagative ones in order to execute the presented inverse problem. The Prony algorithm is applied on $S_{exp}(x_1, \omega)$ for the frequencies $f_{exp}$, the results regarding the estimation of $k_{R_{exp}}^*$ are plotted in Fig. 7.

The theoretical values of $k_R^*$ are calculated from the solution of Eq. (2), using ($\lambda_{theo}^*$, $\mu_{theo}^*$). One can see that the values of $k_R^*$ vary linearly between $52m^{-1}$ and $173m^{-1}$, and the values of $k_R^*$ are varying between $15m^{-1}$ and $49m^{-1}$. Compared to the theoretical values, a constant deviation is to be noticed for the real part of the wavenumber: this leads to the systematic constant relative error of around 5 % that one can be seen in the third graph. Regarding the imaginary part, errors up to 25 % were

<figure><img src="image_5.png" /><figcaption>Fig. 3. Measured components ($u_1, u_3$) of the space-time signal $S_{exp}(x_1, t)$.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Fig. 4. Time signal measured at position x₁ᵢ = 14.2mm.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 5. Polarization parameters ($a_{exp}$, $\chi_{exp}$, $\theta_{exp}$) as a function of frequency and wavenumber of the measured signal $S_{exp}(x_1, t)$.</figcaption></figure>

observed at the lower end of the frequency range studied. However, as the frequency increases, these errors decrease to values below 7 %.

As the four elements mandatory for the execution of the inverse problem are measured experimental, one can run now the optimization process for the estimation of the Lamé parameters as a function of frequency, starting from the polarization and propagative properties of the Rayleigh wave propagating on the block.

The results of the estimation of the complex Lamé coefficients ($\lambda_{exp}^*$, $\mu_{exp}^*$) from the ($\chi_{exp}$, $\theta_{exp}$, $k_R^*$) measurements are plotted in Fig. 8. The deviation relative to the values estimated from the results of the pulse-echo tests and equation (20) ($\lambda_{theo}^*$, $\mu_{theo}^*$) are represented by the error bars: the imaginary parts are negative as they represent losses [15]. In Fig. 8, it is the absolute value of the imaginary part that is presented. One can note that the results of the two methods have the same order of magnitudes. For the real parts, the $\lambda'$ values are around 4GPa and the $\mu'$ values are 1.5GPa. For the imaginary parts, the $\lambda''$ values are varying between [325–520]MPa and the $\mu''$ values are between [136–218]MPa. For the real parts, a good agreement is noticed between the results of the two methods. Despite a small deviation in the upper side of the frequency interval, this deviation is noticeable by the relative error between the two methods, that is around 10 % for these frequencies.

For the imaginary parts, the relative error decreases as the frequency approaches 180 kHz. This frequency corresponds to the central frequency of the temporal signals used to characterize the epoxy block earlier. Since the theoretical values were calculated for this specific frequency and assumed to be constant across frequencies, it is expected to observe such a decrease in error. This highlights the limitation of the prior assumption regarding the constancy of the theoretical values over

a broad frequency range. The interpretation of the relative error is meaningful for frequencies around 180 kHz, but conclusions cannot be drawn for frequencies outside this range. At 180 kHz, the relative error in estimating the imaginary parts of the Lamé coefficients is below 10 % for ($\lambda''$) and ($\mu''$).

These errors are due to the different deviations in the input parameters of the inverse problem algorithm. Fig. 7 has ever shown that an error about 7 % on the estimation of $k_R''$ and around 5 % on $k_R'$. Fig. 6 has illustrated that the errors on $\chi$ and $\theta$ are low, less than 5 %, but they still contribute to the overall deviation on the estimations. In addition, the errors introduced by the optimization process itself, as the problem is nonlinear and any errors introduces on the input may lead the process to a local minimum rather than a global one.

The large errors on the $\lambda''$ estimation can be related to the low sensitivity of the Rayleigh waves to this material parameter. For a 10 % relative variation on the modulus of the Rayleigh wave velocity, the value of the modules of $\lambda^*$ only undergoes a 3.5 % relative variation. Considering $\lambda''$, it only represents 10 % of the modulus, meaning that this has low influences on the wave properties. Despite all of these sources of deviation, the obtained results are still comparable and in good agreement with the results of the pulse-echo method for frequencies around 180 kHz.

The obtained results showed that with the proposed method has the potential to fully characterize viscoelastic materials, by estimating the variation of the complex Lamé coefficients as a function of frequency. From the propagative properties and the polarization properties of a surface Rayleigh wave, monitored in space and time by a 3D laser vibrometer.
<figure><img src="image_8.png" /><figcaption>Fig. 6. Extraction of $\chi$ and $\theta$ as a function of frequency and their comparison with the theoretical values.</figcaption></figure>

## 6. Discussion

In this study, the pulse-echo method was used for the initial characterization of the epoxy block, it serves as a valuable reference for comparing the results obtained from ellipsometry measurements. It is important to note that the viscoelastic properties of epoxy depend on the mixture ratios used during manufacturing and curing [32], which can introduce variability and comparisons with literature results less accurate due to differing fabrication processes. Therefore, characterizing the same block using different methods ensures the consistency and reliability of the values obtained across various characterization techniques.

Despite the sensitivity of the epoxy material properties to fabrication processes, it is still possible to compare the obtained results with those found in the literature. For the real parts of the coefficients, the work of Oral and Ekrem on the speed measurement of ultrasound propagation in epoxy provided the following reference values: $\lambda'_{\text{literature}} = 5.08$ GPa and $\mu'_{\text{literature}} = 1.63$ GPa [34]. In the study proposed in this paper, the density and the velocities of transverse and longitudinal waves are measured during the characterization tests. These experiments confirmed the order of magnitude of the properties proposed by Krautkrämer for epoxy resin [28]. Similar orders of magnitude are also reported by Royer [35]. Therefore, one can conclude that the presented results for these real coefficients ($\lambda'$ and $\mu'$) are in agreement with the literature.

For the imaginary part, a decreasing trend as a function of frequency is observed, which was also noted by Ni et al. over the same frequency range [36]. Ni et al. obtained their results using the time-temperature equivalence method, allowing them to display results over a broader frequency range than that of the excitation used. The obtained curves are known as the “Master curves” of the material. This principle assumes

that the rheological behavior of the polymer material as a function of temperature is proportional to its temporal behavior and inversely proportional to its behavior versus the frequency. By conducting tests at different temperatures within a given frequency range, it is possible to determine the material's rheological behavior (Ni et al., 2020; Royer and Valier-Brasier, 2021). In the results obtained by Ni et al., the trend is decreasing for frequencies around $10^5$ Hz, as observed in our study. The estimated values (136–218 MPa) also fall within the same order of magnitude and the range identified by Ni et al. (100–300 MPa).

In terms of accuracy, the proposed method has a 10 % error on the estimation of the complex Lamé coefficients as function of frequency, this value was based on the relative error obtained for the frequency of 180 kHz. The Master curve technique based on acoustic waves propagation on different temperatures has an error less than 1 % on the real parts of the coefficients and less than 10 % on the imaginary parts [37]. The method based on the reflection of shear waves on the bottom of the simple have shown errors that can go up to 9.2 % for the estimation of the complex shear modulus of epoxy [38], for viscoelastic materials methods based on tensile testing can give rather large errors (~100 %) compared to the methods based on wave propagation, due to the viscous attenuation of the tested materials [39,40]. A method based on the monitoring of Rayleigh wave propagation on the surface a soft material showed results with errors varying from 1 % up to 12 % [7]. The dynamic mechanical analysis method showed results with less than 16 % errors on the estimation of the storage and modulus of viscoelastic materials [41].

In terms of accuracy, the proposed method demonstrates a 10 % error in estimating the complex Lamé coefficients as a function of frequency, based on the relative error observed at 180 kHz. The Master curve technique, which uses acoustic wave propagation at different
<figure><img src="image_9.png" /><figcaption>Fig. 7. Estimation of $k'_R$ using Prony algorithm and comparison with theoretical values.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 8. Estimation of the complex Lamé coefficients and comparison with theoretical values.</figcaption></figure>

temperatures, achieves errors of less than 1 % for the real parts of the coefficients and less than 10 % for the imaginary parts [37]. Methods relying on shear wave reflection at the bottom of samples have shown errors of up to 9.2 % in estimating the complex shear modulus of epoxy [38]. For viscoelastic materials, methods based on tensile testing can yield significantly larger errors (up to ~100 %) compared to wave propagation methods, primarily due to viscous attenuation of the

materials [39,40]. A method based on monitoring Rayleigh wave propagation on the surface of soft materials reported errors ranging from 1 % to 12 % [7]. The dynamic mechanical analysis method has shown errors of less than 16 % in estimating the storage and loss moduli of viscoelastic materials [41].

The proposed method falls within the same range as other methods based on wave propagation. In the literature, the complex shear
modulus is often measured as a function of frequency due to its direct relation to shear wave propagation. Extracting material properties from longitudinal waves requires prior knowledge of the shear modulus, necessitating a dual characterization process using both wave types. The method presented here achieves comprehensive characterization using a single wave by leveraging the polarization of Rayleigh waves, eliminating the need for a second wave type.

Regarding cost, the method described in this study requires access to all three components of the wave, which can be achieved using a 3D laser vibrometer, such as those manufactured by Polytech®. Despite the high cost, the versatility of laser vibrometer acquisition systems allows for efficient scanning of both simple and curved surfaces, significantly saving time in measurements. In terms of numerical cost, the proposed approach is comparable to other methods because the 2D Quaternion Fourier transform is mainly implemented using the FFT algorithm, known for its optimal performance in calculating Fourier transforms.

## 7. Conclusion

In this paper, a method for the characterization of viscoelastic materials was theoretically introduced and validated experimentally. Different application fields of the viscoelastic materials were presented and the classical techniques for the estimation of the different parameters of the materials were described with their limitations. The advantages of using surface Rayleigh waves as investigating tool were listed. Following that, in the theoretical background, the complex notation of the Lamé coefficient, the equation of Rayleigh wave allowing the estimation of the complex wavenumber, and the expressions of the particle displacement from which the polarization parameters are extracted, and the inverse problem is then presented allowing the estimation of the complex lamé coefficients from the surface wave properties.

The quaternion space and the 2DQFT were introduced as a method to estimate the polarization properties of the waves as a function of frequency. The Prony algorithm and its theoretical formulation are presented as tools to estimate the propagative properties of the surface wave.

The experimental validation was started by introducing the results obtained from the pulse-echo method, the technique that uses a shear and pressure wave to characterize the viscoelastic mediums. The results of this technique are used for validation. The experimental setup was presented including its characteristics. The studied material, a block of Epoxy, was also presented, and the extracted bivariant space-time signal. The polarization properties extracted from it using 2DQFT are explained. The comparison of these results with the theoretical values estimated from the results of the pulse-echo tests have been shown in a good agreement (relative error <5 %). The estimation of the propagative properties of the surface wave was then proposed from the application of the Prony algorithm. The comparison of these lasts with the pulse-echo tests results showed a relative error <7 % on the real and imaginary parts of the wavenumber.

Despite these values, the results obtained for the complex Lamé coefficients from the inverse problem were comparable to those obtained from pulse-echo tests, with errors less than 10 %. A brief explanation of this error is provided by studying the sensitivity of Rayleigh waves to the modulus of $\lambda^*$. The proposed method was compared with existing methods in terms of accuracy and cost. In terms of accuracy, the proposed method falls within the same range as other methods based on wave propagation and controlled mechanical deformations.

These results show the ability of the proposed method to estimate the complex Lamé coefficients as a function of frequency. In our case, the method is not completely contact laisse, as the generation of the wave is still done with contact transduces, but this aspect can be changed to be contactless with the use of power lasers to generate the waves.

The presented approach has the advantage of being independent of the rheological model used, meaning that you can run the characterization method, and from the obtained results as a function of frequency

one can choose the adequate model. The propagation of Rayleigh waves on different configurations of materials (multilayers, orthotropic, composites), has been already studied and theoretically modularized, using these models, the proposed method can be adapted to this configuration.

## CRediT authorship contribution statement

**Aziz Bouzzit:** Writing – original draft, Software, Data curation. **Loïc Martinez:** Writing – review & editing, Validation, Supervision. **Andres Arciniegas:** Writing – review & editing, Validation, Supervision. **Stéphane Serfaty:** Writing – review & editing, Validation, Supervision, Project administration. **Nicolas Wilkie-Chancellier:** Writing – review & editing, Validation, Supervision, Project administration.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Data availability

Data will be made available on request.

## Acknowledgements

This work is supported and funded by CY Cergy Paris University.

## References

[1] Viscoelastic LR. Materials. Viscoelastic Materials. Cambridge University Press; 2009. p. 377-440.

[2] Xie X, Zhao G, Zhang C, Tang J, Zhou X, Jian X, et al. An innovative tensile test method to evaluate the effect of the loading rate on viscoelastic interfaces. Eng Fract Mech 2022;276:108872. https://doi.org/10.1016/j. engfracmech.2022.108872.

[3] Trivedi AR, Hawkins N, Siviour CR. Optimising dynamic mechanical analysis experiments on soft rubbers for use in time temperature superposition. MethodsX 2022;9:101831. https://doi.org/10.1016/j.mex.2022.101831.

[4] Walters K, Jones WM. Measurement of Viscosity. Instrumentation Reference Book, Elsevier; 2010, p. 69-75. https://doi.org/10.1016/B978-0-7506-8308-1.00007-3.

[5] Chen S, Urban M, Pislaru C, Kinnick R, Yi Zheng, Aiping Yao, et al. Shearwave dispersion ultrasound vibrometry (SDUV) for measuring tissue elasticity and viscosity. IEEE Trans Ultrason Ferroelectr Freq Control 2009;56:55-62. https:// doi.org/10.1109/TUFFC.2009.1005.

[6] Urban MW, Pislaru C, Nenadic IZ, Kinnick RR, Greenleaf JF. Measurement of viscoelastic properties of in vivo swine myocardium using lamb wave dispersion ultrasound vibrometry (LDUV). IEEE Trans Med Imaging 2013;32:247-61. https:// doi.org/10.1109/TMI.2012.2222656.

[7] Kazemirad SK, Heris H, Mongeau L. Experimental methods for the characterization of the frequency-dependent viscoelastic properties of soft materials. J Acoust Soc Am 2013;133:3186-97. https://doi.org/10.1121/1.4798668.

[8] Royer D, Chenu C. Experimental and theoretical waveforms of Rayleigh waves generated by a thermoelastic laser line source. Ultrasonics 2000;38:891-5. https:// doi.org/10.1016/S0041-624X(00)00022-6.

[9] Apetre N, Ruzzene M, Jacobs LJ, Qu J. Measurement of the Rayleigh wave polarization using 1D Laser vibrometry. NDT E Int 2011;44:247-53. https://doi. org/10.1016/j.ndteint.2010.07.007.

[10] Scholte JG. On Rayleigh waves in visco-elastic media. Physica 1947;13:245-50. https://doi.org/10.1016/0031-8914(47)90083-9.

[11] Sharma MD. Propagation of two Rayleigh waves in viscoelastic medium: Explicit expressions for complex velocities. ZAMM - Journal of Applied Mathematics and Mechanics / Zeitschrift Für Angewandte Mathematik Und Mechanik 2023;103. https://doi.org/10.1002/zamm.202200545.

[12] Boiko AV, Kulik VM, Seoudi BM, Chun HH, Lee I. Measurement method of complex viscoelastic material properties. Int J Solids Struct 2010;47:374-82. https://doi. org/10.1016/j.ijsolstr.2009.09.037.

[13] Royston TJ, Dai Z, Chaunsali R, Liu Y, Peng Y, Magin RL. Estimating material viscoelastic properties based on surface wave measurements: A comparison of techniques and modeling assumptions. J Acoust Soc Am 2011;130:4126-38. https://doi.org/10.1121/1.3655883.

[14] Alfrey T, Doty P. The methods of specifying the properties of viscoelastic materials. J Appl Phys 1945;16:700-13. https://doi.org/10.1063/1.1707524.

[15] Currie PK, Hayes MA, O'Leary PM. Viscoelastic Rayleigh waves. Q Appl Math 1977;35:35-53. https://doi.org/10.1090/qam/99648.
[16] Malischewsky PG, Scherbaum F. Love's formula and H/V-ratio (ellipticity) of Rayleigh waves. Wave Motion 2004;40:57-67. https://doi.org/10.1016/j.wavemoti.2003.12.015.

[17] Rose JL. Ultrasonic Guided Waves in Solid Media. Cambridge University Press; 2014. https://doi.org/10.1017/CBO9781107273610.

[18] Lagarias JC, Reeds JA, Wright MH, Wright PE. Convergence properties of the Nelder-Mead simplex method in low dimensions. SIAM J Optim 1998;9:112-47. https://doi.org/10.1137/S1052623496303470.

[19] Flamant J-M. Une approche générique pour l'analyse et le filtrage des signaux bivariés. (A general approach for the analysis and filtering of bivariate signals), 2018.

[20] Flamant J, Le Bihan N, Chainais P. Time-frequency analysis of bivariate signals. Appl Comput Harmon Anal 2019;46:351-83. https://doi.org/10.1016/j.acha.2017.05.007.

[21] Hamilton WR. Elements of Quaternions. Longmans, Green; 1866.

[22] Jamison JE. Extension of some theorems of complex functional analysis to linear spaces over the quaternions and Cayley numbers. Doctoral Dissertations 1970.

[23] Ell TA. Quaternion-Fourier transforms for analysis of two-dimensional linear time-invariant partial differential systems. In: Proceedings of 32nd IEEE Conference on Decision and Control. IEEE; 1993. p. 1830-41. https://doi.org/10.1109/CDC.1993.325510.

[24] Pei S-C, Ding J-J, Chang J-H. Efficient implementation of quaternion Fourier transform, convolution, and correlation by 2-D complex FFT. IEEE Trans Signal Process 2001;49:2783-97. https://doi.org/10.1109/78.960426.

[25] Loïc M. Nouvelles méthodes d'identification d'ondes de surface - Etude de l'onde A sur une cible courbe. Université du Havre 1998.

[26] Marple Jr SL. Digital spectral analysis. Courier Dover Publications 2019.

[27] Vollmann J, Dual J. High-resolution analysis of the complex wave spectrum in a cylindrical shell containing a viscoelastic medium. Part I. Theory and numerical results. J Acoust Soc Am 1997;102:896-908. https://doi.org/10.1121/1.419956.

[28] Krautkrämer J, Krautkrämer H. Ultrasonic Testing of Materials. Springer Berlin Heidelberg; 1990. https://doi.org/10.1007/978-3-662-10680-8.

[29] Bernard S, Kazemirad S, Cloutier G. A frequency-shift method to measure shear-wave attenuation in soft tissues. IEEE Trans Ultrason Ferroelectr Freq Control 2017;64:514-24. https://doi.org/10.1109/TUFFC.2016.2634329.

[30] Kijanka P, Urban MW. Two-point frequency shift method for shear wave attenuation measurement. IEEE Trans Ultrason Ferroelectr Freq Control 2020;67: 483-96. https://doi.org/10.1109/TUFFC.2019.2945620.

[31] Chunhua Hu, Ning Tu, Wenkai Lu. Seismic attenuation estimation using an improved frequency shift method. IEEE Geosci Remote Sens Lett 2013;10:1026-30. https://doi.org/10.1109/LGRS.2012.2227933.

[32] Rokhlín SI, Lewis DK, Graff KF, Adler L. Real-time study of frequency dependence of attenuation and velocity of ultrasonic waves during the curing reaction of epoxy resin. J Acoust Soc Am 1986;79:1786-93. https://doi.org/10.1121/1.393240.

[33] Adler L, Nagy PB. Measurements of acoustic surface waves on fluid-filled porous rocks. J Geophys Res Solid Earth 1994;99:17863-9. https://doi.org/10.1029/94JB01557.

[34] Oral I, Ekrem M. Measurement of the elastic properties of epoxy resin/polyvinyl alcohol nanocomposites by ultrasonic wave velocities. Express Polym Lett 2022;16: 591-606. https://doi.org/10.3144/expresspolymlett.2022.44.

[35] Royer D, Valier-Brasier T. Ondes élastiques dans les solides 1: Propagation. vol. 1. ISTE Group; 2021.

[36] Ni Y, Song H, Wilcox DA, Medvedev GA, Boudouris BW, Caruthers JM. Rethinking the analysis of the linear viscoelastic behavior of an epoxy polymer near and above the glass transition. Macromolecules 2020;53:1867-80. https://doi.org/10.1021/acs.macromol.9b02634.

[37] Sutherland HJ, Lingle R. An acoustic characterization of polymethyl methacrylate and three epoxy formulations. J Appl Phys 1972;43:4022-6. https://doi.org/10.1063/1.1660868.

[38] Li Y, Chang J, Huang L, Tang Y. Comparative study on viscoelastic evaluation methods of polymer materials based on ultrasonic method. Materials 2019;12: 2948. https://doi.org/10.3390/ma12182948.

[39] Sasmita F, Candra TA, Judawisastra H, Priambodo TA. Young's modulus determination of polyester and epoxy by means of ultrasonic pulse echo testing. IOP Conf Ser Mater Sci Eng 2019;547:012045. https://doi.org/10.1088/1757-899X/547/1/012045.

[40] Shtark A, Grosbein H, Sameach G, Hilton HH. An Alternative Protocol for Determining Viscoelastic Material Properties Based on Tensile Tests Without the Use of Poisson's Ratios. Volume 10: Mechanics of Solids and Structures, Parts A and B, ASMEDC; 2007, p. 437-54. https://doi.org/10.1115/IMECE2007-41068.

[41] Chakravartula A, Komvopoulos K. Viscoelastic properties of polymer surfaces investigated by nanoscale dynamic mechanical analysis. Appl Phys Lett 2006;88. https://doi.org/10.1063/1.2189156.