# Energy velocity of elastic guided waves in immersed plates for complex frequencies and slownesses

Marc Deschamps, Éric Ducasse

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

To cite this version:

Marc Deschamps, Éric Ducasse. Energy velocity of elastic guided waves in immersed plates for complex frequencies and slownesses. Wave Motion, 2025, 138, pp.103565. 10.1016/j.wavemoti.2025.103565. hal-05067658

HAL Id: hal-05067658

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

Submitted on 14 May 2025

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_2.png" /></figure>

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

ELSEVIER

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

# Energy velocity of elastic guided waves in immersed plates for complex frequencies and slownesses

Marc Deschamps<sup>*ID*</sup>, Eric Ducasse

University of Bordeaux, CNRS, Bordeaux INP, I₂M, UMR 5295, F-33400, Talence, France

Arts et Metiers Institute of Technology, CNRS, Bordeaux INP, I₂M, UMR 5295, F-33400 Talence, France

## ARTICLE INFO

**Keywords:**

Plane elastic waveguide

Inhomogeneous plane waves

Complex frequency

Energy velocity

Energy flux

## ABSTRACT

The computation of guided modes in fluid-loaded multilayer plates is generally done by a spatial approach, i.e. solutions are sought for a complex slowness. An alternative approach, less frequently employed, involves seeking solutions for complex frequencies. These frequencies correspond to plate resonances. They denote transient phenomena and the guided modes exhibit non-harmonic behavior. Consequently, conventional methods of averaging over time periods become unsuitable for calculating the means of energy quantities. In other words, the calculation of average fields cannot be reduced to a single average over a time period. To tackle this issue, for a predetermined mode, the average fields are obtained through a single averaging process applied to an arbitrary phase term. This averaging process renders independent the means of all energy quantities from the arbitrary origin phase. As usual, an additional integration across the thickness is conducted to derive total energy quantities. Doing this, the total average fields depend on both time and position on the surface plate. A set of four equations is derived from instantaneous and local energy balance equations. From these averages, the energy velocity can be directly calculated. The equations provide further insights into wave dispersion and damping along the energy flow direction, arising from viscoelastic losses and leakages in fluid.

## 1. Introduction

The objective of this study is to analyze the energy velocity of elastic guided waves with complex frequencies propagating in stratified plates immersed in fluids. From a practical perspective, this study is of significant interest in the field of non-destructive testing of composite materials, as it enables a complementary approach to the analysis of wave propagation within such structures. Indeed, by being simultaneously stratified, anisotropic, dispersive, and viscoelastic, these engineered materials exhibit a combination of characteristics that contribute to a high level of complexity in their behavior (see e.g. [1]).

It is well-established that for harmonic guided waves, i.e. real frequencies, the energy velocity is defined as the average Poynting vector divided by the average total energy, see e.g. [2]. The mean values of the necessary energetic quantities are determined by averaging over both a specific time period and the thickness of the plate. Of course, this is applicable only if the solutions of the dispersion equation of the system fluids/plate are periodic in time. Consequently, for non-harmonic solutions, an alternative description must be involved. The periodicity is lost when seeking solutions for complex frequencies. Such solutions, consistent with experimental observations, arise when the plate is immersed, giving rise to leaky waves in fluids, as illustrated in [3]. From this approach, a mode is characterized by its complex wavevector $k = \omega s$, where $\omega$ and $s$ represent the complex angular and the complex slowness vector, respectively [4–7]. To introduce the energy velocity for non-harmonic guided waves, and more generally for any

<footnote>* Corresponding author at: University of Bordeaux, CNRS, Bordeaux INP, I₂M, UMR 5295, F-33400, Talence, France.</footnote>

<footnote>E-mail addresses: marc.deschamps@u-bordeaux.fr (M. Deschamps), eric.ducasse@u-bordeaux.fr (E. Ducasse).</footnote>
wave with complex frequency and complex wavevector, energetic balance equations are written for such a context. A set of four equations is derived from instantaneous and local energy balance equations. The signal associated with a single mode, defined by a unique pair $(\omega, s)$ that satisfies the dispersion equation, can exhibit various time delays due to arbitrary phase terms. The presence of this random phase can be attributed, experimentally, to inaccuracies in establishing both the time origin and the positions of the transducers. Employing an averaging process to take advantage of the periodic nature offered by this phase, the mean behavior of a single mode is properly characterized. Beyond this definition, important links exist between the complex frequency, the complex slowness vector, and various energy fluxes.<fnref n="1" /> They provide supplementary information for a better understanding of three distinct aspects: (i) wave dispersion, (ii) the disparity between potential and kinetic mean energies, and (iii) damping along the energy flow direction, attributable to viscoelastic losses and leakage.

From this averaged system of equations, the energy velocity is defined, similarly to the standard harmonic case, as the ratio of the energy flux vector to the total energy. In the case of a non-absorbing plate in vacuum, this wave velocity corresponds to the group velocity [2]. By nature, these two velocities are distinct; therefore, it is not surprising that they may differ in general cases. In fact, the group velocity is associated with the propagation of the maximum amplitude of a signal, while the energy velocity is expressed from products of wave fields. If these two velocities are nearly identical for solutions with complex frequencies, they can, however, differ significantly for solutions characterized by complex slownesses, as will be observed. In addition, the energy balance behavior of a plate system under inhomogeneous wave excitation exhibits significant deviations from that observed in the homogeneous case. A notable example is the reflection of inhomogeneous waves by an immersed plate, where the reflection coefficient may exceed unity. This phenomenon can be explained by the existence of interactive energy exchanges between waves in the fluid medium [8].

The novelty of this study lies, first, in proposing an efficient method for calculating the average of non-periodic waves and, second, in utilizing this approach to compare group velocity and energy velocity for complex frequencies or slowness. The paper is organized as follows. Section 2 succinctly introduces guided waves with complex frequency and complex slowness. In Section 3, these general guided waves, which nonetheless exhibit amplitude and front phase planes, are subsequently integrated into the fundamental energy balance equations. This process yields the four instantaneous and local equations. To fully exploit these results, a procedure for calculating the averages of energy quantities is proposed in Section 4. The application of this new definition, tailored to transient waves, and its rich consequences are presented in Section 5. Finally, numerical results are presented and discussed in Section 6.

## 2. Guided waves for complex frequency and complex slowness

In order to establish our notations, let us first describe the problem. The geometry, as defined within the 3D fixed coordinate system $(O, n_x, n_y, n_z)$, comprises a plate with a thickness $\ell$ in the z-direction and infinite extent in the xy-plane. This infinite plate separates two different surrounding fluids, as illustrated in Fig. 1. The observation point, denoted as $M$, is situated in such a way that $\overline{OM} = z = x + z n_z$, where its projection onto the xy-plane is represented by the vector $x = x n_x + y n_y$.

As is customary, the superposition of all partial waves within the plate must satisfy the boundary conditions. Specifically, the normal stress and normal displacements must be continuous at both interfaces. These interfaces serve as the connectors between the guided wave propagating within the plate and the two waves propagating within the adjacent fluid media.

Owing to the fluid coupling, solutions defined by the frequency/wavevector pair $(\omega, k)$ can be both complex [3]. Indeed, real-valued wavevector and frequency parameters correspond exclusively to non-attenuated waves, such as SH guided waves in elastic plates that do not radiate into surrounding fluids. Consequently, in the most comprehensive scenario involving attenuated waves, the boundary conditions can only be satisfied when at least one of $k$ or $\omega$ is complex-valued. The guided wave is assumed to be a summation of homogeneous or inhomogeneous plane waves at the same complex frequency. In accordance with the generalized Snell's law, all of these elementary waves manifest a common complex wavevector at the interface, for which the associated phase front is chosen to always be oriented in the direction of the vector $n_z$ within the xy-plane. The most common scenario for a harmonic leaky guided wave corresponds to a real-valued frequency and a complex wavevector $k = k' + i k''$. The superscripts prime and double prime denote, and will denote, the real and imaginary parts of a complex variable or expression, respectively. To encompass the entire set of solutions, it is imperative to consider complex frequencies. As a result, both the wavevector and frequency are presumed to be complex. The real wavevector, denoted as $k'$, is consistently oriented along the unit vector $n_x$, with no constraints imposed on the imaginary part $k''$ of the wavevector. It is noteworthy that, aside from the symmetry conditions, no specific assumptions are imposed on the coefficients of the elasticity tensor. Consequently, general anisotropy can be thoroughly investigated by simply rotating the crystallographic axes, which correspond to variations in the propagation direction of guided waves (see Fig. 1).

The complex angular frequency is expressed by $\omega = \omega' + i\omega''$. The parameter $\omega''$ explains the transient behavior of the signal associated with the guided waves. It is predominantly positive, although in certain situations, it exhibits negativity [3]. Through this representation, the guided wave can be characterized at any spatial point $M$ and time $t$ by using the subsequent real velocity field, denoted by $v$, which represents the real part of the following complex vector field $v$:

$$v(x, y, z, t) = \hat{v}(z) \exp [i(\omega t - \mathbf{k} \cdot \mathbf{x} + \psi)]. \quad (1)$$

<footnote n="1">¹ One of which is the Poynting vector.</footnote>
<figure><img src="image_5.png" /><figcaption>Fig. 1. The plate in the fixed coordinate system (O, **n**<sub>x</sub>, **n**<sub>y</sub>, **n**<sub>z</sub>), along with the slowness vector **s** and the wavevector **k** associated with the guided wave.</figcaption></figure>

Then, the real velocity field takes the form:

$$
\mathbf{v}(x, y, z, t) = [\hat{\mathbf{v}}'(z) \cos(\phi') - \hat{\mathbf{v}}''(z) \sin(\phi')] \exp(-\phi'') \quad \text{with} \quad \begin{cases} \phi' = \omega' t - k_x' x + \psi, \\ \phi'' = \omega'' t - k_x'' x. \end{cases} \tag{2}
$$

The term $\psi$ lies within the interval $(-\pi, \pi]$. It represents the arbitrary phase of the complex amplitude of the wave. This phase will subsequently be utilized to determine average values, and its nature as well as its impact on real fields will be further discussed in Section 4. To maintain clarity in our text, complex quantities will henceforth be denoted by underlining or by the use of a hat over a symbol. It is important to note that although the phase $\psi$ is arbitrary, it must remain consistent across all fields to satisfy the Snell–Descartes laws.

Likewise, the real strain and stress fields are expressed, respectively, as follows:

$$
\mathcal{E}(x, y, z, t) = \left\{ \hat{\mathcal{E}}(z) \exp \left[ \hat{i} \left( \frac{\omega}{t} - \frac{\mathbf{k}}{\mathbf{k} \cdot \mathbf{x}} \cdot \mathbf{x} + \psi \right) \right] \right\}' \quad \text{and} \quad \mathcal{S}(x, y, z, t) = \left\{ \hat{\mathcal{S}}(z) \exp \left[ \hat{i} \left( \frac{\omega}{t} - \frac{\mathbf{k}}{\mathbf{k} \cdot \mathbf{x}} \cdot \mathbf{x} + \psi \right) \right] \right\}'.
$$

As shown in [3] in the investigation of Lamb waves propagating through isotropic plates immersed in a fluid, there is a noteworthy interest in the introduction of the complex slowness bi-vector, denoted as **s**, which is defined by [4]:

$$
\underline{s} = s' + i s'' = \underline{\tau} \underline{k}, \quad \text{with } \underline{\tau} = \tau' + i \tau'' = \underline{\omega}^{-1}.
$$

In doing so, the slowness vector exhibits spatial wave structures, as it is decoupled from the temporal behavior. This relationship may appear trivial, but its significance lies in the fact that it removes any ambiguity regarding the origins of the imaginary part of the wavevector. In fact, the slowness vector depends solely on quantities related to rheological phenomena, i.e., elasticity or viscoelasticity, and spatial structure, such as illustrated by evanescent surface waves, for example. The frequency solely characterizes temporal effects, despite its influence on the imaginary part of the wavevector. To illustrate this point, consider the case of a wave with $s'' = 0$. In this case, it is evident that the wave exhibits despite everything a complex wavevector $\omega s'$ while any field $\hat{v}(z) \exp[i\underline{\omega}(t-s'\cdot\mathbf{x})+i\underline{\psi}]$ is a function of $(t-s'\cdot\mathbf{x})$ that corresponds to translational motion at velocity $|s'|^{-1}$ in the $s'$ direction.

**3. Basic equations**

In every layer, at a constant complex frequency, the inhomogeneous plane wave described by Eq. (1) complies with the fundamental principle of dynamics, which is given by: $\hat{i} \underline{\omega} \rho \underline{v} = \hat{i} \underline{\omega} [\mathbf{I} \underline{\Sigma}]^T \underline{s} + \partial_z \underline{\sigma}_z$, where $[\bullet]^T$ represents the transposition operation. This equation can be rewritten in the form:

$$
\rho \underline{v} = [\mathbf{I} \underline{\Sigma}]^T \underline{s} - i \underline{\tau} \partial_z \underline{\sigma}_z, \quad (3)
$$

where $\rho$ stands for the mass density. The multiplication by the **I** matrix<sup>2</sup> selects the projection on the $xy$-plane and transforms this projection into a well-oriented Poynting vector in the broad sense [2]. The wave damping resulting from viscoelastic behavior

<footnote>$^2$ $\mathbf{I} = \begin{pmatrix} -1 & 0 & 0 \\ 0 & -1 & 0 \\ 0 & 0 & 0 \end{pmatrix}$.</footnote>
is considered in our analysis. Various models are available to describe absorbing media. One of them characterizes the dynamic behavior of the hysteresis-shaped stress-strain relationship using the complex tensor of viscoelasticity C, defined as follows [2]:

$$ \underline{C} = C' + iC'', \qquad (4) $$

where $C'$ and $C''$ are the elastic and viscosity contributions to this tensor. With this model, each tensor component does not depend on the frequency, which is clearly not the case for general viscoelastic behavior [9]. However, for specific materials and frequency ranges, it is in concordance with experimental observations. In no way, this model can be valid at any frequency. Consequently, to enhance the validity of models, it is imperative to ensure that the imaginary part of each constant varies with respect to frequency. The use of complex frequencies would require additional efforts, beyond the scope of this article. In fact, to our knowledge, analyzing in detail the consequences of such a generalization remains an open problem.

With the behavior described by Eq. (4), the stress tensor is:

$$ \underline{\Sigma} = C' : \underline{\mathcal{E}} + i C'' : \underline{\mathcal{E}}. \qquad (5) $$

This equation can be examined using complex variables, as they establish a connection between two individual fields, $\underline{\Sigma}$ and $\underline{\mathcal{E}}$ in the case of Eq. (5). This does not hold true for the physical quantities associated with the energy equations, as they invariably arise from an interaction between two fields, as will be expressed by Eq. (18). In this physical context, the equations presented in the following sections will generally be real-valued, with only a few exceptional cases warranting complex values. Over the next three subsections, we will establish relationships between energy fluxes and the slowness vector for inhomogeneous waves with complex frequency. It is important to keep in mind that the equations derived throughout Section 3 are dependent on both space and time. To well define the energy velocity for non-harmonic guided waves, and more generally for any complex wavenumbers, the energetic balance equations are obtained for such circumstances.

Anticipating future calculations, let us introduce two complex vectors denoted as $\underline{p}_+ = \underline{I} \underline{\Sigma} \underline{v}^*$ and $\underline{p}_- = \underline{I} \underline{\Sigma} \underline{v}$, whose the real and imaginary parts signify the projection of four energy fluxes along the interfaces, where the superscript * denotes the complex conjugate operator. These two vectors are expressed as follows:

$$ \underline{p}_{\pm} = \underline{p}' \pm \underline{q}'' + i (\underline{p}'' \mp \underline{q}') . \qquad (6) $$

The four real vectors introduced in this definition are such that: $\underline{p}' = \underline{I} \underline{\Sigma}' \underline{v}'$, $\underline{p}'' = \underline{I} \underline{\Sigma}'' \underline{v}'$, $\underline{q}' = \underline{I} \underline{\Sigma}' \underline{v}''$ and $\underline{q}'' = \underline{I} \underline{\Sigma}'' \underline{v}''$. From a physical point of view, although these vectors may exhibit similarities with instantaneous Poynting vectors, the authentic Poynting vector is unambiguously issued from the vector $\underline{p}'$.

### 3.1. Instantaneous and local kinetic energy

First, let us examine the instantaneous volume density of kinetic energy. By definition, it is expressed as follows:

$$ e_k = \frac{1}{2} \rho \underline{v}' \cdot \underline{v}'. \qquad (7) $$

By virtue of Newton's second law (3), the volumetric density of kinetic energy can be redefined by expressing one of the velocity vectors $\underline{v}$ as a function of the stress tensor. This gives rise to $e_k = \frac{1}{2} \{\underline{s} \cdot \underline{I} \underline{\Sigma}\}' \cdot \underline{v}' + \frac{1}{2} \{\underline{\tau} \cdot \underline{\partial}_z \underline{\sigma}_z\}'' \cdot \underline{v}'.$ From this equation and employing the definitions provided in Eq. (6), the solution presented in Eq. (2) leads to the following expression for its kinetic energy:

$$ e_k = \frac{1}{2} (\underline{s}' \cdot \underline{p}' - \underline{s}'' \cdot \underline{p}'') + \tau' a_k + \tau'' b_k, \qquad (8) $$

where $a_k = \frac{1}{2} (\underline{\partial}_z \underline{\sigma}_z'') \cdot \underline{v}'$ and $b_k = \frac{1}{2} (\underline{\partial}_z \underline{\sigma}_z') \cdot \underline{v}'.$ To extract additional information from Eq. (3), an insightful statement can be obtained by calculating from this equation the complex scalar product $\rho \underline{v} \cdot \underline{v}'.$ So, it is obtained:

$$ \rho \underline{v} \cdot \underline{v}' = (\underline{I} \underline{\Sigma})^T \underline{s} - i \underline{\tau} \cdot \underline{\partial}_z \underline{\sigma}_z \cdot \underline{v}' . \qquad (9) $$

It is worth emphasizing that due to the complex nature of frequency, this term does not depict the averaged kinetic energy over a period, as typically considered in the context of harmonic waves. This is for the simple reason that the field under consideration is by choose non-periodic. At this stage, it is merely a mathematical expression; however, an adequate physical interpretation will be provided later on. The expression on the left-hand side of Eq. (9) is inherently real-valued by design. It is immediately deduced that $\{((\underline{I} \underline{\Sigma})^T \underline{s} - i \underline{\tau} \cdot \underline{\partial}_z \underline{\sigma}_z) \cdot \underline{v}'\}'' = 0$. By using the flux quantities introduced in Eq. (6) leads to:

$$ s'' \cdot \underline{p}' + s' \cdot \underline{p}'' - s' \cdot \underline{q}' + s'' \cdot \underline{q}'' + 4 \tau' c_k + 4 \tau'' d_k = 0, \qquad (10) $$

with $c_k = -\frac{1}{4} [(\underline{\partial}_z \underline{\sigma}_z') \cdot \underline{v}' + (\underline{\partial}_z \underline{\sigma}_z'') \cdot \underline{v}'']$ and $d_k = \frac{1}{4} [(\underline{\partial}_z \underline{\sigma}_z'') \cdot \underline{v}' - (\underline{\partial}_z \underline{\sigma}_z') \cdot \underline{v}'']$. The apparently artificial introduction of the factor $\frac{1}{4}$ is guided by the fact that the balance equations presented later by Eq. (25) must contain directly interpretable average energy terms with physical significance.
## 3.2. Instantaneous and local potential energy

The instantaneous volume density of total potential energy $e_t$ is expressed as:

$$e_t = \frac{1}{2} \Sigma' : \mathcal{E}' = e_p - e_q, \qquad (11)$$

where the newly introduced real terms are such that: $e_p = \frac{1}{2} \mathcal{E}' : C' : \mathcal{E}'$ and $e_q = \frac{1}{2} \mathcal{E}'' : C'' : \mathcal{E}'$. They represent, respectively, the elastic potential energy and the dissipative energy losses. For further developments, let us introduce the following property:

$$\mathbf{A} : \mathcal{E} = (\mathbf{I} \mathbf{A}) : (\mathbf{s} \otimes \mathbf{v}) - i \underline{\tau} (\mathbf{A} \mathbf{n}_z) \cdot \partial_z \mathbf{v}, \qquad (12)$$

which connects any second-order symmetric tensor A with the strain tensor for the solution given by Eq. (2). Exploiting this relation to the total potential energy expressed by Eq. (11), it is obtained:

$$e_t = \frac{1}{2} (\mathbf{I} \Sigma') : (\mathbf{s}' \otimes \mathbf{v}' - \mathbf{s}'' \otimes \mathbf{v}'') + \tau' a_p + \tau'' b_p,$$

where $a_p = \frac{1}{2} \boldsymbol{\sigma}'_z \cdot (\partial_z \mathbf{v}'')$ and $b_p = \frac{1}{2} \boldsymbol{\sigma}'_z \cdot (\partial_z \mathbf{v}')$. Expanding this last equation for the solution given by Eq. (2) and using relations (6) leads to:

$$e_p = \frac{1}{2} (\mathbf{s}' \cdot \mathbf{p}' - \mathbf{s}'' \cdot \mathbf{q}') + \tau' a_p + \tau'' b_p + e_q. \qquad (13)$$

An interesting equation can be obtained from the calculation of the real-valued scalar expression $\underline{\mathcal{E}} : C' : \underline{\mathcal{E}}^*$. Using the stress expression from Eq. (5), it is then derived the supplementary equation:

$$\underline{\mathcal{E}} : C' : \underline{\mathcal{E}}^* = \underline{\mathbf{S}}^* : \underline{\mathcal{E}} + 4i e_d, \qquad (14)$$

where $e_d = \frac{1}{4} \underline{\mathcal{E}} : C'' : \underline{\mathcal{E}}^*$. Since the scalar $e_d$ is a real and positive quantity [10], the imaginary part of the right-hand side of Eq. (14) is necessarily zero. Therefore, using both Eq. (12) and the variables given in Eq. (6), the subsequent equation must be satisfied:

$$\mathbf{s}'' \cdot \mathbf{p}' - \mathbf{s}' \cdot \mathbf{p}'' + \mathbf{s}' \cdot \mathbf{q}' + \mathbf{s}'' \cdot \mathbf{q}'' + 4 \tau' c_p + 4 \tau'' d_p + 4 e_d = 0, \qquad (15)$$

$$\text{with } c_p = -\frac{1}{4} [\boldsymbol{\sigma}'_z \cdot (\partial_z \mathbf{v}') + \boldsymbol{\sigma}''_z \cdot (\partial_z \mathbf{v}'')] \text{ and } d_p = \frac{1}{4} [\boldsymbol{\sigma}'_z \cdot (\partial_z \mathbf{v}'') - \boldsymbol{\sigma}''_z \cdot (\partial_z \mathbf{v}')].$$

## 3.3. Instantaneous and local balance energy equations

We will now proceed to employ the equations outlined in the preceding two subsections. By performing algebraic manipulations, the following set of four equations can easily be established:

$$\begin{align}
e_k + e_p &= \mathbf{s}' \cdot \mathbf{p}' - \frac{1}{2} \mathbf{s}'' \cdot \mathbf{p}'' + \tau' a_+ + \tau'' b_+ + e_q, \tag{a} \\
e_k - e_p &= -\frac{1}{2} \mathbf{s}'' \cdot \mathbf{p}''_+ + \tau' a_- + \tau'' b_- - e_q, \tag{b} \\
0 &= \frac{1}{2} \mathbf{s}'' \cdot \mathbf{p}''_+ + \tau' c_+ + \tau'' d_+ + e_d, \tag{c} \\
0 &= \frac{1}{2} \mathbf{s}' \cdot \mathbf{p}''_+ + \tau' c_- + \tau'' d_- - e_d. \tag{d}
\end{align}$$

The coefficients $a_\pm$, $b_\pm$, $c_\pm$ and $d_\pm$ are defined with respect to vertical flux quantities, i.e., oriented along the z-direction, and are such that: $a_\pm = a_k \pm a_p$, $b_\pm = b_k \pm b_p$, $c_\pm = c_k \pm c_p$ and $d_\pm = d_k \pm d_p$, with in particular the important following factorizations: $b_+ = \frac{1}{2} \partial_z (\boldsymbol{\sigma}'_z \cdot \mathbf{v}')$, $c_+ = -\frac{1}{4} \partial_z \{\boldsymbol{\sigma}_z \cdot \mathbf{v}^*\}'$ and $d_- = -\frac{1}{4} \partial_z \{\boldsymbol{\sigma}_z \cdot \mathbf{v}^*\}''$. The last three relations are of interest because, as will be shown, they will offer a means to simplify the calculations of averages introduced latter. The first two equations are derived by adding and subtracting Eqs. (8) and (13), respectively. The last two equations are obtained by performing summation and subtraction operations on Eqs. (10) and (15), respectively. All of these equations depend on both time and space. They represent instantaneous and local relationships between energy terms. Strictly speaking, the energy balance associated with the guided wave is determined by Eq. (16a), in which the volume density of total energy, comprising the summation of volume densities of kinetic and potential energies, is compared to various wave fluxes. The second equation, denoted as Eq. (16b), offers insights into the relative magnitudes of the kinetic and potential terms. The final two equations, Eqs. (16c) and (16d), establish a connection between the real and imaginary components of the slowness vector $\underline{s}$ and the flux quantities in the $xy$-plane.

All of these relationships pertain to local and instantaneous quantities. In order to effectively leverage them, it becomes imperative to carry out spatial and temporal averaging. This will be the focal point of the upcoming section.

## 4. Average value of the energy quantities

Indeed, owing to the complex-valued nature of frequency, it is evident that not all fields manifest both temporal and spatial harmonicity. This leads to the immediate fact that the computation of average fluxes cannot be streamlined into a singular time-averaged value, as is commonly undertaken for harmonic plane waves in conventional practices. In other words, considering the spatial and temporal characteristics of complex harmonic waves, it becomes essential to propose an alternative definition for averaging that adequately accommodates the intricacies inherent in this complexity. The phase term $\psi$ introduced in Eq. (1) plays a crucial role in defining averaged quantities. As previously mentioned, this term arises from the arbitrary complex amplitude of
the wave. From an experimental point of view, its value is influenced, for instance, by the choice of the time origin, the positions of possible transducers, and the amplitude of the source. In any case, this term is dependent on wave propagation. This becomes evident when considering that the solution $(\underline{k}, \omega)$ of the dispersion equation for the plate system is independent of the complex arbitrary amplitude and, consequently, of $\psi$. However, the instantaneous real physical quantities do depend on $\psi$, see, for example, Eq. (1). This dependence does not manifest as a simple arbitrary factor that can be factored out of these quantities, which would be necessary for its elimination in the energy balance equations. This fact necessitates the precise definition of a unique mode, characterized by a specific pair $(\underline{k}, \omega)$, while still allowing $\psi$ to take any value within the interval $(-\pi, \pi]$. This association between a single solution and arbitrary phases constitutes a unique mode. The next question is how to describe the intrinsic propagation of this mode. Inspired by signal processing techniques commonly used to analyze the stochastic nature of various phenomena (see, e.g., [11]), the proposed approach to tackle this challenge is as follows: (i) Represent a mode as a family of functions associated with physical fields, corresponding to the infinite set of possible values for $\psi$; (ii) Consider this phase term as a random variable, even though it is not strictly so, by assuming it follows a uniform distribution. Each realization corresponds to assigning an arbitrary value to $\psi$; and (iii) Derive the physical quantities for a single mode by computing their probabilistic average, specifically by performing an average over a period $T = 2\pi$. This approach extends the concept of the average value, usually defined in textbooks for harmonic waves, to inhomogeneous waves characterized by both complex frequencies and wave numbers. Naturally, this new definition is compatible with the traditional one when dealing with harmonic waves, as the phenomenon is ergodic in this case, see e.g. [11]. In other words, the expectation at any time $t_0$, defined as the probabilistic average over the phase $\psi$, is equivalent to the time average for any fixed phase value $\psi_0$. It is important to note that this average is independent of both $t_0$ and $\psi_0$.

Finally, let us express the average field $\mathbb{F}(x, y, t)$ at the position $(x, y)$ and time $t$ for any field $f(x, y, z, t)$ as follows:

$$ \mathbb{F}(x, y, t) = \langle \mathbb{f}(x, y, z, t) \rangle_{\ell} = \frac{1}{\ell} \int_{0}^{\ell} \mathbb{f}(x, y, z, t) dz, \quad \text{where} \\ \mathbb{f}(x, y, z, t) = \langle f(x, y, z, t) \rangle_{\psi} = \frac{1}{2\pi} \int_{-\pi}^{\pi} f(x, y, z, t) d\psi. \tag{17} $$

In the context of this average and to enhance readability, the double-bar capital notation will be subsequently employed to represent the average values associated with both the phase $\psi$ and the vertical $z$-position. Conversely, double-bar lowercase notation will be exclusively reserved for average values related to the phase $\psi$. The energy flux, potential energy, and kinetic energy, all denoted at this stage as $w(x, y, z, t)$, can be expressed as functions of two physical fields. These fields correspond to the real parts of two complex fields $\underline{a}(x, y, z, t)$ and $\underline{b}(x, y, z, t)$ satisfying the following relationship:

$$ w(x, y, z, t) = \{ \underline{a}(x, y, z, t) \}' \cdot \{ \underline{b}(x, y, z, t) \}', \tag{18} $$

where the bilinear operator $\cdot$ depends upon the specific energetic quantity under examination. Its specific definition will be explicitly provided as needed later on. In agreement with Eq. (2), the two fields are of the form: $\underline{a}(x, y, z, t) = \hat{a}(z) \exp(i \underline{\phi})$ and $\underline{b}(x, y, z, t) = \hat{b}(z) \exp(i \underline{\phi})$.

To calculate the integrals described by Eq. (17), it proves advantageous to decouple the dependence on the $z$-coordinate from that on $x, y$ and $t$. This can be achieved as follows:

$$ w(x, y, z, t) = g(x, y, t) \left\{ \left\{ \hat{a}(z) \cdot \hat{b}(z)^* \right\}' + \left\{ \left( \hat{a}(z) \cdot \hat{b}(z) \right) \exp(2i \underline{\phi}') \right\}' \right\}, \tag{19} $$

where $g(x, y, t) = \frac{1}{2} \exp(-2 \underline{\phi}')$. The complex phase depends only on $x, y$ and $t$, while the fields $\hat{a}$ and $\hat{b}$ are exclusively functions of $z$. Applying averages as defined in Eq. (17) to the generic fields of the form Eq. (19) gives rise to:

$$ \mathbb{W}(x, y, t) = \langle w \rangle_{\ell}, \quad \text{where} \quad w = g(x, y, t) \{ \hat{a}(z) \cdot \hat{b}(z)^* \}', \tag{20} $$

because $\underline{\phi}''$ does not depend on $\psi$ while the following average is zero:

$$ \int_{-\pi}^{\pi} \exp(2i \underline{\phi}') d\underline{\psi} = \exp[2i(\underline{\omega}'t - k_{\bar{x}}'\bar{x})] \int_{-\pi}^{\pi} \exp(2i \underline{\psi}) d\underline{\psi} = 0. $$

It is important to note that the averaged quantities $\mathbb{W}(x, y, t)$ depend on the observation point on the surface and time. This dependence is solely encapsulated by the $g$-factor. The $g$-factor is, in fact, the square of the wave amplitude at the location $(x, y)$ and time $t$. To simplify the notations, the dependence on 3D-space and time will be omitted until this is possible without compromising clarity.

Let us observe that the average value $w$ given in Eq. (20) is strictly equal to $w = \frac{1}{2} \{\underline{a} \cdot \underline{b}^*\}'$. This equality arises from the inherent nature of the $g$ factor being a real value, and the shared phase term between both vectors. In passing, let us emphasize that expressions (7) and (11) now attain a physical interpretation through averaging over $\psi$, but not over a hypothetical space or time period.

The fields that require averaging do not always conform directly with the representation established by Eq. (18). As an example, when averaging the term $p'' = I (\underline{\Sigma}' v'' + \underline{\Sigma}'' v')$ from Eq. (6), it is necessary to address terms in the forms of $\underline{a}' \cdot \underline{b}''$ and $\underline{a}'' \cdot \underline{b}'$. From Eq. (20), it is a simple matter to average these quantities and all other combinations of the real parts of vectors $\underline{a}$ and $\underline{b}$. This can be obtained by using the obvious equation $\langle \underline{a} \cdot \underline{b} \rangle_{\psi} = 0$. In fact, elementary algebraic manipulations demonstrate that:

$$ \begin{aligned} \langle \underline{a}'' \cdot \underline{b}' \rangle_{\psi} &= - \langle \underline{a}' \cdot \underline{b}'' \rangle_{\psi} = - \langle \underline{a}' \cdot \{ -i \underline{b} \}' \rangle_{\psi} = \frac{1}{2} \{ \underline{a} \cdot \underline{b}^* \}'' \\ \langle \underline{a}'' \cdot \underline{b}'' \rangle_{\psi} &= \langle \underline{a}' \cdot \underline{b}' \rangle_{\psi} = \frac{1}{2} \{ \underline{a} \cdot \underline{b}^* \}' \end{aligned} \tag{21} $$
and consequently:

$$ \begin{aligned} \langle \{a \cdot b^*\}' \rangle_{\psi} &= \langle a' \cdot b' + a'' \cdot b'' \rangle_{\psi} = 2 \langle a' \cdot b' \rangle_{\psi}, \\ \langle \{a \cdot b^*\}'' \rangle_{\psi} &= \langle a'' \cdot b' - a' \cdot b'' \rangle_{\psi} = 2 \langle a'' \cdot b' \rangle_{\psi}. \end{aligned} \qquad (22) $$

It is noteworthy that, for harmonic guided waves only, the averages over $\psi$ are strictly equivalent to averages over a time period, since the ergodicity property defined above is satisfied. In the same way, for waves periodic in space, i.e., of wavevector with real-valued components, the average in space over a wavelength coincides with the stochastic average over $\psi$. Consequently, the stochastic averaging defined here is more general than the usual time or space averaging, but is in accordance with the literature for harmonic or spatially-periodic case. Note that the issue of space averaging applies not only to waves but also to the mechanical properties of periodic media in determining equivalent homogeneous quantities [12].

## 5. Averaged energy equations

The average definition summarized by Eq. (17) will be used in this section to obtain averaged expressions of equations relative to instantaneous fields obtained in Section 3 for general complex guided waves. To achieve this, let us express the averages over $\psi$ of the variables involved in the flux vectors defined by Eq. (6). In this case, the indefinite fields $\mathbf{a}$ and $\mathbf{b}$ are replaced by $\mathbf{I}\Sigma$ and $\mathbf{v}$ (or $\mathbf{v}^*$), respectively, and the linear operator $\mathbf{\cdot}$ by the matrix-vector product. So, employing the relationships provided in Eq. (21) results in the following average over $\psi$ of these four flux vectors:

$$ \mathbf{P}' = \mathbf{q}'' = \frac{1}{2} \mathbf{P}'_+ = g \mathbf{I} \{ \hat{\mathbf{\Sigma}} \hat{\mathbf{v}}^* \}' , \quad \mathbf{P}'' = -\mathbf{q}' = \frac{1}{2} \mathbf{P}''_+ = g \mathbf{I} \{ \hat{\mathbf{\Sigma}} \hat{\mathbf{v}}^* \}'' , \qquad (23) $$

and, consequently, $\mathbf{P}' = 0$ and $\mathbf{P}'' = 0$. Let us note that the average over $\psi$ of the viscoelastic contribution $e_q$, introduced in Eq. (13), is zero, because:

$$ \mathbf{\varepsilon}_q = \frac{1}{2} \langle \mathbf{\varepsilon}'' : \mathbf{C}'' : \mathbf{\varepsilon}' \rangle_{\psi} = \frac{1}{4} \{ \mathbf{\varepsilon} : \mathbf{C}'' : \mathbf{\varepsilon}^* \}'' = 0. $$

For its part, the loss of energy $e_d$, introduced in Eq. (14), does not depend on the phase $\psi$. As a consequence, this energy loss can be averaged as follows:

$$ \mathbf{\varepsilon}_d = \frac{1}{4} \langle \mathbf{\varepsilon} : \mathbf{C}'' : \mathbf{\varepsilon}^* \rangle_{\psi} = \frac{g}{2} \langle \hat{\mathbf{\varepsilon}} : \mathbf{C}'' : \hat{\mathbf{\varepsilon}}^* \rangle_{\psi} . \qquad (24) $$

By employing the double-bar notation as introduced in Eqs. (17) to define mean quantities, and utilizing the basic expressions for specific averages provided by Eqs. (21) and (22), along with the assistance of Eqs. (23) through (24), the balance equations outlined in Eq. (16) can also be integrated. In doing so, the associated relationships can be reformulated in terms of scalar products between the complex slowness vector and the complex mean Poynting vector, as follows:

$$ \begin{aligned} s' \cdot \mathbf{P}' &= \mathbf{E}_{tot} - \tau' \mathbf{F}_{in}'' - \tau'' \mathbf{F}_{out}', & (a) \\ s'' \cdot \mathbf{P}'' &= \mathbf{E}_{p-k} + \tau' \mathbf{F}_{out}'' + \tau'' \mathbf{F}_{in}', & (b) \end{aligned} \qquad (25) $$

$$ s'' \cdot \mathbf{P}' = \tau' \mathbf{F}_{out}'' - \tau'' \mathbf{F}_{in}'' - \mathbf{E}_d, \qquad (c) $$

$$ s' \cdot \mathbf{P}'' = \tau' \mathbf{F}_{in}'' - \tau'' \mathbf{F}_{out}'' + \mathbf{E}_d, \qquad (d) $$

where $\mathbf{E}_{tot} = \mathbf{E}_p + \mathbf{E}_k$ and $\mathbf{E}_{p-k} = \mathbf{E}_p - \mathbf{E}_k$ represent, respectively, the mean summation and the difference between the potential and kinetic mean energies. The quantities $\mathbf{E}_k$, $\mathbf{E}_p$ and $\mathbf{E}_d$ are the mean kinetic and potential energies and the mean loss of energy, respectively. Reintroducing the z-dependence, they take the form: $\mathbf{E}_p = \frac{g}{2} \langle \hat{\mathbf{\varepsilon}}(z) : \mathbf{C}'(z) : \hat{\mathbf{\varepsilon}}^*(z) \rangle_{\ell}$, $\mathbf{E}_k = \frac{g}{2} \langle \rho(z) \hat{\mathbf{v}}(z) \cdot \hat{\mathbf{v}}^*(z) \rangle_{\ell}$ and $\mathbf{E}_d = \frac{g}{2} \langle \hat{\mathbf{\varepsilon}}(z) : \mathbf{C}''(z) : \hat{\mathbf{\varepsilon}}^*(z) \rangle_{\ell}$. The average of the eight coefficients $a_{\pm}$, $b_{\pm}$, $c_{\pm}$ and $d_{\pm}$ defined in Eq. (16), when combined with the relationships provided by Eq. (21), results in only two independent variables $\mathbf{F}_{in}$ and $\mathbf{F}_{out}$ such that:

$$ \underline{\mathbf{E}}_{in} = \frac{g}{4\ell} \langle (\partial_z \hat{\mathbf{\sigma}}_z) \cdot \hat{\mathbf{v}}^* - \hat{\mathbf{\sigma}}_z \cdot (\partial_z \hat{\mathbf{v}}^*) \rangle_{\ell} \quad \text{and} \quad \underline{\mathbf{E}}_{out} = \frac{1}{2\ell} [\underline{\mathbf{P}}_{-z}(\ell) - \underline{\mathbf{P}}_{-z}(0)], $$

where the vector $\underline{\mathbf{P}}_{-z}(z) = \langle -\frac{1}{2} \hat{\mathbf{\sigma}}_z(z) \cdot \underline{\mathbf{v}}^*(z) \rangle_{\psi} = -\frac{g}{2} \hat{\mathbf{\sigma}}_z(z) \cdot \underline{\mathbf{v}}^*(z)$ denotes the vertical complex flux averaged over $\psi$ at the z-position. Indeed, $\mathbf{A}_+ = \mathbf{D}_+ = \mathbf{F}''_{in}$, $\mathbf{B}_+ = -\mathbf{C}_+ = \mathbf{F}'_{out}$, $\mathbf{A}_- = \mathbf{D}_- = \mathbf{F}''_{out}$ and $\mathbf{B}_- = -\mathbf{C}_- = \mathbf{F}'_{in}$. These quantities provide two types of contributions. For $\mathbf{F}_{in}$, the calculation of the average over the thickness, achieved through numerical integration, requires knowledge of the waveforms at every vertical position in the plate. This term is associated to the energy fluxes in the z-direction within the plate. In contrast, the coefficient $\mathbf{F}_{out}$ depends solely on the outgoing energy fluxes at the two boundaries, as the integral over the thickness can be analytically obtained. These coefficients are associated with the loss of energy resulting from leakage in the two external fluids.

Before discussing Eqs. (25), it is pertinent to examine the capabilities of the function $\hat{\mathbf{v}}(z)$ introduced in Eq. (1). This function does not affect the structure of the four averaged relations in its form, as its inclusion does not require any specific assumption for deriving this system. It represents the dependence of the velocity vector field on the z-direction, and its expression depends on the nature of the plate, i.e., whether it is isotropic, anisotropic, multilayered, or continuously variable. The dispersion curves and their associated fields are numerically obtained in [13] for all these cases, except for the continuously variable plate, for which additional numerical efforts are required to fully address this problem.
To interpret physically Eqs. (25a) and (25c), they can be rewritten as follows:

$$ 
\begin{aligned}
s' \cdot c_e &= 1 - D' = 1 - (\tau' F''_{in} + \tau'' F''_{out}) E_{tot}^{-1}, && (a) \\
-s'' \cdot c_e &= D'' = (\tau' F''_{out} + \tau'' F''_{in} + E_d) E_{tot}^{-1}, && (b)
\end{aligned}
\quad (26) $$

where $c_e$ denotes the energy velocity vector of guided modes in the plate plane. It is given by:

$$ c_e = \frac{P'}{E_{tot}} = c_e \mathbf{n}_e, \quad (27) $$

where the two vectors $P'$ and $\mathbf{n}_e$ denote the mean value of the power flow vector $p'$ and the unit vector in the Poynting direction, respectively. The newly introduced equations in Eq. (26) require further elucidation.

### Dispersion term $D'$, Eq. (26a)

First, let us analyze Eq. (26a). The term $D'$ can be interpreted as a contribution responsible for frequency dispersion, i.e., the phase and energy velocities depend on the frequency. This phenomenon arises from the interference of partial waves along the z-direction. To reinforce this interpretation, let us note the phase velocity projection on the $\mathbf{n}_e$-direction as $c'_e = (s' \cdot \mathbf{n}_e)^{-1}$. From Eq. (26a), it is a simple matter to obtain the following expression:

$$ c_e = (1 - D') c'_e. \quad (28) $$

For non-dispersive waves, such as bulk waves, the phase velocity in the direction of the Poynting vector is equal to the energy velocity, i.e., $D' = 0$. Otherwise, this term is never zero and guided waves are always dispersive.

### Damping term $D''$, Eq. (26b)

This term provides insightful information regarding the damping along the energy flow direction, i.e., along the direction of energy velocity. In accordance with the wave definition given by Eq. (2), the negative sign preceding Eq. (26b) indicates that the wave amplitude decreases in the energy direction for damping vectors $s''$ oriented oppositely. The component of the vector $-s''$ on the energy direction $\mathbf{n}_e$, noted $s''_e$, can be expressed as follows: $s''_e = -s'' \cdot \mathbf{n}_e$. Considering only the contribution of real frequencies, the attenuation factor along the energy flow direction due to this damping is given by $k''_e = \omega' s''_e$. The distance $\lambda_e$ covered during a time period at the energy velocity is given by $\lambda_e = 2\pi c_e (\omega')^{-1}$. This distance is not the wavelength in the $\mathbf{n}_e$-direction, since by virtue of Eq. (28) the energy velocity is not associated to the slowness vector in this direction. Finally, let us calculate the attenuation factor for the distance $\lambda_e$, which is given by: $k''_e \lambda_e = 2\pi D''$. The variable $D''$ clearly appears as an intrinsic quantity that contains information on different damping along the energy flow direction. Three factors contribute to the loss of amplitude, including viscoelastic effects, fluid leakages, and the coupling between time and space, inherent to the transient behavior.

For harmonic leaky waves, it must be positive. To gain a deeper insight into this assertion, let us analyze more precisely the vertical complex flux averaged over $\psi$, denoted by $\mathcal{P}_z(z_f)$ and calculated at the z-position of the two interfaces plate/fluid, i.e., located at $z_f = 0, \ell'$. At these interfaces, the normal stress and the normal displacement vector are continuous. Consequently, the vertical complex fluxes are also continuous at these two interfaces. Subsequently, the fluxes within the solid and the fluid adhere to the subsequent relation:

$$ \frac{\mathcal{P}_z(z_f)}{\ell} = \pm \frac{g}{2 \rho_f} |\hat{p}_f|^2 s_{zf}^* \quad (29) $$

where $\hat{p}_f, \rho_f$ and $s_{zf}$ are the variables associated with the inhomogeneous planes waves in the two fluids. They are, respectively, pressure amplitude terms, fluid densities and complex slownesses of the waves in the z-direction. In the context of leaky waves, the selection between the solutions $+$ or $-$ is contingent upon ensuring that the normal flux in the upper fluid is negative, whereas it is positive in the lower fluid. As a result, for harmonic wave, i.e., $\tau'' = 0, D''$ is positive since $F''_{out}$ and $E_d$ are both positive. The quantity $D''$ can be zero, as will be reported in Table 1. In such situation, the damping vector $s''$ is not necessary zero, it can be orthogonal to the energy flow direction. This provides an infinite number of possibilities because, to satisfy Eq. (26c), it is sufficient to add a component to the damping vector $s''_\perp$ that belongs to the plane orthogonal to the $\mathbf{n}_e$-direction. This additional component is entirely arbitrary, unlike the component along the energy flow direction.

This reinforces the fact that $D''$ is effectively an intrinsic quantity in a sense that this damping wave does not depend on the potential additional attenuation along $s''_\perp$. From that point of view, for the simplest case of harmonic bulk waves and Lamb waves in free single layer, some complementary discussions can be found in [14,15], respectively.

To conclude this part, Table 1 provides synthetic conditions under which the three variables $D', E_{p-k},^3$ and $D'',$ as defined in Eq. (26), can be zero. For each variable, the selected conditions must be regarded as both necessary and cumulative. The analysis encompasses scenarios involving plates in a vacuum as well as bulk waves. If one is under consideration, the other is, naturally, not applicable. If a particular condition is not necessary, it is referred to as being indifferent. For immersed plates, all the quantities $D', E_{p-k}$, and $D''$ are never zero, except for some very particular circumstances or for SH modes. This case is not documented on Table 1.

<footnote n="3">³ Deduced directly from Eq. (25b).</footnote>
<figure><img src="image_6.png" /><figcaption>Fig. 2. The phase velocity $c_φ$ and the attenuation $k''$ as functions of frequency $f = \omega/c_φ^2$, for solutions involving: both real slowness and frequency (—), real frequency (····), real slowness (— —) and a free plate (—). (a) Phase velocity $c_φ$ (b) Attenuation $k''$.</figcaption></figure>

with the fluid velocity straight line as the frequencies increase. This occurrence, analyzed in [17] for real frequencies, should merit additional investigations in the general case. However, it falls outside the scope of the present paper.

Before analyzing the energy and group velocity curves, let us comment on the definition of the average by Eq. (17). The remark concerns the inhomogeneity of waves in fluids, i.e., when $c_φ < c_F$. In this case, because the solution is real, the component of the damping vector lies along the z-axis. For physical reasons, in the specific case under consideration, it is oriented in the direction of negative z to tend to zero as z tends to infinity. The main consequence of this fact is that the energy flux of this single inhomogeneous plane wave is oriented along the interface. Consequently, this does not constitute a leakage. Nevertheless, this contribution must be included in the energy balance of the guided wave. To achieve this, the integration over the thickness defined in relation (17) must be extended to an integration over $[0, \infty)$, as illustrated in [18], which is denoted as:

$$ \mathbb{F}(x, y, t) = \langle \mathbb{F}(x, y, z, t) \rangle_{\infty} . \qquad (30) $$

This extension is feasible because the function to be integrated tends to zero at infinity. However, this assumption does not hold for the standard leaky wave, as in this case, the field diverges to infinity at infinity. This is not significant, as the energy of the guided wave is confined to the plate, rendering the integral extension unnecessary. This is further substantiated by Eq. (26b), where the contribution of the waves in the fluid is manifested as energy losses.

Taking these remarks into consideration, Fig. 3.b presents the energy and group velocities as functions of frequency for the physical solutions to our problem. The group velocity $c_g$ is defined from the phase variations with respect to frequency, variations which are directly connected to the group slowness $s_g$, such that: $s_g = \frac{\partial k'}{\partial \omega'}$, see e.g. [19], from which the group velocity is calculated. For both real frequency and slowness, it is well established that the group and energy velocities are identical [2]. To facilitate the interpretation of curves shown in Fig. 3.b, the following plotting conventions are employed: (i) real solutions are represented by solid lines, while complex solutions are represented by dashed lines; (ii) short-dashed and long-dashed lines indicate, respectively, complex slowness and complex frequency solutions, irrespective of color; (iii) when they can be differentiated, the red and blue colors represent energy and group velocities, respectively; (iv) when the two velocities cannot be discerned because they are equal for real solutions, or because they are numerically very close in the case of complex frequency solutions, they are plotted in purple, either dashed or solid in accordance with item (i); (v) the curves depicted with solid red lines and dot points represent the energy velocity without taking into account the waves in the fluids for real solutions.
<figure><img src="image_7.png" /><figcaption>Fig. 3. Group and energy velocities corresponding only to the physical solutions derived from the phase velocity of Fig. 2.a for a frequency range extended. (a) Phase velocity $c_φ$ for solutions involving: real slowness and frequency (—), real frequency (···) and real slowness (---). (b) Energy velocity (•••), group velocity (••••), and points where the group velocity equals the energy velocity (—, calculated using Eq. (30)) or approximately equals the energy velocity (—, calculated using Eq. (17)). Energy velocity, for real solutions, calculated using the average from Eq. (17) (---). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

Clearly when the solutions are complex, i.e., when $c_φ > c_F$, the group (blue dashed line) and energy (red dashed line) velocities differ drastically for complex slowness solutions. This is not surprising since, for inhomogeneous waves, these two velocities do not coincide, as demonstrated in [20]. This is not the case for complex frequency solutions, where, within the investigated phase velocity range, these two velocities, although not formally equal, are very close and difficult to distinguish (purple dashed line). On the contrary, for the *Scholte 1* and *Scholte 2* modes, the solutions are real, and both the group and energy velocities are rigorously identical as demonstrated for real solutions in [2]. To obtain this coincidence, it is important to note that the energy velocity must be defined according to Eq. (30) (solid purple line), rather than relying on the initial expression provided by Eq. (17) (solid red lines with dot points). The undesirable consequence is the emergence of discontinuity.

## 7. Conclusion

A set of four equations governing the instantaneous and local energy balance has been derived for a general inhomogeneous guided wave with a complex frequency propagating in an immersed multilayered plate. To derive the intrinsic energy balance equations, these equations were averaged over the phase term associated with the arbitrary amplitude phase of the guided waves. By proceeding in this way, the resulting set of equations depends explicitly on both time and space, enabling a clear expression of energy velocity for non-harmonic solutions of the dispersion equation, regardless of whether the guided waves exhibit complex wavenumbers and frequencies. The energy velocity, as obtained, has been compared with the group velocity. While these two velocities are nearly identical for solutions with complex frequencies, they can differ significantly for solutions with complex slownesses, as illustrated by numerical results.
CRediT authorship contribution statement

Marc Deschamps: Writing – original draft, Conceptualization. Eric Ducasse: Software, Methodology, Conceptualization.

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

No data was used for the research described in the article.

References

[1] M. Bruneau, C. Potel, Materials and Acoustics Handbook, John Wiley & Sons, 2013.

[2] B.A. Auld, Acoustic Fields and Waves in Solids, Vol. 2, Wiley-Interscience, 1973.

[3] O. Poncelet, M. Deschamps, Lamb waves generated by complex harmonic inhomogeneous plane waves, J. Acoust. Soc. Am. 102 (1997) 292–300.

[4] M. Hayes, Inhomogeneous plane waves, Arch. Ration. Mech. Anal. 85 (1984) 41–79.

[5] M. Hayes, Inhomogeneous electromagnetic plane waves in crystals, Arch. Ration. Mech. Anal. 97 (1987) 221–260.

[6] P. Boulanger, M. Hayes, Inhomogeneous plane waves in viscous fluids, Contin. Mech. Thermodyn. 2 (1990) 1–16.

[7] G. Caviglia, A. Morro, E. Pagani, Inhomogeneous waves in viscoelastic media, Wave Motion 12 (2) (1990) 143–159.

[8] M. Deschamps, P. Chevée, Reflection and refraction of a heterogeneous plane wave by a solid layer, Wave Motion 15 (1) (1992) 61–75.

[9] B. Hosten, M. Deschamps, B.R. Tittmann, Inhomogeneous wave generation and propagation in lossy anisotropic solids. Application to the viscoelastic characterization of composite materials, J. Acoust. Soc. Am. 82 (5) (1987) 1763–1770.

[10] D. Royer, E. Dieulesaint, Ondes élastiques dans les solides, Masson, Paris, 1996.

[11] Hwei P. Hsu, Theory and Problems of Probability, Random Variables, and Random Processes, McGraw-Hill, 1997.

[12] L.M. Brekhovskikh, O.A. Godin, Acoustic waves in absorbing anisotropic media, in: Acoustics of Layered Media I, Springer, 1990, pp. 158–164.

[13] E. Ducasse, M. Deschamps, Mode computation of immersed multilayer plates by solving an eigenvalue problem, Wave Motion 112 (2022) 102962.

[14] M. Deschamps, F. Assouline, Attenuation along the Poynting vector direction of inhomogeneous plane waves in absorbing and anisotropic solids, Acta Acust. United Acust. 86 (2000) 295–302.

[15] G. Neau, M.J.S. Lowe, M. Deschamps, Propagation of Lamb waves in anisotropic and absorbing plates: Theoretical derivation and experiments, AIP Conf. Proc. 615 (1) (2002) 1062–1069.

[16] A.L. Shuvalov, Theory of plane subsonic elastic waves in fluid-loaded anisotropic plates, Proc. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci. 458 (2028) (2002) 1323–1352.

[17] A.L. Shuvalov, O. Poncelet, M. Deschamps, Analysis of the dispersion spectrum of fluid-loaded anisotropic plates: flexural-type branches and real-valued loops, J. Sound Vib. 290 (3) (2006) 1175–1201.

[18] I.A. Nedospasov, V.G. Mozhaev, I.E. Kuznetsova, Unusual energy properties of leaky backward Lamb waves in a submerged plate, Ultrasonics 77 (2017) 95–99.

[19] M.J. Lighthill, Group velocity, Inst. Maths Applies 1 (1965) 1–28.

[20] M. Hayes, M.J.P. Musgrave, On energy flux and group velocity, Wave Motion 1 (1) (1979) 75–82.