# An inverse method for characterization of dynamic response of 2D structures under stochastic conditions

Xuefeng Li, Abdelmalek Zine, Mohamed Najib Ichchou, Noureddine Bouhaddi, Pascal Fossat

► To cite this version:

Xuefeng Li, Abdelmalek Zine, Mohamed Najib Ichchou, Noureddine Bouhaddi, Pascal Fossat. An inverse method for characterization of dynamic response of 2D structures under stochastic conditions. Chinese Journal of Aeronautics, 2024, 37 (3), pp.440 - 455. 10.1016/j.cja.2024.01.007. hal-04625670

HAL Id: hal-04625670

https://hal.science/hal-04625670

Submitted on 26 Jun 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
# An inverse method for characterization of dynamic response of 2D structures under stochastic conditions

Xuefeng Li, Abdelmalek Zine, Mohamed Najib Ichchou, Noureddine Bouhaddi, Pascal Fossat

## ► To cite this version:

Xuefeng Li, Abdelmalek Zine, Mohamed Najib Ichchou, Noureddine Bouhaddi, Pascal Fossat. An inverse method for characterization of dynamic response of 2D structures under stochastic conditions. Chinese Journal of Aeronautics, 2024, 37 (3), pp.440 - 455. 10.1016/j.cja.2024.01.007. hal-04625670

HAL Id: hal-04625670

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

Submitted on 26 Jun 2024

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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>

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

FULL LENGTH ARTICLE

# An inverse method for characterization of dynamic response of 2D structures under stochastic conditions

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

Xuefeng LI<sup>a</sup>, Abdelmalek ZINE<sup>b</sup>, Mohamed ICHCHOU<sup>a,*</sup>, Noureddine BOUHADDI<sup>c</sup>, Pascal FOSSAT<sup>a</sup>

<sup>a</sup> Laboratory of Tribology and Systems Dynamic, Central School of Lyon, Ecully 69134, France

<sup>b</sup> Institut Camille Jordan, Central School of Lyon, Ecully 69134, France

<sup>c</sup> Department of Applied Mechanics, University of Burgundy Franche-Comté, Besançon 25000, France

Received 17 March 2023; revised 8 June 2023; accepted 2 July 2023

Available online 10 January 2024

KEYWORDS

Inverse method;
Dispersion relation;
Wavenumber space;
Periodic plates;
Stochastic conditions;
Wave propagation
characterization

**Abstract** The reliable estimation of the wavenumber space (k-space) of the plates remains a long-term concern for acoustic modeling and structural dynamic behavior characterization. Most current analyses of wavenumber identification methods are based on the deterministic hypothesis. To this end, an inverse method is proposed for identifying wave propagation characteristics of two-dimensional structures under stochastic conditions, such as wavenumber space, dispersion curves, and band gaps. The proposed method is developed based on an algebraic identification scheme in the polar coordinate system framework, thus named Algebraic K-Space Identification (AKSI) technique. Additionally, a model order estimation strategy and a wavenumber filter are proposed to ensure that AKSI is successfully applied. The main benefit of AKSI is that it is a reliable and fast method under four stochastic conditions: (A) High level of signal noise; (B) Small perturbation caused by uncertainties in measurement points' coordinates; (C) Non-periodic sampling; (D) Unknown structural periodicity. To validate the proposed method, we numerically benchmark AKSI and three other inverse methods to extract dispersion curves on three plates under stochastic conditions. One experiment is then performed on an isotropic steel plate. These investigations demonstrate that AKSI is a good in-situ k-space estimator under stochastic conditions.

© 2024 Production and hosting by Elsevier Ltd. on behalf of Chinese Society of Aeronautics and Astronautics. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).

E-mail address: mohamed.ichchou@ec-lyon.fr (M. ICHCHOU).
Peer review under responsibility of Editorial Committee of CJA.

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

Production and hosting by Elsevier

1. Introduction

In the field of the transportation industry, characterizing the vibroacoustic dynamic behavior of structures is essential for structural optimization,¹⁻⁴ noise control,⁵,⁶ and damage detection.⁷⁻¹¹ In the last few decades, inverse methods on wavenum-

<footnote>* Corresponding author.</footnote>
<figure><img src="image_2.png" /></figure>

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

FULL LENGTH ARTICLE

# An inverse method for characterization of dynamic response of 2D structures under stochastic conditions

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

Xuefeng LI<sup>a</sup>, Abdelmalek ZINE<sup>b</sup>, Mohamed ICHCHOU<sup>a,*</sup>, Noureddine BOUHADDI<sup>c</sup>, Pascal FOSSAT<sup>a</sup>

<sup>a</sup> Laboratory of Tribology and Systems Dynamic, Central School of Lyon, Ecully 69134, France

<sup>b</sup> Institut Camille Jordan, Central School of Lyon, Ecully 69134, France

<sup>c</sup> Department of Applied Mechanics, University of Burgundy Franche-Comté, Besançon 25000, France

Received 17 March 2023; revised 8 June 2023; accepted 2 July 2023

Available online 10 January 2024

KEYWORDS

Inverse method;
Dispersion relation;
Wavenumber space;
Periodic plates;
Stochastic conditions;
Wave propagation
characterization

**Abstract** The reliable estimation of the wavenumber space (k-space) of the plates remains a long-term concern for acoustic modeling and structural dynamic behavior characterization. Most current analyses of wavenumber identification methods are based on the deterministic hypothesis. To this end, an inverse method is proposed for identifying wave propagation characteristics of two-dimensional structures under stochastic conditions, such as wavenumber space, dispersion curves, and band gaps. The proposed method is developed based on an algebraic identification scheme in the polar coordinate system framework, thus named Algebraic K-Space Identification (AKSI) technique. Additionally, a model order estimation strategy and a wavenumber filter are proposed to ensure that AKSI is successfully applied. The main benefit of AKSI is that it is a reliable and fast method under four stochastic conditions: (A) High level of signal noise; (B) Small perturbation caused by uncertainties in measurement points' coordinates; (C) Non-periodic sampling; (D) Unknown structural periodicity. To validate the proposed method, we numerically benchmark AKSI and three other inverse methods to extract dispersion curves on three plates under stochastic conditions. One experiment is then performed on an isotropic steel plate. These investigations demonstrate that AKSI is a good in-situ k-space estimator under stochastic conditions.

© 2024 Production and hosting by Elsevier Ltd. on behalf of Chinese Society of Aeronautics and Astronautics. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).

E-mail address: mohamed.ichchou@ec-lyon.fr (M. ICHCHOU).
Peer review under responsibility of Editorial Committee of CJA.

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

Production and hosting by Elsevier

1. Introduction

In the field of the transportation industry, characterizing the vibroacoustic dynamic behavior of structures is essential for structural optimization,¹⁻⁴ noise control,⁵,⁶ and damage detection.⁷⁻¹¹ In the last few decades, inverse methods on wavenum-

<footnote>* Corresponding author.</footnote>
ber identification of 2D structures have attracted increasing attention due to their ability to achieve characterization of structural dynamic behavior, such as dispersion relation, wavenumber space (k-space), band gaps, and damping loss factors. By solving an inverse problem, the extracted wavenumber can be used to determine the mechanical parameters of structures.¹²⁻¹⁴ Numerous 2D inverse methods are available in the literature. These inverse methods use only the responses of the structure as input parameters to extract wavenumber, leading to their wide range of applications in industrial practice. In general, these methods can be divided into non-linear inverse methods and linear Prony family inverse methods.

In the non-linear family inverse methods. The Inhomogeneous Wave Correlation (IWC) method is the most popular inverse method. This method was proposed by Berthaut et al.,¹⁵ which aims to find complex wavenumbers that maximize the correlation between an inhomogeneous wave mode and measured wave field at a fixed wave propagation angle θ. At a given frequency, the k-space can be formed by repeating the IWC procedure at each angle θ. In the past years, this method has been applied to extract k-space on a series of complex structures, such as locally resonant metamaterial,¹⁶ honeycomb sandwich plate,¹² and curved structures.⁵ These works show that IWC is a reliable inverse method in the medium and high-frequency domain, even in the presence of measurement error.¹⁷,¹⁸ Furthermore, this method is not limited by periodic sampling. However, as the IWC is close to the Fourier transform, the extracted wavenumber is less accurate in the low-frequency range. In addition, this method involves a non-linear iterative process, leading to a high computational cost. A very similar method, called Spatial Laplace Transform for Complex Wavenumber recovery (SLaTCoW) method, has been proposed in Ref. 19 which differs from IWC in that it applies the Laplace transform to the signal. This method has been used to extract dispersion curves for a polycrystalline aluminum plate.²⁰ In addition to the same limitations as IWC, another difficulty with this method is the choice of a suitable cost function.

In the Prony family inverse method, the classical Prony method is well known for high-resolution parameter identification in noiseless conditions. In the context of vibroacoustic, the classical Prony method has been applied to extract dispersion curves of a constrained-layer damping sandwich plate²¹ and the k-space of a wooden plate.²² Margerit et al.²³ proposed the High-Resolution Wavevector Analysis (HRWA) method, allowing the automatic k-space identification without the need to the selection of angles in prior. However, its practical implementation is challenging due to several factors, including the estimation of the total number of waves propagating in 2D structures and the use of 2D Finite Impulse Response (FIR) Filter. These factors limit the ability of HRWA to extract a smooth and full k-space even in perfect condition. Recently, Boukadia et al.²⁴ proposed INverse CONvolution MEthod

(INCOME) to overcome this problem. This method draws the inspiration of the Wave Finite Element Method (WFEM) to rebuild the formula of the Prony method in the convolution framework, allowing us to extract a full and accurate k-space by a single 2D convolution kernel in the noiseless condition. Unfortunately, this method is limited to periodic sampling and sensitive to signal noise. In terms of the application of the Prony method on periodic structures, Ribeiro et al.²⁵ proposed a Bloch Wavenumber Identification method, called “BWI”, to extract the Bloch wavenumber of periodic structures. An important conclusion can be obtained from BWI: the Prony method requires that the periodicity of the sampling must be an integer multiple of the structural periodicity. In practice, however, the structural periodicity is often not accurately known due to structural deformation or unknown when periodic structures are surrounded by other materials, resulting in the failure of wavenumber extraction using linear Prony family inverse methods.

In Ref. 26, the authors proposed the Algebraic Wavenumber Identification (AWI) method to extract wavenumbers of 1D structures and this method has been validated to have a good performance under different stochastic conditions.²⁷ However, it is worth noting that AWI is designed for 1D structures and it is not able to extract the k-space of 2D structures alone. In addition, the methodology of AWI does not take into account direction parameters in 2D structures. Based on this context, the present work aims to extend AWI to extract the k-space of 2D structures, overcoming the main limitations of other 2D inverse methods mentioned in Table 1. To differentiate AWI, the method proposed is named Algebraic K-Space Identification (AKSI) technique.

Similar to AWI, the methodology of AKSI is also developed in the algebraic parameter estimator framework,²⁸ involving the transformation of signals: spatial domain → wavenumber domain → spatial domain. However, the parameter estimation of multidimensional systems (such as 2D structures) by the algebraic parameter estimator needs to solve the inverse Laplace transforms for complex higher-order partial differential equations, which is more problematic than ordinary differential equations in unidimensional systems (such as 1D structures). To address this problem, we established AKSI with the direction parameter in the polar coordinate system based on the conversion of the wavenumber and the geometric coordinates of the measurement points between the Cartesian and polar coordinate systems. This results in that the essence of the proposed method lies in estimating the unknown parameters of an ordinary differential equation, avoiding the introduction of partial differential equations, which allows us to extract the wavenumber in each direction of 2D structures under stochastic conditions. In addition, we proposed a model order estimation and a wavenumber filtering strategy adapted to AKSI to guarantee the successful application of general properties between popular 2D inverse methods for dispersion curves extraction of 2D structures.

<table><caption>**Table 1** Comparison of general properties between popular 2D inverse methods for dispersion curves extraction of 2D structures.</caption><thead><tr><td>Method</td><td>Solver</td><td>Sampling limit</td><td>Robust to perturbation</td><td>Unknown periodicity</td><td>K-space</td></tr></thead><tbody><tr><td>IWC</td><td>Nonlinear</td><td>None</td><td>Yes</td><td>Yes</td><td>Yes</td></tr><tr><td>SLaTCoW</td><td>Nonlinear</td><td>None</td><td>Yes</td><td>Yes</td><td>Yes</td></tr><tr><td>HRWA</td><td>Linear</td><td>Periodicity</td><td>No</td><td>No</td><td>Yes</td></tr><tr><td>INCOME</td><td>Linear</td><td>Periodicity</td><td>No</td><td>No</td><td>Yes</td></tr><tr><td>BWI</td><td>Linear</td><td>Periodicity</td><td>Yes</td><td>No</td><td>No</td></tr></tbody></table>
ber identification of 2D structures have attracted increasing attention due to their ability to achieve characterization of structural dynamic behavior, such as dispersion relation, wavenumber space (k-space), band gaps, and damping loss factors. By solving an inverse problem, the extracted wavenumber can be used to determine the mechanical parameters of structures.¹²⁻¹⁴ Numerous 2D inverse methods are available in the literature. These inverse methods use only the responses of the structure as input parameters to extract wavenumber, leading to their wide range of applications in industrial practice. In general, these methods can be divided into non-linear inverse methods and linear Prony family inverse methods.

In the non-linear family inverse methods. The Inhomogeneous Wave Correlation (IWC) method is the most popular inverse method. This method was proposed by Berthaut et al.,¹⁵ which aims to find complex wavenumbers that maximize the correlation between an inhomogeneous wave mode and measured wave field at a fixed wave propagation angle θ. At a given frequency, the k-space can be formed by repeating the IWC procedure at each angle θ. In the past years, this method has been applied to extract k-space on a series of complex structures, such as locally resonant metamaterial,¹⁶ honeycomb sandwich plate,¹² and curved structures.⁵ These works show that IWC is a reliable inverse method in the medium and high-frequency domain, even in the presence of measurement error.¹⁷,¹⁸ Furthermore, this method is not limited by periodic sampling. However, as the IWC is close to the Fourier transform, the extracted wavenumber is less accurate in the low-frequency range. In addition, this method involves a non-linear iterative process, leading to a high computational cost. A very similar method, called Spatial Laplace Transform for Complex Wavenumber recovery (SLaTCoW) method, has been proposed in Ref. 19 which differs from IWC in that it applies the Laplace transform to the signal. This method has been used to extract dispersion curves for a polycrystalline aluminum plate.²⁰ In addition to the same limitations as IWC, another difficulty with this method is the choice of a suitable cost function.

In the Prony family inverse method, the classical Prony method is well known for high-resolution parameter identification in noiseless conditions. In the context of vibroacoustic, the classical Prony method has been applied to extract dispersion curves of a constrained-layer damping sandwich plate²¹ and the k-space of a wooden plate.²² Margerit et al.²³ proposed the High-Resolution Wavevector Analysis (HRWA) method, allowing the automatic k-space identification without the need to the selection of angles in prior. However, its practical implementation is challenging due to several factors, including the estimation of the total number of waves propagating in 2D structures and the use of 2D Finite Impulse Response (FIR) Filter. These factors limit the ability of HRWA to extract a smooth and full k-space even in perfect condition. Recently, Boukadia et al.²⁴ proposed INverse CONvolution MEthod

(INCOME) to overcome this problem. This method draws the inspiration of the Wave Finite Element Method (WFEM) to rebuild the formula of the Prony method in the convolution framework, allowing us to extract a full and accurate k-space by a single 2D convolution kernel in the noiseless condition. Unfortunately, this method is limited to periodic sampling and sensitive to signal noise. In terms of the application of the Prony method on periodic structures, Ribeiro et al.²⁵ proposed a Bloch Wavenumber Identification method, called “BWI”, to extract the Bloch wavenumber of periodic structures. An important conclusion can be obtained from BWI: the Prony method requires that the periodicity of the sampling must be an integer multiple of the structural periodicity. In practice, however, the structural periodicity is often not accurately known due to structural deformation or unknown when periodic structures are surrounded by other materials, resulting in the failure of wavenumber extraction using linear Prony family inverse methods.

In Ref. 26, the authors proposed the Algebraic Wavenumber Identification (AWI) method to extract wavenumbers of 1D structures and this method has been validated to have a good performance under different stochastic conditions.²⁷ However, it is worth noting that AWI is designed for 1D structures and it is not able to extract the k-space of 2D structures alone. In addition, the methodology of AWI does not take into account direction parameters in 2D structures. Based on this context, the present work aims to extend AWI to extract the k-space of 2D structures, overcoming the main limitations of other 2D inverse methods mentioned in Table 1. To differentiate AWI, the method proposed is named Algebraic K-Space Identification (AKSI) technique.

Similar to AWI, the methodology of AKSI is also developed in the algebraic parameter estimator framework,²⁸ involving the transformation of signals: spatial domain → wavenumber domain → spatial domain. However, the parameter estimation of multidimensional systems (such as 2D structures) by the algebraic parameter estimator needs to solve the inverse Laplace transforms for complex higher-order partial differential equations, which is more problematic than ordinary differential equations in unidimensional systems (such as 1D structures). To address this problem, we established AKSI with the direction parameter in the polar coordinate system based on the conversion of the wavenumber and the geometric coordinates of the measurement points between the Cartesian and polar coordinate systems. This results in that the essence of the proposed method lies in estimating the unknown parameters of an ordinary differential equation, avoiding the introduction of partial differential equations, which allows us to extract the wavenumber in each direction of 2D structures under stochastic conditions. In addition, we proposed a model order estimation and a wavenumber filtering strategy adapted to AKSI to guarantee the successful application.

<table><caption>**Table 1** Comparison of general properties between popular 2D inverse methods for dispersion curves extraction of 2D structures.</caption><thead><tr><td>Method</td><td>Solver</td><td>Sampling limit</td><td>Robust to perturbation</td><td>Unknown periodicity</td><td>K-space</td></tr></thead><tbody><tr><td>IWC</td><td>Nonlinear</td><td>None</td><td>Yes</td><td>Yes</td><td>Yes</td></tr><tr><td>SLaTCoW</td><td>Nonlinear</td><td>None</td><td>Yes</td><td>Yes</td><td>Yes</td></tr><tr><td>HRWA</td><td>Linear</td><td>Periodicity</td><td>No</td><td>No</td><td>Yes</td></tr><tr><td>INCOME</td><td>Linear</td><td>Periodicity</td><td>No</td><td>No</td><td>Yes</td></tr><tr><td>BWI</td><td>Linear</td><td>Periodicity</td><td>Yes</td><td>No</td><td>No</td></tr></tbody></table>
tion of AKSI in the k-space extraction of 2D structures under stochastic conditions. The main contributions are the following:

* Extend AWI to extract the k-space of 2D structures under stochastic conditions. Moreover, the proposed strategy for model order estimation and wavenumber filtering is necessary for the proposed method to be implemented successfully.

* As the AKSI inherits the advantages of an algebraic estimator, it exhibits good robustness under perturbation conditions. In addition, AKSI only needs to solve several linear algebraic equations, which can reduce computational costs compared to nonlinear methods such as Inhomogeneous Wave Correlation (IWC).

* Consider a more realistic signal model that takes into account external signal noise and small perturbations. In addition, the performance of the proposed method is investigated when the structural periodicity of periodic structures is unknown.

* Compare the proposed method with the other three popular 2D inverse methods numerically and experimentally, validating that the proposed method can overcome the main limitations of other popular inverse methods listed in Table 1. The investigations provide a benchmark for developing inverse methods in the future.

First, briefly review the tested inverse methods, including IWC, BWI, and INCOME. Secondly, a stochastic signal model of the plate is presented, and then the formula of AKSI is given. Thirdly, benchmark AKSI and other inverse methods in three numerical studies. Then, an experimental result of identifying the dispersion curve and k-space using AKSI and IWC for a thin isotropic plate is presented. Finally, conclusions and open issues to be addressed in the future are provided.

## 2. Review of tested inverse methods

### 2.1. Inhomogeneous Wave Correlation (IWC) method

The principle of IWC is similar to the Modal Assurance Criterion (MAC).²⁹ At a fixed frequency, assuming that the displacements $U(x_j, y_j)$ of the N measurement points are known at a fixed wave propagation angle $\theta$, the complex wavenumber can be obtained by maximizing the IWC function:

$$ \text{IWC}(k, \theta) = \frac{\left| \sum_{j=1}^{N} U(x_j, y_j) \hat{\sigma}_{k, \theta}^{*}(x_j, y_j) \Omega_j \right|}{\sqrt{\sum_{j=1}^{N} |U(x_j, y_j)|^2 \Omega_j \sum_{i=1}^{N} |\hat{\sigma}_{k, \theta}(x_j, y_j)|^2 \Omega_j}} \quad (1) $$

where inhomogeneous wave model $\hat{\sigma}_{k,\theta}^*(x_j, y_j)$ is defined as $e^{-jk(\theta)(x_j \cos \theta + y_j \sin \theta)}$; the superscript * denotes the complex conjugate; $\Omega_j$ is an estimation of surface around point at $(x_j, y_j)$; k is the complex wavenumber.

### 2.2. Bloch Wavenumber Identification (BWI)

The Prony method is extended to extract the Bloch wavenumber of periodic structures.²⁵ The core of this method is that the measurement points first need to be grouped. The measure-

ment points in each group need to satisfy two requirements: (A) The measurement points have the same sampling interval; (B) The sampling interval is an integer multiple of the structural periodicity. Then the wavenumber in x and y directions can be estimated along the line of the Prony method. When only one wave type propagates in an infinite plate, two groups of measured data $U$ and $U'$ need to be first arranged, then the coefficients of the characteristic polynomial $a\lambda_x^2 + b\lambda_x + c = 0$ can be estimated using the least-squares method. The whole process can be expressed in matrix format:

$$ \begin{bmatrix} U_1 & U_2 & U_3 \\ U_2 & U_3 & U_4 \\ \vdots & \vdots & \vdots \\ U_{N-2} & U_{N-1} & U_N \end{bmatrix} \begin{bmatrix} a \\ b \\ c \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ 0 \end{bmatrix} \quad (2) $$

Based on the obtained coefficients vector $[a \ b \ c]^T$, propagation constants $\lambda_x$ can be calculated by the characteristic polynomial. Finally, one can obtain the wavenumber by $\lambda_x = e^{-jk_x \Delta x}$ with sampling interval $\Delta x$. The wavenumber $k_y$ can be done in the same way.

### 2.3. Inverse Convolution Method (INCOME)

The main contribution of INCOME is that all properties of the k-space of 2D structures can be estimated by only a 2D convolution kernel. Moreover, the k-space can be extracted with high numerical precision in perfect conditions. To simplify the explanation of INCOME, only one wave type is assumed to propagate in the 2D structure. The wavenumber can be estimated in three steps. The first step is to design the 2D convolutional kernel as

$$ \begin{bmatrix} S_{1,3} & S_{2,3} & S_{3,3} \\ S_{1,2} & S_{2,2} & S_{3,2} \\ S_{1,1} & S_{2,1} & S_{3,1} \end{bmatrix} = \begin{bmatrix} e & h & f \\ g & d & g \\ f & h & e \end{bmatrix} \quad (3) $$

The second step is to estimate the elements of this 2D convolutional kernel by calculating the convolution of displacements $U(x_j, y_j)_{(i,j) \in [1,N] \times [1,M]}$ with the 2D convolution kernel $S$ as

$$ V = U * S = 0 \quad (4) $$

where “*” is the convolution symbolic.

The explicit form of Eq. (4) can be expressed as

$$ \begin{cases} \forall(p, q) \in [1, N-3] \times [1, M-3] \\ V_{p,q} = \sum_{(i,j) \in [1,3]^2} S_{ij} U_{p+3-i,q+3-j} = 0 \end{cases} \quad (5) $$

Once elements of the 2D convolution kernel are calculated, the propagation constants can be obtained by

$$ \begin{aligned} & a + b \left( \hat{\lambda}_x + \frac{1}{\hat{\lambda}_x} \right) + c \left( \hat{\lambda}_y + \frac{1}{\hat{\lambda}_y} \right) + d \left( \frac{\hat{\lambda}_x}{\hat{\lambda}_y} + \frac{\hat{\lambda}_y}{\hat{\lambda}_x} \right) \\ & + e \left( \hat{\lambda}_x \hat{\lambda}_y + \frac{1}{\hat{\lambda}_x \hat{\lambda}_y} \right) = 0 \end{aligned} \quad (6) $$
tion of AKSI in the k-space extraction of 2D structures under stochastic conditions. The main contributions are the following:

* Extend AWI to extract the k-space of 2D structures under stochastic conditions. Moreover, the proposed strategy for model order estimation and wavenumber filtering is necessary for the proposed method to be implemented successfully.

* As the AKSI inherits the advantages of an algebraic estimator, it exhibits good robustness under perturbation conditions. In addition, AKSI only needs to solve several linear algebraic equations, which can reduce computational costs compared to nonlinear methods such as Inhomogeneous Wave Correlation (IWC).

* Consider a more realistic signal model that takes into account external signal noise and small perturbations. In addition, the performance of the proposed method is investigated when the structural periodicity of periodic structures is unknown.

* Compare the proposed method with the other three popular 2D inverse methods numerically and experimentally, validating that the proposed method can overcome the main limitations of other popular inverse methods listed in Table 1. The investigations provide a benchmark for developing inverse methods in the future.

First, briefly review the tested inverse methods, including IWC, BWI, and INCOME. Secondly, a stochastic signal model of the plate is presented, and then the formula of AKSI is given. Thirdly, benchmark AKSI and other inverse methods in three numerical studies. Then, an experimental result of identifying the dispersion curve and k-space using AKSI and IWC for a thin isotropic plate is presented. Finally, conclusions and open issues to be addressed in the future are provided.

## 2. Review of tested inverse methods

### 2.1. Inhomogeneous Wave Correlation (IWC) method

The principle of IWC is similar to the Modal Assurance Criterion (MAC).²⁹ At a fixed frequency, assuming that the displacements $U(x_j, y_j)$ of the N measurement points are known at a fixed wave propagation angle $\theta$, the complex wavenumber can be obtained by maximizing the IWC function:

$$ \text{IWC}(k, \theta) = \frac{\left| \sum_{j=1}^{N} U(x_j, y_j) \hat{\sigma}_{k, \theta}^{*}(x_j, y_j) \Omega_j \right|}{\sqrt{\sum_{j=1}^{N} |U(x_j, y_j)|^2 \Omega_j \sum_{i=1}^{N} |\hat{\sigma}_{k, \theta}(x_j, y_j)|^2 \Omega_j}} \quad (1) $$

where inhomogeneous wave model $\hat{\sigma}_{k,\theta}^*(x_j, y_j)$ is defined as $e^{-jk(\theta)(x_j \cos \theta + y_j \sin \theta)}$; the superscript * denotes the complex conjugate; $\Omega_j$ is an estimation of surface around point at $(x_j, y_j)$; k is the complex wavenumber.

### 2.2. Bloch Wavenumber Identification (BWI)

The Prony method is extended to extract the Bloch wavenumber of periodic structures.²⁵ The core of this method is that the measurement points first need to be grouped. The measure-

ment points in each group need to satisfy two requirements: (A) The measurement points have the same sampling interval; (B) The sampling interval is an integer multiple of the structural periodicity. Then the wavenumber in x and y directions can be estimated along the line of the Prony method. When only one wave type propagates in an infinite plate, two groups of measured data $U$ and $U'$ need to be first arranged, then the coefficients of the characteristic polynomial $a\lambda_x^2 + b\lambda_x + c = 0$ can be estimated using the least-squares method. The whole process can be expressed in matrix format:

$$ \begin{bmatrix} \begin{bmatrix} U_1 & U_2 & U_3 \\ U_2 & U_3 & U_4 \\ \vdots & \vdots & \vdots \\ U_{N-2} & U_{N-1} & U_N \end{bmatrix} & \begin{bmatrix} a \\ b \\ c \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ 0 \end{bmatrix} \\ \begin{bmatrix} U'_1 & U'_2 & U'_3 \\ U'_2 & U'_3 & U'_4 \\ \vdots & \vdots & \vdots \\ U'_{N-2} & U'_{N-1} & U'_N \end{bmatrix} & \begin{bmatrix} e & h & f \\ g & d & g \\ f & h & e \end{bmatrix} \end{bmatrix} \quad (2) $$

Based on the obtained coefficients vector $[a \ b \ c]^T$, propagation constants $\lambda_x$ can be calculated by the characteristic polynomial. Finally, one can obtain the wavenumber by $\lambda_x = e^{-jk_x \Delta x}$ with sampling interval $\Delta x$. The wavenumber $k_y$ can be done in the same way.

### 2.3. Inverse Convolution Method (INCOME)

The main contribution of INCOME is that all properties of the k-space of 2D structures can be estimated by only a 2D convolution kernel. Moreover, the k-space can be extracted with high numerical precision in perfect conditions. To simplify the explanation of INCOME, only one wave type is assumed to propagate in the 2D structure. The wavenumber can be estimated in three steps. The first step is to design the 2D convolutional kernel as

$$ \begin{bmatrix} S_{1,3} & S_{2,3} & S_{3,3} \\ S_{1,2} & S_{2,2} & S_{3,2} \\ S_{1,1} & S_{2,1} & S_{3,1} \end{bmatrix} = \begin{bmatrix} e & h & f \\ g & d & g \\ f & h & e \end{bmatrix} \quad (3) $$

The second step is to estimate the elements of this 2D convolutional kernel by calculating the convolution of displacements $U(x_j, y_j)_{(i,j) \in [1,N] \times [1,M]}$ with the 2D convolution kernel $S$ as

$$ \mathbf{V} = \mathbf{U} * \mathbf{S} = \mathbf{0} \quad (4) $$

where “*” is the convolution symbolic.

The explicit form of Eq. (4) can be expressed as

$$ \begin{cases} \forall(p, q) \in [1, N-3] \times [1, M-3] \\ V_{p,q} = \sum_{(i,j) \in [1,3]^2} S_{ij} U_{p+3-i,q+3-j} = 0 \end{cases} \quad (5) $$

Once elements of the 2D convolution kernel are calculated, the propagation constants can be obtained by

$$ \begin{aligned} & a + b \left( \lambda_x + \frac{1}{\lambda_x} \right) + c \left( \lambda_y + \frac{1}{\lambda_y} \right) + d \left( \frac{\lambda_x}{\lambda_y} + \frac{\lambda_y}{\lambda_x} \right) \\ & + e \left( \lambda_x \lambda_y + \frac{1}{\lambda_x \lambda_y} \right) = 0 \end{aligned} \quad (6) $$
The final step is to calculate the wavenumber using the relationship:

$$ \left\{ \begin{aligned} \lambda_x &= e^{-jk_x \Delta x} = e^{-jk \cos \theta \Delta x} \\ \lambda_y &= e^{-jk_y \Delta y} = e^{-jk \sin \theta \Delta y} \end{aligned} \right. \qquad (7) $$

where $\Delta x$ and $\Delta y$ are the sampling interval along x and y directions. The essence of INCOME is the Prony method because the elements of the 2D convolution kernel are calculated along the way of the Prony method. The difference is that the Prony method can only extract the wavenumber in one direction using the characteristic polynomial with a single variable, as explained in Section 2.2. In contrast, INCOME proposes a characteristic polynomial containing two variables $\lambda_x$ and $\lambda_y$ in the WFEM framework, allowing a full k-space to be extracted by changing the angle $\theta$.

### 3. Algebraic K-Space Identification (AKSI)

Concerned with introducing the theory of AKSI in polar coordinate systems. First, Section 3.1 establishes a stochastic signal model considering the influence of signal noise and the small perturbation problem. Then, Section 3.2 introduces the formulation of AKSI in polar coordinates, the essence of which lies in the identification of the unknown parameters of an Ordinary Differential Equation (ODE). Next, Section 3.3 describes two essential signal processing tools for implementing the AKSI algorithm. The first aims to estimate model order and the second concerns the wavenumbers filtering. Finally, Section 3.4 presents the implementation procedure of AKSI.

#### 3.1. Stochastic signal model

Under perfect conditions, the harmonic displacement at any measured point $(x, y)$ of a plate can be modeled by $n_w$ plane waves as

$$ U(x, y) = \sum_{m=1}^{n_w} A_m e^{ik_{x,m}x} e^{ik_{y,m}y} \qquad (8) $$

where $A_m$ is the complex amplitude of m th wave; $k_{x,m}$ and $k_{y,m}$ are the components of the m th wavenumber $k_m$ in x and y directions, respectively. Fig. 1(a) illustrates the displacement field of the isotropic aluminum plate, studied in Section 4.1, at 150 Hz. As shown in Fig. 1(a), when the wave propagation angle $\theta$ is given, the displacement $U(x_n, y_n)$ of any Point P in this direction can also be expressed as

$$ U(\theta, x_n, y_n) = \sum_{m=1}^{n_w} A_m e^{ik_{\theta,x_n,m}x_n} e^{ik_{\theta,y_n,m}y_n} \qquad (9) $$

In practice, the structural responses will inevitably be affected by various disturbances. To create a more realistic signal model, two disturbance factors are considered: (A) Small perturbation due to the geometrical coordinates' uncertainty; (B) External signal noise. Small perturbations are common in practice and can lead to mesh distortion, as shown in Fig. 1(b). The direct consequence is that the geometric coordinates of the measured points do not match the actual structural response. For example, under the influence of the small perturbation, Point P moves to Point P' with an unknown and random coordinate, leading to a mismatch between the

<figure><img src="image_6.png" /><figcaption>(a) Transformation of displacement and coordinates of measurement point P in Cartesian and polar coordinate systems</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(b) Comparison of a perfect 2D grid and a 2D distortion grid due to small perturbation problem</figcaption></figure>

<caption>**Fig. 1** Coordinate system transformation and 2D grid distortion.</caption>

actual response of the structure (the response of Point P') and the coordinates of the measured point (coordinate of Point P). Furthermore, in contrast to the 1D case, the small perturbation of the 2D measurement points is more complex, because the coordinates of the measurement points can be perturbed randomly in any direction. To solve this problem, in a probabilistic framework, the small perturbation problem can be quantified by introducing two random variables into the coordinate of each point.

$$ \left\{ \begin{aligned} \hat{x}_n &= x_n \pm \xi_n \Delta x_n \\ \hat{y}_n &= y_n \pm \zeta_n \Delta y_n \end{aligned} \right. \qquad n = 1, 2, \dots, N \qquad (10) $$

where random variables $\xi_n$ and $\zeta_n$ are defined as the small perturbation ratio which can obey different distributions, such as uniform distribution and gradient distribution. Therefore, the signal model in the presence of small perturbation can be expressed as

$$ U(\theta, \hat{x}_n, \hat{y}_n) = \sum_{m=1}^{n_w} A_m e^{ik_{\theta,x_n,m}(x_n \pm \xi_n \Delta x_n)} e^{ik_{\theta,y_n,m}(y_n \pm \zeta_n \Delta y_n)} \qquad (11) $$

where the small perturbation part can be expanded by a first-order Taylor series for illustration of measurement error caused by the small perturbation problem:
The final step is to calculate the wavenumber using the relationship:

$$ \begin{cases} \lambda_x = e^{-jk_x \Delta x} = e^{-jk \cos \theta \Delta x} \\ \lambda_y = e^{-jk_y \Delta y} = e^{-jk \sin \theta \Delta y} \end{cases} \qquad (7) $$

where $\Delta x$ and $\Delta y$ are the sampling interval along x and y directions. The essence of INCOME is the Prony method because the elements of the 2D convolution kernel are calculated along the way of the Prony method. The difference is that the Prony method can only extract the wavenumber in one direction using the characteristic polynomial with a single variable, as explained in Section 2.2. In contrast, INCOME proposes a characteristic polynomial containing two variables $\lambda_x$ and $\lambda_y$ in the WFEM framework, allowing a full k-space to be extracted by changing the angle $\theta$.

### 3. Algebraic K-Space Identification (AKSI)

Concerned with introducing the theory of AKSI in polar coordinate systems. First, Section 3.1 establishes a stochastic signal model considering the influence of signal noise and the small perturbation problem. Then, Section 3.2 introduces the formulation of AKSI in polar coordinates, the essence of which lies in the identification of the unknown parameters of an Ordinary Differential Equation (ODE). Next, Section 3.3 describes two essential signal processing tools for implementing the AKSI algorithm. The first aims to estimate model order and the second concerns the wavenumbers filtering. Finally, Section 3.4 presents the implementation procedure of AKSI.

#### 3.1. Stochastic signal model

Under perfect conditions, the harmonic displacement at any measured point $(x, y)$ of a plate can be modeled by $n_w$ plane waves as

$$ U(x, y) = \sum_{m=1}^{n_w} A_m e^{ik_{0,m}x} e^{ik_{y,m}y} \qquad (8) $$

where $A_m$ is the complex amplitude of m th wave; $k_{0,m}$ and $k_{y,m}$ are the components of the m th wavenumber $k_m$ in x and y directions, respectively. Fig. 1(a) illustrates the displacement field of the isotropic aluminum plate, studied in Section 4.1, at 150 Hz. As shown in Fig. 1(a), when the wave propagation angle $\theta$ is given, the displacement $U(x_n, y_n)$ of any Point P in this direction can also be expressed as

$$ U(\theta, x_n, y_n) = \sum_{m=1}^{n_w} A_m e^{ik_{0,m}x_n} e^{ik_{y,m}y_n} \qquad (9) $$

In practice, the structural responses will inevitably be affected by various disturbances. To create a more realistic signal model, two disturbance factors are considered: (A) Small perturbation due to the geometrical coordinates' uncertainty; (B) External signal noise. Small perturbations are common in practice and can lead to mesh distortion, as shown in Fig. 1(b). The direct consequence is that the geometric coordinates of the measured points do not match the actual structural response. For example, under the influence of the small perturbation, Point P moves to Point P' with an unknown and random coordinate, leading to a mismatch between the

<figure><img src="image_6.png" /><figcaption>(a) Transformation of displacement and coordinates of measurement point P in Cartesian and polar coordinate systems</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(b) Comparison of a perfect 2D grid and a 2D distortion grid due to small perturbation problem</figcaption></figure>

<caption>**Fig. 1** Coordinate system transformation and 2D grid distortion.</caption>

actual response of the structure (the response of Point P') and the coordinates of the measured point (coordinate of Point P). Furthermore, in contrast to the 1D case, the small perturbation of the 2D measurement points is more complex, because the coordinates of the measurement points can be perturbed randomly in any direction. To solve this problem, in a probabilistic framework, the small perturbation problem can be quantified by introducing two random variables into the coordinate of each point.

$$ \begin{cases} \hat{x}_n = x_n \pm \xi_n \Delta x_n \\ \hat{y}_n = y_n \pm \zeta_n \Delta y_n \end{cases} \qquad n = 1, 2, \dots, N \qquad (10) $$

where random variables $\xi_n$ and $\zeta_n$ are defined as the small perturbation ratio which can obey different distributions, such as uniform distribution and gradient distribution. Therefore, the signal model in the presence of small perturbation can be expressed as

$$ U(\theta, \hat{x}_n, \hat{y}_n) = \sum_{m=1}^{n_w} A_m e^{ik_{0,m}(x_n \pm \xi_n \Delta x_n)} e^{ik_{y,m}(y_n \pm \zeta_n \Delta y_n)} \qquad (11) $$

where the small perturbation part can be expanded by a first-order Taylor series for illustration of measurement error caused by the small perturbation problem:
$$
\begin{equation}
\begin{aligned}
U(\theta, \hat{x}_n, \hat{y}_n) = & \sum_{m=1}^{n_w} A_m e^{jk_{\theta,x,m} \xi_n \Delta x_n} (1 \pm j k_{\theta,x,m} \xi_n \Delta y_n + O(\xi_n^2 \Delta x_n^2)) e^{jk_{\theta,y,m} y_n} \\
& \times (1 \pm j k_{\theta,y,m} \xi_n \Delta y_n + O(\xi_n^2 \Delta y_n^2)) \approx \sum_{m=1}^{n_w} A_m e^{jk_{\theta,x,m} \xi_n \Delta x_n} e^{jk_{\theta,y,m} y_n} \\
& \pm A_m e^{jk_{\theta,x,m} \xi_n \Delta x_n} e^{jk_{\theta,y,m} y_n} (1 + j k_{\theta,x,m} \xi_n \Delta x_n + j k_{\theta,y,m} \xi_n \Delta y_n \\
& + k_{\theta,x,m} \xi_n \Delta x_n k_{\theta,y,m} \xi_n \Delta y_n)
\end{aligned}
\tag{12}
\end{equation}
$$

where the first part is the unperturbed part, which is the same as the model in perfect condition, as shown in Eq. (9). The second part is the perturbed part which is not additive noise but depends on the degree of uncertainty of the geometric coordinates, the properties of the structure, and the frequency. It should be noted that in practice, the random variables $\xi_n$ and $\xi_n$ are unknown, and therefore, the random coordinate of each measurement point can only be considered as the unperturbed coordinate. Based on this, the stochastic signal model in the presence of small perturbation can be represented as

$$ \hat{U}(\theta, x_n, y_n) = U(\theta, \hat{x}_n, \hat{y}_n) \quad (13) $$

Then additive noise, such as external signal noise $\hat{\epsilon}(\theta, x_n, y_n)$, is considered. Thus the stochastic signal model is

$$ S(\theta, x_n, y_n) = \hat{U}(\theta, x_n, y_n) + \hat{\epsilon}(\theta, x_n, y_n) \quad (14) $$

The signal model described in Eq. (14) is built in the Cartesian coordinate system. As shown in Fig. 1(a), according to the transformation of the coordinate system, the coordinate of the measurement Point P can be expressed in the polar coordinate system:

$$ \begin{cases} x_n = r_n \cos \theta \\ y_n = r_n \sin \theta \end{cases} \quad n = 1, 2, \dots, N \quad (15) $$

where $r_n$ is the distance between the excitation Points O and P. The relationship between the wavenumber and its components in x and y directions are given by

$$ \begin{cases} k_{\theta,x,m} = k_{\theta,m} \cos \theta \\ k_{\theta,y,m} = k_{\theta,m} \sin \theta \end{cases} \quad n = 1, 2, \dots, N \quad (16) $$

The displacement of Point P in the polar coordinate system can be obtained by substituting Eqs. (15) and (16) into Eq. (9):

$$ U(\theta, r_n) = \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} r_n} \quad (17) $$

where the model is the signal model in perfect condition. The corresponding stochastic signal model is

$$ S(\theta, r_n) = \hat{U}(\theta, r_n) + \hat{\epsilon}(\theta, r_n) \quad (18) $$

where $\hat{\epsilon}(\theta, r_n)$ is the additional noise part; $\hat{U}(\theta, r_n)$ is the part affected by the small perturbation problem which can be expressed as

$$
\begin{align}
\hat{U}(\theta, r_n) &= \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} \hat{r}_n} \nonumber \\
&= \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} \cos \theta (r_n \cos \theta \pm \xi_n \Delta x_n)} e^{jk_{\theta,m} \sin \theta (r_n \sin \theta \pm \xi_n \Delta y_n)} \nonumber \\
&= \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} (r_n \pm \xi_n \Delta x_n \cos \theta \pm \xi_n \Delta y_n \sin \theta)} \tag{19}
\end{align}
$$

where it is clear that $\hat{r}_n$ is related to $\xi_n$ and $\xi_n$.

## 3.2. AKSI formula in polar coordinate system

Introduce the formula of AKSI on wavenumber identification of a plate in the polar coordinate system framework. To facilitate the expression of the general formula for AKSI, the signal noise part $\hat{\epsilon}(\theta, r)$ is ignored for the moment. Thus, Eq. (18) can be rewritten as

$$ S(\theta, r_n) = \sum_{m=1}^{n_w} A_m e^{R_m r_n} \quad (20) $$

where $R_m$ is equal to $jk_{\theta,m}$. In the wavenumber domain, Eq. (20) can be written in the wavenumber domain by Laplace transform:

$$ S(\theta, s) = \frac{A_1}{s - R_1} + \frac{A_2}{s - R_2} + \dots + \frac{A_{n_w}}{s - R_{n_w}} \quad (21) $$

where $s$ is the variable in the wavenumber domain.

It is obvious that Eq. (20) can be regarded as a solution to an ordinary differential equation. Thus, in the wavenumber domain, the characteristic polynomial of this ODE can be defined as

$$ \Psi(s) = \prod_{m=1}^{n_w} (s - R_m) = \sum_{i=0}^{n_w} \gamma(n_w - i) s^i \quad (22) $$

where $\gamma(n_w - i)_{[i \in [0, n_w]]}$ are the unknown coefficients of the characteristic polynomial. It is to be noted that if these unknown coefficients are estimated, the wavenumber can be calculated by solving this characteristic polynomial. Therefore, the following paragraphs start with the estimation of the coefficients $\gamma(n_w - i)_{[i \in [0, n_w]]}$ in three steps based on the algebraic identification scheme: (A) By Laplace transform, an ordinary differential equation with unknown coefficients is established in the wavenumber domain; (B) Inverse Laplace transform is applied to the differential equation to obtain the corresponding expression in the spatial domain; (C) Least-squares method is used to estimate the unknown coefficients, which are closely related to the wavenumber.

The first step is establishing an ODE in the wavenumber domain. A polynomial formula can be obtained by multiplying Eqs. (21) and (22) as

$$
\begin{align}
S(\theta, s) \Psi(s) &= \left( \frac{A_1}{s - R_1} + \frac{A_2}{s - R_2} + \dots + \frac{A_{n_w}}{s - R_{n_w}} \right) \prod_{m=1}^{n_w} (s - R_m) \nonumber \\
&= \sum_{m=1}^{n_w} A_m \prod_{i=1, i \neq m}^{n_w} (s - R_i) \tag{23}
\end{align}
$$

It is easy to see that Eq. (23) is a polynomial of order $n_w - 1$. Thus, a differential equation can be built by taking the derivative of Eq. (23) with respect to the variable $sn_w$ times:

$$ \frac{d^{n_w} S(\theta, s) \Psi(s)}{ds^{n_w}} = \frac{d^{n_w} [S(\theta, s) \sum_{i=0}^{n_w} \gamma(n_w - i) s^i]}{ds^{n_w}} \quad (24) $$

To calculate Eq. (24), two formulas are presented as

$$ \frac{d^{n_w} S(\theta, s) \Psi(s)}{ds^{n_w}} = \sum_{j=0}^{n_w} \binom{n_w}{j} \frac{d^{n_w-j} S(\theta, s)}{ds^{n_w-j}} \cdot \frac{d^j \Psi(s)}{ds^j} \quad (25) $$

$$ \frac{d^{n_w} (s^j)}{ds^{n_w}} = \frac{j!}{(j - n_w)!} s^{j - n_w} \quad (26) $$

where $\binom{n_w}{j} = \frac{n_w!}{j!(n_w-j)!}$ with $n_w!$ is the factorial of $n_w$. Eq. (25) is a Leibniz formula and Eq. (26) is a formula for higher-order algebraic derivatives. Substituting these two formulas into Eq. (24), Eq. (27) can be obtained:

$$ \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! s^{j-i} \frac{d^i S(\theta, s)}{ds^i} \gamma(i) = 0 \quad (27) $$
$$
\begin{equation}
\begin{aligned}
U(\theta, \hat{x}_n, \hat{y}_n) = & \sum_{m=1}^{n_w} A_m e^{jk_{\theta,x,m} \xi_n \Delta x_n} (1 \pm j k_{\theta,x,m} \xi_n \Delta y_n + O(\xi_n^2 \Delta x_n^2)) e^{jk_{\theta,y,m} y_n} \\
& \times (1 \pm j k_{\theta,y,m} \xi_n \Delta y_n + O(\xi_n^2 \Delta y_n^2)) \approx \sum_{m=1}^{n_w} A_m e^{jk_{\theta,x,m} \xi_n \Delta x_n} e^{jk_{\theta,y,m} y_n} \\
& \pm A_m e^{jk_{\theta,x,m} \xi_n \Delta x_n} e^{jk_{\theta,y,m} y_n} (1 + j k_{\theta,x,m} \xi_n \Delta x_n + j k_{\theta,y,m} \xi_n \Delta y_n + k_{\theta,x,m} \xi_n \Delta x_n k_{\theta,y,m} \xi_n \Delta y_n)
\end{aligned}
\tag{12}
\end{equation}
$$

where the first part is the unperturbed part, which is the same as the model in perfect condition, as shown in Eq. (9). The second part is the perturbed part which is not additive noise but depends on the degree of uncertainty of the geometric coordinates, the properties of the structure, and the frequency. It should be noted that in practice, the random variables $\xi_n$ and $\xi_n$ are unknown, and therefore, the random coordinate of each measurement point can only be considered as the unperturbed coordinate. Based on this, the stochastic signal model in the presence of small perturbation can be represented as

$$ \hat{U}(\theta, x_n, y_n) = U(\theta, \hat{x}_n, \hat{y}_n) \quad (13) $$

Then additive noise, such as external signal noise $\varepsilon(\theta, x_n, y_n)$, is considered. Thus the stochastic signal model is

$$ S(\theta, x_n, y_n) = \hat{U}(\theta, x_n, y_n) + \varepsilon(\theta, x_n, y_n) \quad (14) $$

The signal model described in Eq. (14) is built in the Cartesian coordinate system. As shown in Fig. 1(a), according to the transformation of the coordinate system, the coordinate of the measurement Point P can be expressed in the polar coordinate system:

$$ \begin{cases} x_n = r_n \cos \theta \\ y_n = r_n \sin \theta \end{cases} \quad n = 1, 2, \dots, N \quad (15) $$

where $r_n$ is the distance between the excitation Points O and P. The relationship between the wavenumber and its components in x and y directions are given by

$$ \begin{cases} k_{\theta,x,m} = k_{\theta,m} \cos \theta \\ k_{\theta,y,m} = k_{\theta,m} \sin \theta \end{cases} \quad n = 1, 2, \dots, N \quad (16) $$

The displacement of Point P in the polar coordinate system can be obtained by substituting Eqs. (15) and (16) into Eq. (9):

$$ U(\theta, r_n) = \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} r_n} \quad (17) $$

where the model is the signal model in perfect condition. The corresponding stochastic signal model is

$$ S(\theta, r_n) = \hat{U}(\theta, r_n) + \varepsilon(\theta, r_n) \quad (18) $$

where $\varepsilon(\theta, r_n)$ is the additional noise part; $\hat{U}(\theta, r_n)$ is the part affected by the small perturbation problem which can be expressed as

$$
\begin{align}
\hat{U}(\theta, r_n) &= \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} \hat{r}_n} \nonumber \\
&= \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} \cos \theta (r_n \cos \theta \pm \xi_n \Delta x_n)} e^{jk_{\theta,m} \sin \theta (r_n \sin \theta \pm \xi_n \Delta y_n)} \nonumber \\
&= \sum_{m=1}^{n_w} A_m e^{jk_{\theta,m} (r_n \pm \xi_n \Delta x_n \cos \theta \pm \xi_n \Delta y_n \sin \theta)} \tag{19}
\end{align}
$$

where it is clear that $\hat{r}_n$ is related to $\xi_n$ and $\xi_n$.

## 3.2. AKSI formula in polar coordinate system

Introduce the formula of AKSI on wavenumber identification of a plate in the polar coordinate system framework. To facilitate the expression of the general formula for AKSI, the signal noise part $\varepsilon(\theta, r)$ is ignored for the moment. Thus, Eq. (18) can be rewritten as

$$ S(\theta, r_n) = \sum_{m=1}^{n_w} A_m e^{R_m r_n} \quad (20) $$

where $R_m$ is equal to $jk_{\theta,m}$. In the wavenumber domain, Eq. (20) can be written in the wavenumber domain by Laplace transform:

$$ S(\theta, s) = \frac{A_1}{s - R_1} + \frac{A_2}{s - R_2} + \cdots + \frac{A_{n_w}}{s - R_{n_w}} \quad (21) $$

where $s$ is the variable in the wavenumber domain.

It is obvious that Eq. (20) can be regarded as a solution to an ordinary differential equation. Thus, in the wavenumber domain, the characteristic polynomial of this ODE can be defined as

$$ \Psi(s) = \prod_{m=1}^{n_w} (s - R_m) = \sum_{i=0}^{n_w} \gamma(n_w - i) s^i \quad (22) $$

where $\gamma(n_w - i)_{[i \in [0, n_w]]}$ are the unknown coefficients of the characteristic polynomial. It is to be noted that if these unknown coefficients are estimated, the wavenumber can be calculated by solving this characteristic polynomial. Therefore, the following paragraphs start with the estimation of the coefficients $\gamma(n_w - i)_{[i \in [0, n_w]]}$ in three steps based on the algebraic identification scheme: (A) By Laplace transform, an ordinary differential equation with unknown coefficients is established in the wavenumber domain; (B) Inverse Laplace transform is applied to the differential equation to obtain the corresponding expression in the spatial domain; (C) Least-squares method is used to estimate the unknown coefficients, which are closely related to the wavenumber.

The first step is establishing an ODE in the wavenumber domain. A polynomial formula can be obtained by multiplying Eqs. (21) and (22) as

$$
\begin{align}
S(\theta, s) \Psi(s) &= \left( \frac{A_1}{s - R_1} + \frac{A_2}{s - R_2} + \cdots + \frac{A_{n_w}}{s - R_{n_w}} \right) \prod_{m=1}^{n_w} (s - R_m) \nonumber \\
&= \sum_{m=1}^{n_w} A_m \prod_{i=1, i \neq m}^{n_w} (s - R_i) \tag{23}
\end{align}
$$

It is easy to see that Eq. (23) is a polynomial of order $n_w - 1$. Thus, a differential equation can be built by taking the derivative of Eq. (23) with respect to the variable $sn_w$ times:

$$ \frac{d^{n_w} S(\theta, s) \Psi(s)}{ds^{n_w}} = \frac{d^{n_w} [S(\theta, s) \sum_{i=0}^{n_w} \gamma(n_w - i) s^i]}{ds^{n_w}} \quad (24) $$

To calculate Eq. (24), two formulas are presented as

$$ \frac{d^{n_w} S(\theta, s) \Psi(s)}{ds^{n_w}} = \sum_{j=0}^{n_w} \binom{n_w}{j} \frac{d^{n_w-j} S(\theta, s)}{ds^{n_w-j}} \cdot \frac{d^j \Psi(s)}{ds^j} \quad (25) $$

$$ \frac{d^{n_w} (s^j)}{ds^{n_w}} = \frac{j!}{(j - n_w)!} s^{j - n_w} \quad (26) $$

where $\binom{n_w}{j} = \frac{n_w!}{j!(n_w-j)!}$ with $n_w!$ is the factorial of $n_w$. Eq. (25) is a Leibniz formula and Eq. (26) is a formula for higher-order algebraic derivatives. Substituting these two formulas into Eq. (24), Eq. (27) can be obtained:

$$ \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! s^{j-i} \frac{d^i S(\theta, s)}{ds^i} \gamma(i) = 0 \quad (27) $$
Now, the second step needs to proceed, which aims to return Eq. (27) into the spatial domain using the inverse Laplace transform. For this purpose, the following inverse Laplace transform is introduced:

$$ \mathcal{L}^{-1}\left(\frac{1}{s^l} \cdot \frac{d^l S(\theta, s)}{ds^l}\right) = \frac{1}{(I-1)!} \int_0^{r_n} v_{I-1, J}(\tau) S(\theta, \tau) d\tau \quad (28) $$

with

$$ v_{I,J}(\tau) = (r_n - \tau)^l (-\tau)^J \quad (29) $$

where $\tau$ is integral variable; $J$ and $I$ are two indices, representing the order of the derivative and the exponent of the variable $s$, respectively.

To satisfy inverse Laplace transform Eq. (28), Eq. (27) requires first dividing by $s^{n_w+1}$:

$$ \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! \frac{1}{s^{n_w+1+i-j}} \cdot \frac{d^l S(\theta, s)}{ds^l} \gamma^l(i) = 0 \quad (30) $$

Then, after applying the inverse Laplace transform Eq. (28) to Eq. (30), one can obtain Eq. (31) in the spatial domain:

$$ \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! \frac{1}{(n_w + i - j)!} \times \int_0^{r_n} (r_n - \tau)^{n_w + i - j} (-\tau)^j S(\theta, \tau) d\tau \gamma^j(i) = 0 \quad (31) $$

Eq. (31) can also be expressed as

$$ \sum_{i=0}^{n_w} \phi(i, \theta, r_n) \gamma^l(i) = 0 \quad (32) $$

with

$$ \phi(i, \theta, r_n) = \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! \frac{1}{(n_w + i - j)!} \times \int_0^{r_n} (r_n - \tau)^{n_w + i - j} (-\tau)^j S(\theta, \tau) d\tau \quad (33) $$

where the integrals can be calculated using numerical integration, such as the trapezoidal rule.

Finally, the third step is to estimate $\gamma^l(i)$ of Eq. (32) using the least-squares method. Eq. (32) can be expressed in the matrix format as

$$ HX = M \quad (34) $$

with

$$ \left\{ \begin{array}{l} H = \begin{bmatrix} \phi(n_w, \theta, r_1) & \phi(n_w - 1, \theta, r_1) & \cdots & \phi(0, \theta, r_1) \\ \phi(n_w, \theta, r_2) & \phi(n_w - 1, \theta, r_2) & \cdots & \phi(0, \theta, r_2) \\ \vdots & \vdots & \ddots & \vdots \\ \phi(n_w, \theta, r_N) & \phi(n_w - 1, \theta, r_N) & \cdots & \phi(0, \theta, r_N) \end{bmatrix} \\ X = \begin{bmatrix} \gamma(n_w) \\ \gamma(n_w - 1) \\ \vdots \\ \gamma(0) \end{bmatrix} \\ M = \begin{bmatrix} 0 \\ 0 \\ \vdots \\ 0 \end{bmatrix} \end{array} \right. \quad (35) $$

where $X$ is the eigenvector corresponding to the smallest eigenvalue of $H^T H$. Once coefficients vector $X$ of the characteristic polynomial are obtained, the $R_m$ can be calculated by Eq. (22). Then the wavenumber $k_{0,m}$ can be obtained by $k_{0,m} = -jR_m$.

It is worth noting that the integrals of the AKSI formula resemble low-pass filter,³⁰,³¹ which can reduce the influence of the measurement error, such as signal noise and the small perturbation problem. Moreover, the AKSI treats the signals as a continuous function. Therefore, it is not limited to periodic sampling. On the other hand, AKSI requires only simple linear problems to be solved, thus reducing the computational cost. In addition, this method is an $\theta$-dependent method that requires only a small number of samples on the line grid as input parameters, thus reducing the time needed to arrange the 2D measurement points in the experimental situation. The sampling characteristics of AKSI increase its applicability in periodic structures with unknown periodicity.

### 3.3. Model order estimation and wavenumber filtering

In the theory of AKSI in Section 3.2, the model order $n_w$ (the number of waves) is assumed to be known a priori. However, this parameter is unknown in most practical cases. Because incorrect estimation of the model order affects the accuracy of the extracted wavenumber, the signal requires pre-processing (model order estimation). Many model order estimation methods are available. In general terms, they can be divided into two categories: information-theoretic criteri,³²,³³ and subspace-based methods.³⁴⁻³⁶ The first category of methods can be applied to signals obtained by any sampling way, but it is sensitive to signal noise. Comparatively, the second category of methods is robust to signal noise but only works with periodic signals. Maximum Description Length (MDL) and ESTimation of Error (ESTER), the representative of these two categories of model order estimation methods, are used to constitute a model order strategy. As the theory of MDL³³ and ESTER³⁴ has been well established, focused on introducing the model order estimation strategy, which is suitable for different stochastic conditions:

* Periodic samples: ESTER is preferred because of its robustness to noise. The interested reader can be referred to Refs. 23,37 where ESTER has been used in the HRWA method for estimating the number of waves in one- and two-dimensional structures.

* Non-periodic samples less affected by perturbations: MDL is used to estimate the model order. MDL has been used to estimate the number of waves.³⁸

* Non-periodic samples affected by perturbations in a bigger degree: the interpolation method is first used to generate the periodic samples from the non-periodic samples, and then ESTER is used to estimate the model order.

* Non-periodic samples or periodic samples strongly affected by perturbations: the most effective way is first to divide the frequency band into several sub-frequency bands, and then manually adjust the model order in each sub-frequency band until a better dispersion curve is obtained, and finally repeat the second step in the other frequency bands.

The application of ESTER on model order estimation is illustrated by using the extraction of the wavenumber of the
Now, the second step needs to proceed, which aims to return Eq. (27) into the spatial domain using the inverse Laplace transform. For this purpose, the following inverse Laplace transform is introduced:

$$ \mathcal{L}^{-1}\left(\frac{1}{s^l} \cdot \frac{d^l S(\theta, s)}{ds^l}\right) = \frac{1}{(I-1)!} \int_0^{r_n} v_{I-1, J}(\tau) S(\theta, \tau) d\tau \quad (28) $$

with

$$ v_{I,J}(\tau) = (r_n - \tau)^l (-\tau)^J \quad (29) $$

where $\tau$ is integral variable; $J$ and $I$ are two indices, representing the order of the derivative and the exponent of the variable $s$, respectively.

To satisfy inverse Laplace transform Eq. (28), Eq. (27) requires first dividing by $s^{n_w+1}$:

$$ \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! \frac{1}{s^{n_w+1+i-j}} \cdot \frac{d^l S(\theta, s)}{ds^l} \gamma^l(i) = 0 \quad (30) $$

Then, after applying the inverse Laplace transform Eq. (28) to Eq. (30), one can obtain Eq. (31) in the spatial domain:

$$ \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! \frac{1}{(n_w + i - j)!} \times \int_0^{r_n} (r_n - \tau)^{n_w + i - j} (-\tau)^j S(\theta, \tau) d\tau \gamma^j(i) = 0 \quad (31) $$

Eq. (31) can also be expressed as

$$ \sum_{i=0}^{n_w} \phi(i, \theta, r_n) \gamma^l(i) = 0 \quad (32) $$

with

$$ \phi(i, \theta, r_n) = \sum_{i=0}^{n_w} \sum_{j=i}^{n_w} \binom{n_w}{j} \binom{n_w - i}{n_w - j} (n_w - j)! \frac{1}{(n_w + i - j)!} \times \int_0^{r_n} (r_n - \tau)^{n_w + i - j} (-\tau)^j S(\theta, \tau) d\tau \quad (33) $$

where the integrals can be calculated using numerical integration, such as the trapezoidal rule.

Finally, the third step is to estimate $\gamma^l(i)$ of Eq. (32) using the least-squares method. Eq. (32) can be expressed in the matrix format as

$$ HX = M \quad (34) $$

with

$$ \left\{ \begin{array}{l} H = \begin{bmatrix} \phi(n_w, \theta, r_1) & \phi(n_w - 1, \theta, r_1) & \cdots & \phi(0, \theta, r_1) \\ \phi(n_w, \theta, r_2) & \phi(n_w - 1, \theta, r_2) & \cdots & \phi(0, \theta, r_2) \\ \vdots & \vdots & \ddots & \vdots \\ \phi(n_w, \theta, r_N) & \phi(n_w - 1, \theta, r_N) & \cdots & \phi(0, \theta, r_N) \end{bmatrix} \\ X = \begin{bmatrix} \gamma(n_w) \\ \gamma(n_w - 1) \\ \vdots \\ \gamma(0) \end{bmatrix} \\ M = \begin{bmatrix} 0 \\ 0 \\ \vdots \\ 0 \end{bmatrix} \end{array} \right. \quad (35) $$

where $X$ is the eigenvector corresponding to the smallest eigenvalue of $H^T H$. Once coefficients vector $X$ of the characteristic polynomial are obtained, the $R_m$ can be calculated by Eq. (22). Then the wavenumber $k_{0,m}$ can be obtained by $k_{0,m} = -jR_m$.

It is worth noting that the integrals of the AKSI formula resemble low-pass filter,³⁰,³¹ which can reduce the influence of the measurement error, such as signal noise and the small perturbation problem. Moreover, the AKSI treats the signals as a continuous function. Therefore, it is not limited to periodic sampling. On the other hand, AKSI requires only simple linear problems to be solved, thus reducing the computational cost. In addition, this method is an $\theta$-dependent method that requires only a small number of samples on the line grid as input parameters, thus reducing the time needed to arrange the 2D measurement points in the experimental situation. The sampling characteristics of AKSI increase its applicability in periodic structures with unknown periodicity.

### 3.3. Model order estimation and wavenumber filtering

In the theory of AKSI in Section 3.2, the model order $n_w$ (the number of waves) is assumed to be known a priori. However, this parameter is unknown in most practical cases. Because incorrect estimation of the model order affects the accuracy of the extracted wavenumber, the signal requires pre-processing (model order estimation). Many model order estimation methods are available. In general terms, they can be divided into two categories: information-theoretic criteri,³²,³³ and subspace-based methods.³⁴⁻³⁶ The first category of methods can be applied to signals obtained by any sampling way, but it is sensitive to signal noise. Comparatively, the second category of methods is robust to signal noise but only works with periodic signals. Maximum Description Length (MDL) and ESTimation of Error (ESTER), the representative of these two categories of model order estimation methods, are used to constitute a model order strategy. As the theory of MDL³³ and ESTER³⁴ has been well established, focused on introducing the model order estimation strategy, which is suitable for different stochastic conditions:

* Periodic samples: ESTER is preferred because of its robustness to noise. The interested reader can be referred to Refs. 23,37 where ESTER has been used in the HRWA method for estimating the number of waves in one- and two-dimensional structures.

* Non-periodic samples less affected by perturbations: MDL is used to estimate the model order. MDL has been used to estimate the number of waves.³⁸

* Non-periodic samples affected by perturbations in a bigger degree: the interpolation method is first used to generate the periodic samples from the non-periodic samples, and then ESTER is used to estimate the model order.

* Non-periodic samples or periodic samples strongly affected by perturbations: the most effective way is first to divide the frequency band into several sub-frequency bands, and then manually adjust the model order in each sub-frequency band until a better dispersion curve is obtained, and finally repeat the second step in the other frequency bands.

The application of ESTER on model order estimation is illustrated by using the extraction of the wavenumber of the
isotropic plate studied in Section 4.1 as an example. Blue asterisks in Fig. 2 show the estimation results when ESTER is applied to the 50 periodic displacements along x direction at 8.5 kHz when $n_w \in [1,20]$. It can be seen that the ESTER value reaches the minimum when $n_w = 8$. Thus the model order at this frequency is taken as 8. The zoomed sub-pictures show that the dispersion curve between 8.0 kHz and 9.0 kHz obtained by AKSI with $n_w = 8$ is accurate, while AKSI fails to estimate wavenumber when $n_w = 4$, which also illustrates that ESTER can provide reliable model order estimation for AKSI.

The value of model order may vary with frequency, resulting in the number of extracted wavenumbers may also vary at different frequencies. Therefore, the dominant wavenumber $k_{0,m}$ of propagating waves needs to be selected through wavenumber filtering:

$$k_{0,m} = \begin{cases} \mathrm{Re}(k_{0,m}) > 0 \\ \mathrm{Im}(k_{0,m}) > 0 \\ |\mathrm{Re}(k_{0,m})| > |\mathrm{Im}(k_{0,m})| \end{cases} \quad (36)$$

where the first condition reflects the characteristic of the propagating wave, and the second condition aims to filter out evanescent waves.

### 3.4. Implementation of AKSI procedure

The flowchart for extracting dispersion curves and k-space using the AKSI procedure is shown in Fig. 3. The good performance of the AKSI algorithm under stochastic conditions is attributed to the nature of the algebraic identification estimator, which has been proven to be a fast algebraic estimation process with high robustness to signal noise and structured perturbations.28,30 It is worth noting that AKSI has only been able to identify the imaginary part of the wavenumbers of plane waves or approximated plane waves. This limitation is because the signal model used in AKSI is a plane wave model, as shown in Eq. (8). A characteristic of the plane wave is that the wavefront is a straight line. However, the wave field of the 2D structure often presents a non-straight wavefront, especially for low-damped thin plates. As a result, the 2D non-plane wave signal cannot be represented by a sum of damped exponential functions, resulting in the fact that the imaginary part of the wavenumbers of 2D structures remains the most

challenging parameter to be identified for the plane wave model-based inverse methods. This effect comes mainly from the fact that the imaginary part of the wavenumber is related to the decay of wave propagation, which is influenced by a variety of factors, such as interference from evanescent waves, energy transformation in different media, and the geometry of the sources, leading to the challenge of accurately extracting the imaginary part of wavenumber. Assessing the impact of these factors on the extraction of the imaginary part is outside the scope of this paper. For this reason, AKSI is only applied to extract the imaginary part of the wavenumber of the Levy-Type (LT) periodic plate studied in Section 4.3. Due to the high damping effect of LT periodic structures, the displacement field has an approximately straight wavefront.

## 4. Numerical study and parametric survey

The performance of AKSI in extracting dispersion curves and k-space of 2D structures under stochastic conditions is verified by a series of numerical cases. The stochastic conditions here include: (A) High signal noise; (B) Small perturbation problems; (C) Non-periodic sampling; (D) Periodic structures with approximately known and unknown periodicity. In addition, we compare AKSI with three other popular inverse methods, including IWC, INCOME, and BWI, to show that AKSI overcomes some main limitations of these inverse methods. We organized various case studies in Fig. 4 in which the type of plates tested, the purpose of the tests, the stochastic conditions tested, and the methods compared are all mentioned.

To facilitate obtaining disturbed 2D displacement fields under signal noise and small perturbation conditions, the transverse displacements of the plates in the first two numerical cases are calculated by analytical methods. In contrast, in the third numerical case, the transverse displacements of the LT periodic plate are obtained with the finite element software COMSOL for ease of modeling. The inverse methods are compared objectively based on the same data set from the 2D grid in each case. The k-space is formed from the extracted wavenumber of 24 angles in the first two case studies, as shown in Fig. 5, where these angles ensure that AKSI can use the displacements measured on the 2D grid as input parameters. These angles in the first quadrant are in the following order: $\theta_1 = 0^\circ$, $\theta_2 = 18.4349^\circ$, $\theta_3 = 26.5651^\circ$, $\theta_4 = 45^\circ$, $\theta_5 = 63.4349^\circ$, $\theta_6 = 71.5651^\circ$, and $\theta_7 = 90^\circ$. The angles in the

<figure><img src="image_8.png" /><figcaption>Fig. 2 Model order estimation.</figcaption></figure>
isotropic plate studied in Section 4.1 as an example. Blue asterisks in Fig. 2 show the estimation results when ESTER is applied to the 50 periodic displacements along x direction at 8.5 kHz when $n_w \in [1,20]$. It can be seen that the ESTER value reaches the minimum when $n_w = 8$. Thus the model order at this frequency is taken as 8. The zoomed sub-pictures show that the dispersion curve between 8.0 kHz and 9.0 kHz obtained by AKSI with $n_w = 8$ is accurate, while AKSI fails to estimate wavenumber when $n_w = 4$, which also illustrates that ESTER can provide reliable model order estimation for AKSI.

The value of model order may vary with frequency, resulting in the number of extracted wavenumbers may also vary at different frequencies. Therefore, the dominant wavenumber $k_{0,m}$ of propagating waves needs to be selected through wavenumber filtering:

$$k_{0,m} = \begin{cases} \mathrm{Re}(k_{0,m}) > 0 \\ \mathrm{Im}(k_{0,m}) > 0 \\ |\mathrm{Re}(k_{0,m})| > |\mathrm{Im}(k_{0,m})| \end{cases} \quad (36)$$

where the first condition reflects the characteristic of the propagating wave, and the second condition aims to filter out evanescent waves.

### 3.4. Implementation of AKSI procedure

The flowchart for extracting dispersion curves and k-space using the AKSI procedure is shown in Fig. 3. The good performance of the AKSI algorithm under stochastic conditions is attributed to the nature of the algebraic identification estimator, which has been proven to be a fast algebraic estimation process with high robustness to signal noise and structured perturbations.28,30 It is worth noting that AKSI has only been able to identify the imaginary part of the wavenumbers of plane waves or approximated plane waves. This limitation is because the signal model used in AKSI is a plane wave model, as shown in Eq. (8). A characteristic of the plane wave is that the wavefront is a straight line. However, the wave field of the 2D structure often presents a non-straight wavefront, especially for low-damped thin plates. As a result, the 2D non-plane wave signal cannot be represented by a sum of damped exponential functions, resulting in the fact that the imaginary part of the wavenumbers of 2D structures remains the most

challenging parameter to be identified for the plane wave model-based inverse methods. This effect comes mainly from the fact that the imaginary part of the wavenumber is related to the decay of wave propagation, which is influenced by a variety of factors, such as interference from evanescent waves, energy transformation in different media, and the geometry of the sources, leading to the challenge of accurately extracting the imaginary part of wavenumber. Assessing the impact of these factors on the extraction of the imaginary part is outside the scope of this paper. For this reason, AKSI is only applied to extract the imaginary part of the wavenumber of the Levy-Type (LT) periodic plate studied in Section 4.3. Due to the high damping effect of LT periodic structures, the displacement field has an approximately straight wavefront.

## 4. Numerical study and parametric survey

The performance of AKSI in extracting dispersion curves and k-space of 2D structures under stochastic conditions is verified by a series of numerical cases. The stochastic conditions here include: (A) High signal noise; (B) Small perturbation problems; (C) Non-periodic sampling; (D) Periodic structures with approximately known and unknown periodicity. In addition, we compare AKSI with three other popular inverse methods, including IWC, INCOME, and BWI, to show that AKSI overcomes some main limitations of these inverse methods. We organized various case studies in Fig. 4 in which the type of plates tested, the purpose of the tests, the stochastic conditions tested, and the methods compared are all mentioned.

To facilitate obtaining disturbed 2D displacement fields under signal noise and small perturbation conditions, the transverse displacements of the plates in the first two numerical cases are calculated by analytical methods. In contrast, in the third numerical case, the transverse displacements of the LT periodic plate are obtained with the finite element software COMSOL for ease of modeling. The inverse methods are compared objectively based on the same data set from the 2D grid in each case. The k-space is formed from the extracted wavenumber of 24 angles in the first two case studies, as shown in Fig. 5, where these angles ensure that AKSI can use the displacements measured on the 2D grid as input parameters. These angles in the first quadrant are in the following order: $\theta_1 = 0^\circ$, $\theta_2 = 18.4349^\circ$, $\theta_3 = 26.5651^\circ$, $\theta_4 = 45^\circ$, $\theta_5 = 63.4349^\circ$, $\theta_6 = 71.5651^\circ$, and $\theta_7 = 90^\circ$. The angles in the

<figure><img src="image_8.png" /><figcaption>Fig. 2 Model order estimation.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>Fig. 3 Flowchart for dispersion curve and k-space identification using AKSI procedure.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 4 Tree diagram of numerical study.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>Fig. 3 Flowchart for dispersion curve and k-space identification using AKSI procedure.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 4 Tree diagram of numerical study.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>Fig. 5 Schematic representation of angle θ distribution of k-space extracted in numerical study.</figcaption></figure>

other quadrants can be obtained symmetrically from those in the first quadrant.

## 4.1. Numerical study 1: An isotropic plate

AKSI is used to extract dispersion curves and k-space of a simply supported thin isotropic plate subjected to a point load under perfect and signal noise conditions. Section 4.1.1 compares the performance of AKSI, INCOME, and IWC under perfect conditions, whereas Section 4.1.2 compares the performance of these three inverse methods under signal noise conditions.

The structure is a 200 cm × 200 cm aluminum plate (elastic modulus E = 79 GPa, Poisson's ratio v = 0.3, density of plate ρ = 2700 kg/m³, damping loss factor η = 0.005) with the thickness of 0.17 cm. The point load is put on the center of the plate. The transverse displacement U(x, y, ω) of this Kirchhoff plate can be calculated using³⁹

$$U(x, y, \omega) = \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{F_0 \Phi_{nm}(x_0, y_0) \Phi_{nm}(x, y)}{[\pi^4 D (\frac{n_1^2}{l_1^2} + \frac{n_2^2}{l_2^2})^2 - \rho h \omega^2] \int_0^{l_1} \int_0^{l_2} \Phi_{nm}^2(x, y) dx dy} \quad (37)$$

with

$$D = \frac{E(1 + i\eta)h^3}{12(1 - \nu^2)} \quad (38)$$

$$\Phi_{nm}(x, y) = \sin \frac{n\pi x}{l_1} \sin \frac{m\pi y}{l_2} \quad (39)$$

where F₀ the force and it is defined as 10 N in this case; l₁, l₂, and h are the length, width, and thickness of the plate, respectively; ω is angular frequency; D is bending stiffness; Φₘₘ(x, y) is eigen-modes. The transverse displacements calculated on a 501 × 501 mesh grid with the 4 mm sampling interval along x and y directions are used as the input parameters of INCOME and IWC. AKSI was performed on 50 measurement points distributed along each direction.

<figure><img src="image_12.png" /><figcaption>Fig. 6 Comparison of dispersion curve and k-space at 500 Hz obtained by IWC, INCOME, and AKSI for isotropic aluminum plate under perfect condition.</figcaption></figure>

### 4.1.1. Comparison of inverse methods under perfect condition

Under the perfect condition, the dispersion curves in x direction extracted by INCOME, IWC, and AKSI are shown in Fig. 6(a). INCOME and AKSI can provide accurate dispersion curves in the full frequency range, while IWC performs better at higher frequencies. This is because that displacement cannot contain many wavelengths in the low-frequency range, resulting in fluctuations in the dispersion curve obtained by IWC. The k-space comparison at 500 Hz is provided in Fig. 6(b), where all three tested inverse methods have a good performance on k-space identification, although the k-space obtained by IWC is not as smooth as that extracted from INCOME and AKSI. It is noted that INCOME is a method that can provide wavenumber with high numerical precision under the perfect condition²⁴ but has a poor performance in noisy environments, especially when the signal is disturbed by the high level of signal noise.

### 4.1.2. Comparison of inverse methods under signal noise condition

Fig. 7 shows the wavenumbers extracted by AKSI when considering displacements affected by different signal noise levels as input parameters. The additive noise here is Gaussian white noise with the Signal-to-Noise Ratio (SNR) generated by the AWGN function in MATLAB and it can be measured by SNR. Fig. 7 shows that the dispersion curve obtained by AKSI is reliable when SNR is higher than 20, having a good agreement with the analytical solution.
<figure><img src="image_11.png" /><figcaption>Fig. 5 Schematic representation of angle θ distribution of k-space extracted in numerical study.</figcaption></figure>

other quadrants can be obtained symmetrically from those in the first quadrant.

## 4.1. Numerical study 1: An isotropic plate

AKSI is used to extract dispersion curves and k-space of a simply supported thin isotropic plate subjected to a point load under perfect and signal noise conditions. Section 4.1.1 compares the performance of AKSI, INCOME, and IWC under perfect conditions, whereas Section 4.1.2 compares the performance of these three inverse methods under signal noise conditions.

The structure is a 200 cm × 200 cm aluminum plate (elastic modulus E = 79 GPa, Poisson's ratio v = 0.3, density of plate ρ = 2700 kg/m³, damping loss factor η = 0.005) with the thickness of 0.17 cm. The point load is put on the center of the plate. The transverse displacement U(x, y, ω) of this Kirchhoff plate can be calculated using³⁹

$$U(x, y, \omega) = \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{F_0 \Phi_{nm}(x_0, y_0) \Phi_{nm}(x, y)}{[\pi^4 D (\frac{n_1^2}{l_1^2} + \frac{n_2^2}{l_2^2})^2 - \rho h \omega^2] \int_0^{l_1} \int_0^{l_2} \Phi_{nm}^2(x, y) dx dy} \quad (37)$$

with

$$D = \frac{E(1 + i\eta)h^3}{12(1 - \nu^2)} \quad (38)$$

$$\Phi_{nm}(x, y) = \sin \frac{n\pi x}{l_1} \sin \frac{m\pi y}{l_2} \quad (39)$$

where F₀ the force and it is defined as 10 N in this case; l₁, l₂, and h are the length, width, and thickness of the plate, respectively; ω is angular frequency; D is bending stiffness; Φₘₘ(x, y) is eigen-modes. The transverse displacements calculated on a 501 × 501 mesh grid with the 4 mm sampling interval along x and y directions are used as the input parameters of INCOME and IWC. AKSI was performed on 50 measurement points distributed along each direction.

<figure><img src="image_12.png" /><figcaption>Fig. 6 Comparison of dispersion curve and k-space at 500 Hz obtained by IWC, INCOME, and AKSI for isotropic aluminum plate under perfect condition.</figcaption></figure>

### 4.1.1. Comparison of inverse methods under perfect condition

Under the perfect condition, the dispersion curves in x direction extracted by INCOME, IWC, and AKSI are shown in Fig. 6(a). INCOME and AKSI can provide accurate dispersion curves in the full frequency range, while IWC performs better at higher frequencies. This is because that displacement cannot contain many wavelengths in the low-frequency range, resulting in fluctuations in the dispersion curve obtained by IWC. The k-space comparison at 500 Hz is provided in Fig. 6(b), where all three tested inverse methods have a good performance on k-space identification, although the k-space obtained by IWC is not as smooth as that extracted from INCOME and AKSI. It is noted that INCOME is a method that can provide wavenumber with high numerical precision under the perfect condition²⁴ but has a poor performance in noisy environments, especially when the signal is disturbed by the high level of signal noise.

### 4.1.2. Comparison of inverse methods under signal noise condition

Fig. 7 shows the wavenumbers extracted by AKSI when considering displacements affected by different signal noise levels as input parameters. The additive noise here is Gaussian white noise with the Signal-to-Noise Ratio (SNR) generated by the AWGN function in MATLAB and it can be measured by SNR. Fig. 7 shows that the dispersion curve obtained by AKSI is reliable when SNR is higher than 20, having a good agreement with the analytical solution.
<figure><img src="image_13.png" /><figcaption>Fig. 7 Dispersion curve obtained by AKSI for isotropic aluminum plate under different levels of signal noise conditions in different frequency ranges.</figcaption></figure>

Figs. 8(a) and (b) show the comparison of the dispersion curves along x direction and the k-space extracted by AKSI, IWC, and INCOME when SNR = 25, respectively. The k-space is formed by estimated wavenumbers along 24 angles of interest at 300 Hz. From Figs. 8(a) and (b), one can observe that AKSI results are in good agreement with the analytical solution even only using 25 measurement points, which indicates that AKSI is a reliable method to extract dispersion curves and k-space in noisy environments. The dispersion curve and k-space extracted by IWC are also acceptable even though the dispersion curve is full of fluctuation at low frequencies. Comparatively, the wavenumbers extracted by INCOME deviate significantly from the analytical wavenumbers at some frequencies in extracted dispersion curve and some angles in extracted k-space. This is mainly because the nature of INCOME is Prony which suffers from a severe ill-conditioned problem when the input parameters are affected by undesirable noise level.

## 4.2. Numerical study 2: An orthotropic plate

Compare the performance of AKSI, INCOME, and IWC in extracting bending wavenumbers of a strongly orthotropic thin plate under the small perturbation and the non-periodic sampling conditions. The structure is a 180 cm × 180 cm orthotropic plate of thickness $h = 10^{-3}$ m, $E_1 = 120$ GPa, $E_2 = 10$ GPa, shear modulus $G_{12} = 4.9$ GPa, $v_{12} = 0.3$, $\rho = 1510$ kg/m$^3$, $\eta_1 = 0.005$, and $\eta_2 = 0.005$. This plate has the simply supported boundary condition and is excited with a point load on the center of the plate. The transverse displacement $U(x, y, \omega)$ of this orthotropic plate is given by Ref. 40 where the same plate has been tested:

<figure><img src="image_14.png" /><figcaption>(a) Dispersion curve comparison</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(b) K-space comparison</figcaption></figure>

<caption>Fig. 8 Comparison of dispersion curve and k-space at 300 Hz obtained by INCOME and AKSI for isotropic aluminum plate under signal noise condition (SNR = 25).</caption>

$$U(x, y, \omega) = \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{F_o \Phi_{nm}(x_0, y_0) \Phi_{nm}(x, y)}{\left\{ \pi^4 \left[ D_1 \left(\frac{n}{z_1}\right)^4 + D_2 \left(\frac{m}{z_2}\right)^4 + 2(D_{12} + D_{66}) \left(\frac{n}{z_1}\right)^2 \left(\frac{m}{z_2}\right)^2 \right] - \rho h \omega^2 \right\} \int_0^{z_1} \int_0^{z_2} \Phi_{nm}^2(x, y) dx dy} \quad (40)$$
<figure><img src="image_13.png" /><figcaption>Fig. 7 Dispersion curve obtained by AKSI for isotropic aluminum plate under different levels of signal noise conditions in different frequency ranges.</figcaption></figure>

Figs. 8(a) and (b) show the comparison of the dispersion curves along x direction and the k-space extracted by AKSI, IWC, and INCOME when SNR = 25, respectively. The k-space is formed by estimated wavenumbers along 24 angles of interest at 300 Hz. From Figs. 8(a) and (b), one can observe that AKSI results are in good agreement with the analytical solution even only using 25 measurement points, which indicates that AKSI is a reliable method to extract dispersion curves and k-space in noisy environments. The dispersion curve and k-space extracted by IWC are also acceptable even though the dispersion curve is full of fluctuation at low frequencies. Comparatively, the wavenumbers extracted by INCOME deviate significantly from the analytical wavenumbers at some frequencies in extracted dispersion curve and some angles in extracted k-space. This is mainly because the nature of INCOME is Prony which suffers from a severe ill-conditioned problem when the input parameters are affected by undesirable noise level.

## 4.2. Numerical study 2: An orthotropic plate

Compare the performance of AKSI, INCOME, and IWC in extracting bending wavenumbers of a strongly orthotropic thin plate under the small perturbation and the non-periodic sampling conditions. The structure is a 180 cm × 180 cm orthotropic plate of thickness $h = 10^{-3}$ m, $E_1 = 120$ GPa, $E_2 = 10$ GPa, shear modulus $G_{12} = 4.9$ GPa, $v_{12} = 0.3$, $\rho = 1510$ kg/m$^3$, $\eta_1 = 0.005$, and $\eta_2 = 0.005$. This plate has the simply supported boundary condition and is excited with a point load on the center of the plate. The transverse displacement $U(x, y, \omega)$ of this orthotropic plate is given by Ref. 40 where the same plate has been tested:

<figure><img src="image_14.png" /><figcaption>(a) Dispersion curve comparison</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(b) K-space comparison</figcaption></figure>

<caption>Fig. 8 Comparison of dispersion curve and k-space at 300 Hz obtained by INCOME and AKSI for isotropic aluminum plate under signal noise condition (SNR = 25).</caption>

$$U(x, y, \omega) = \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{F_o \Phi_{nm}(x_0, y_0) \Phi_{nm}(x, y)}{\left\{ \pi^4 \left[ D_1 \left(\frac{n}{z_1}\right)^4 + D_2 \left(\frac{m}{z_2}\right)^4 + 2(D_{12} + D_{66}) \left(\frac{n}{z_1}\right)^2 \left(\frac{m}{z_2}\right)^2 \right] - \rho h \omega^2 \right\} \int_0^{z_1} \int_0^{z_2} \Phi_{nm}^2(x, y) dx dy} \quad (40)$$
<figure><img src="image_16.png" /><figcaption>Fig. 9 Comparison of perfect displacement field and noisy displacement field of the isotropic plate at 200 Hz.</figcaption></figure>

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

with

$$D_1 = \frac{E_1(1 + i\eta_1)h^3}{12(1 - v_{12}v_{21})} \quad (41)$$

$$D_2 = \frac{E_2(1 + i\eta_2)h^3}{12(1 - v_{12}v_{21})} \quad (42)$$

$$\Phi_{nm}(x, y) = \sin \frac{n\pi x}{z_1} \sin \frac{m\pi y}{z_2} \quad (43)$$

where $D_1$, $D_2$, and $D_{12}$ are the bending stiffness; $D_{66}$ is the torsional stiffness; $v_{12}$ and $v_{21}$ are the major and minor Poisson's ratios. Their relationship is explained in Refs. 41,42:

$$D_3 = D_{12} + D_{66} = v_{12}D_2 + \frac{G_{12}h^3}{6} \quad (44)$$

$$v_{12}E_2 = v_{21}E_1 \quad (45)$$

The transverse displacements on a 901 × 901 mesh grid are calculated from Eq. (40) and used as the input parameters of INCOME. It should be noted that smaller sampling intervals can improve the accuracy of the numerical integration in Eq. (33) and hence the robustness of the AKSI to perturbations. Based on this, the sampling intervals in x and y directions are chosen as 2 mm. The AKSI is carried out using 100 samples at each angle. In Section 4.2.1, the comparison of INCOME and AKSI on wavenumber identification is provided when periodic samples are affected by the small perturbation problem. Section 4.2.2 compares the performance of IWC and AKSI using non-periodic samples as input parameters under perfect conditions.

### 4.2.1. Comparison of inverse methods under small perturbation condition

To create samples in the small perturbation condition, two small perturbation ratios $\xi$ and $\zeta$ are introduced into the coordinates of eigen-modes at each measured point $(x, y)$ in Eq. (41) as

$$\Phi_{nm}(\hat{x}, \hat{y}) = \sin \frac{n\pi\hat{x}}{a} \sin \frac{m\pi\hat{y}}{b} \quad (46)$$

with

$$\left\{ \begin{array}{l} \hat{x}_n = x_n \pm \xi_n \Delta x_n \\ \hat{y}_n = y_n \pm \zeta_n \Delta y_n \end{array} \right. \quad n = 1, 2, \dots, N \quad (47)$$

where $\xi$ and $\zeta$ obey the uniform distribution, and their maximum value is 5. Fig. 9 compares the displacement field at 200 Hz under perfect and small perturbation conditions. It can be seen that the displacement field is perturbed under the small perturbation condition. The effect of the small perturbation problem on transverse displacements can also be

<figure><img src="image_18.png" /><figcaption>Fig. 10 Comparison of displacement curves at 200 Hz along 0° direction under perfect and small perturbation conditions.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. 11 Dispersion curves at five different directions obtained by AKSI under small perturbation condition.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>Fig. 9 Comparison of perfect displacement field and noisy displacement field of the isotropic plate at 200 Hz.</figcaption></figure>

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

with

$$D_1 = \frac{E_1(1 + i\eta_1)h^3}{12(1 - v_{12}v_{21})} \qquad (41)$$

$$D_2 = \frac{E_2(1 + i\eta_2)h^3}{12(1 - v_{12}v_{21})} \qquad (42)$$

$$\Phi_{nm}(x, y) = \sin \frac{n\pi x}{z_1} \sin \frac{m\pi y}{z_2} \qquad (43)$$

where $D_1$, $D_2$, and $D_{12}$ are the bending stiffness; $D_{66}$ is the torsional stiffness; $v_{12}$ and $v_{21}$ are the major and minor Poisson's ratios. Their relationship is explained in Refs. 41,42:

$$D_3 = D_{12} + D_{66} = v_{12}D_2 + \frac{G_{12}h^3}{6} \qquad (44)$$

$$v_{12}E_2 = v_{21}E_1 \qquad (45)$$

The transverse displacements on a 901 × 901 mesh grid are calculated from Eq. (40) and used as the input parameters of INCOME. It should be noted that smaller sampling intervals can improve the accuracy of the numerical integration in Eq. (33) and hence the robustness of the AKSI to perturbations. Based on this, the sampling intervals in x and y directions are chosen as 2 mm. The AKSI is carried out using 100 samples at each angle. In Section 4.2.1, the comparison of INCOME and AKSI on wavenumber identification is provided when periodic samples are affected by the small perturbation problem. Section 4.2.2 compares the performance of IWC and AKSI using non-periodic samples as input parameters under perfect conditions.

### 4.2.1. Comparison of inverse methods under small perturbation condition

To create samples in the small perturbation condition, two small perturbation ratios $\xi$ and $\zeta$ are introduced into the coordinates of eigen-modes at each measured point $(x, y)$ in Eq. (41) as

$$\Phi_{nm}(\hat{x}, \hat{y}) = \sin \frac{n\pi\hat{x}}{a} \sin \frac{m\pi\hat{y}}{b} \qquad (46)$$

with

$$\begin{cases} \hat{x}_n = x_n \pm \xi_n \Delta x_n \\ \hat{y}_n = y_n \pm \zeta_n \Delta y_n \end{cases} \qquad n = 1, 2, \dots, N \qquad (47)$$

where $\xi$ and $\zeta$ obey the uniform distribution, and their maximum value is 5. Fig. 9 compares the displacement field at 200 Hz under perfect and small perturbation conditions. It can be seen that the displacement field is perturbed under the small perturbation condition. The effect of the small perturbation problem on transverse displacements can also be

<figure><img src="image_18.png" /><figcaption>Fig. 10 Comparison of displacement curves at 200 Hz along 0° direction under perfect and small perturbation conditions.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. 11 Dispersion curves at five different directions obtained by AKSI under small perturbation condition.</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>**Fig. 12** Comparison of dispersion curve at 0° direction and k-space at 200 Hz obtained by INCOME and AKSI for orthotropic plate under small perturbation condition.</figcaption></figure>

observed by the displacement curve at 0° direction, as shown in Fig. 10. Under this condition, the dispersion curves obtained by AKSI at different angles, including $\theta_1, \theta_3, \theta_4, \theta_5$, and $\theta_7$ of Fig. 5, are plotted in Fig. 11. Fig. 11 shows that the acceptable dispersion curve can be extracted from AKSI under perturbation conditions. Figs. 12(a) and (b) present the comparison of the dispersion curve along x direction and the k-space at 200 Hz obtained by IWC, INCOME, and AKSI, respectively. Figs. 12(a) and (b) show that on the one hand, AKSI is a good candidate for identifying the dispersion curve and k-space under small perturbation conditions due to its remarkable filtering capabilities. On the other hand, the profile of the

<figure><img src="image_21.png" /><figcaption>**Fig. 13** Comparison of dispersion curves obtained by AKSI and IWC based on 20 non-periodic displacements at 90° direction under perfect condition.</figcaption></figure>

k-space from AKSI illustrates that the tested plate has orthotropic dynamic behavior. For IWC, the extracted dispersion curve at low frequencies is full of fluctuations and the dispersion curve becomes smoother with frequency increasing. The expected k-space can also be used to approximately identify the orthotropic dynamic behavior of the plate. In contrast, the results of the dispersion curve and k-space extracted by INCOME in Figs. 12(a) and (b) illustrate that this method is not able to extract reliable wavenumber under small perturbation conditions due to the nature of the Prony method.

### 4.2.2. Comparison of inverse methods under non-periodic sampling condition

To compare the performance of AKSI and IWC on the dispersion curve identification under the non-periodic sampling condition, 20 non-periodic displacements are created based on the 36 periodic unperturbed displacements with the sampling interval of 1 cm at 90° direction. Fig. 13(a) shows the corresponding displacement curves at 100 Hz. The resulting dispersion curves are shown in Fig. 13(b), indicating that AKSI and IWC are not limited to periodic sampling. Furthermore, as was described in Section 4.1.1, the dispersion curve extracted from AKSI is more accurate than that extracted from IWC in the low-frequency range. The ability of AKSI to be free from periodic sampling comes from the fact that AKSI treats the signal as a continuous signal function.

### 4.3. Numerical study 3: A periodic LT plate

To study the sensitivity of AKSI, IWC, and BWI to the structural periodicity, they are applied to extract the dispersion curve of a periodic LT plate under two stochastic conditions: (A) Structural periodicity is precisely known or approximately known; (B) Structural periodicity is unknown. These two stochastic cases are studied in Sections 4.3.1 and 4.3.2, respectively. The structure is modeled by the solid element of COMSOL Multiphysics. A boundary force excitation is applied on the left side of the periodic structure along x direction, and the boundary condition is set to be free-free. The overall model is made of 5 unit cells and has a total of 742548 degrees of freedom. Fig. 14 shows the structural geometry, the composition of the unit cell, the distribution of measurement points, and the displacement field at 700 Hz. It can be seen that unit cell is made of two parts: steel part ($E = 200$ GPa, $v = 0.3$, $\rho = 7800$ kg/m$^3$, $\eta = 0.0001$) and nylon part ($E = 2$ GPa,

<figure><img src="image_22.png" /></figure>
<figure><img src="image_20.png" /><figcaption>**Fig. 12** Comparison of dispersion curve at 0° direction and k-space at 200 Hz obtained by INCOME and AKSI for orthotropic plate under small perturbation condition.</figcaption></figure>

observed by the displacement curve at 0° direction, as shown in Fig. 10. Under this condition, the dispersion curves obtained by AKSI at different angles, including $θ_1, θ_3, θ_4, θ_5$, and $θ_7$ of Fig. 5, are plotted in Fig. 11. Fig. 11 shows that the acceptable dispersion curve can be extracted from AKSI under perturbation conditions. Figs. 12(a) and (b) present the comparison of the dispersion curve along x direction and the k-space at 200 Hz obtained by IWC, INCOME, and AKSI, respectively. Figs. 12(a) and (b) show that on the one hand, AKSI is a good candidate for identifying the dispersion curve and k-space under small perturbation conditions due to its remarkable filtering capabilities. On the other hand, the profile of the

<figure><img src="image_21.png" /><figcaption>**Fig. 13** Comparison of dispersion curves obtained by AKSI and IWC based on 20 non-periodic displacements at 90° direction under perfect condition.</figcaption></figure>

k-space from AKSI illustrates that the tested plate has orthotropic dynamic behavior. For IWC, the extracted dispersion curve at low frequencies is full of fluctuations and the dispersion curve becomes smoother with frequency increasing. The expected k-space can also be used to approximately identify the orthotropic dynamic behavior of the plate. In contrast, the results of the dispersion curve and k-space extracted by INCOME in Figs. 12(a) and (b) illustrate that this method is not able to extract reliable wavenumber under small perturbation conditions due to the nature of the Prony method.

### 4.2.2. Comparison of inverse methods under non-periodic sampling condition

To compare the performance of AKSI and IWC on the dispersion curve identification under the non-periodic sampling condition, 20 non-periodic displacements are created based on the 36 periodic unperturbed displacements with the sampling interval of 1 cm at 90° direction. Fig. 13(a) shows the corresponding displacement curves at 100 Hz. The resulting dispersion curves are shown in Fig. 13(b), indicating that AKSI and IWC are not limited to periodic sampling. Furthermore, as was described in Section 4.1.1, the dispersion curve extracted from AKSI is more accurate than that extracted from IWC in the low-frequency range. The ability of AKSI to be free from periodic sampling comes from the fact that AKSI treats the signal as a continuous signal function.

### 4.3. Numerical study 3: A periodic LT plate

To study the sensitivity of AKSI, IWC, and BWI to the structural periodicity, they are applied to extract the dispersion curve of a periodic LT plate under two stochastic conditions: (A) Structural periodicity is precisely known or approximately known; (B) Structural periodicity is unknown. These two stochastic cases are studied in Sections 4.3.1 and 4.3.2, respectively. The structure is modeled by the solid element of COMSOL Multiphysics. A boundary force excitation is applied on the left side of the periodic structure along x direction, and the boundary condition is set to be free-free. The overall model is made of 5 unit cells and has a total of 742548 degrees of freedom. Fig. 14 shows the structural geometry, the composition of the unit cell, the distribution of measurement points, and the displacement field at 700 Hz. It can be seen that unit cell is made of two parts: steel part ($E = 200$ GPa, $ν = 0.3$, $ρ = 7800$ kg/m$^3$, $η = 0.0001$) and nylon part ($E = 2$ GPa,

<figure><img src="image_22.png" /></figure>
<figure><img src="image_23.png" /><figcaption>Fig. 14 5-UC LT periodic structure with a 2D grid and its displacement field at 700 Hz.</figcaption></figure>

v = 0.33, ρ = 1150 kg/m³, η = 0.0001). Both parts have the same lengths of 5 cm and the same width of 10 cm, while the thickness of the steel and nylon parts are 3 mm and 2 mm. The simulated transverse displacements are measured on a regular 50 × 11 grid at the frequency range from 10 Hz to 1500 Hz with the step of 10 Hz, which is the input parameter of BWI and IWC, while the displacements on the midline of the plate along x direction (as shown by the blue dot in Fig. 14(a)) are used for AKSI. In addition, Figs. 14(c) and (d) presents the displacement field at 700 Hz and the corresponding displacement curve on the midline of the plate along x direction.

### 4.3.1. Structural periodicity is precisely known or approximately known

For the following extracted results, we use $k_{AWI}$, $k_{BWI}$, $k_{IWC}$, and $k_{WFEM}$ as the wavenumber from AWI, BWI, IWC, and WFEM, respectively. To extract the dispersion curve of periodic structures and identify band gaps accurately, BWI requires that the periodicity of the sampling interval must be an integer multiple of the structural periodicity. As shown in Fig. 15, BWI can provide an accurate dispersion curve when the periodicity of the structure is precisely known under the perfect condition. From extracted dispersion curve, one can clearly distinguish three band gaps displayed as green patches in Fig. 15. The band gap refers to the frequency range where the wave cannot propagate in the medium. The displacement

<figure><img src="image_24.png" /><figcaption>Fig. 15 Complex wavenumbers obtained by BWI for LT periodic structure when structural periodicity (length of a unit cell) is exactly known.</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>Fig. 16 Complex wavenumbers obtained by BWI for LT periodic structure when structural periodicity (length of a unit cell) is approximately known.</figcaption></figure>

field at 700 Hz shown in Fig. 14(c) can also illustrate the effect of the band gap on wave attenuation. Furthermore, the wavenumbers extracted by BWI are restricted to the first Brillouin zone, where the wavelengths are longer than the length of the unit cell. Such wavenumbers are known as Bloch wavenumbers of periodic structures.

However, due to structural deformation or damage, the structural periodicity cannot be determined precisely in reality.

<figure><img src="image_26.png" /><figcaption>Fig. 17 Comparison of complex wavenumbers obtained by IWC, BWI, and AKSI for LT periodic structure when structural periodicity (length of a unit cell) is unknown.</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>Fig. 14 5-UC LT periodic structure with a 2D grid and its displacement field at 700 Hz.</figcaption></figure>

v = 0.33, ρ = 1150 kg/m³, η = 0.0001). Both parts have the same lengths of 5 cm and the same width of 10 cm, while the thickness of the steel and nylon parts are 3 mm and 2 mm. The simulated transverse displacements are measured on a regular 50 × 11 grid at the frequency range from 10 Hz to 1500 Hz with the step of 10 Hz, which is the input parameter of BWI and IWC, while the displacements on the midline of the plate along x direction (as shown by the blue dot in Fig. 14(a)) are used for AKSI. In addition, Figs. 14(c) and (d) presents the displacement field at 700 Hz and the corresponding displacement curve on the midline of the plate along x direction.

### 4.3.1. Structural periodicity is precisely known or approximately known

For the following extracted results, we use $k_{AWI}$, $k_{BWI}$, $k_{IWC}$, and $k_{WFEM}$ as the wavenumber from AWI, BWI, IWC, and WFEM, respectively. To extract the dispersion curve of periodic structures and identify band gaps accurately, BWI requires that the periodicity of the sampling interval must be an integer multiple of the structural periodicity. As shown in Fig. 15, BWI can provide an accurate dispersion curve when the periodicity of the structure is precisely known under the perfect condition. From extracted dispersion curve, one can clearly distinguish three band gaps displayed as green patches in Fig. 15. The band gap refers to the frequency range where the wave cannot propagate in the medium. The displacement

<figure><img src="image_24.png" /><figcaption>Fig. 15 Complex wavenumbers obtained by BWI for LT periodic structure when structural periodicity (length of a unit cell) is exactly known.</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>Fig. 16 Complex wavenumbers obtained by BWI for LT periodic structure when structural periodicity (length of a unit cell) is approximately known.</figcaption></figure>

field at 700 Hz shown in Fig. 14(c) can also illustrate the effect of the band gap on wave attenuation. Furthermore, the wavenumbers extracted by BWI are restricted to the first Brillouin zone, where the wavelengths are longer than the length of the unit cell. Such wavenumbers are known as Bloch wavenumbers of periodic structures.

However, due to structural deformation or damage, the structural periodicity cannot be determined precisely in reality.

<figure><img src="image_26.png" /><figcaption>Fig. 17 Comparison of complex wavenumbers obtained by IWC, BWI, and AKSI for LT periodic structure when structural periodicity (length of a unit cell) is unknown.</figcaption></figure>
Considering an approximate unit-cell length (9 cm), the corre-
sponding dispersion curve obtained by BWI is shown in
Fig. 16. The dispersion curve shows significant scatter, result-
ing in the incorrect identification of band gaps.

4.3.2. *Structural periodicity is unknown*

In some cases, the periodicity of periodic structures may be
completely unknown, such as when the periodic structures
are covered by other materials. To investigate this condition,
the sampling interval is set to 1 cm instead of 10 cm (the struc-
tural periodicity).

Fig. 17 shows the comparison of dispersion curves extracted by IWC, BWI, and AKSI. Two phenomena can be observed from Fig. 17: (A) The location and width of the band gaps can be identified by the dispersion curves extracted by AKSI and IWC under this stochastic condition. This is because these two methods have no sampling requirement, and their performance is not related to the periodicity of periodic structures. Comparatively, the BWI requires that the periodicity of the sampling interval must be an integer multiple of the structural periodicity. Therefore, it is not able to provide the accurate wavenumber when structural periodicity is unknown. (B) The dispersion curve identified by AKSI is limited to the first Brillouin zone, where the wavelengths are longer than one period of the structure. That means that the real part of the wavenumber is always within [0, 31.4 rad/m] (31.4 rad/m is calculated by $\pi/\Delta$, where $\Delta$ is the length of the unit cell which is 10 cm in this case). This is because the signal model in AKSI is defined as a sum of exponential functions, as shown in Eq. (8), which is also the displacement model of the Bloch wave theorem<sup>43</sup> for periodic structures. Therefore, the wavenumber obtained from AKSI can agree well with that of WFEM, having the characteristics of the Bloch wavenumber. Comparatively, IWC considers the signal model as an inhomogeneous wave propagating in the whole structure without the characteristics of the Bloch wavenumber, resulting in the dispersion curve that can be bigger than 31.4 rad/m. It is worth noting that both formats of dispersion

<figure><img src="image_27.png" /><figcaption>Fig. 18 Experimental set-up of a thin isotropic steel plate.</figcaption></figure>

<figure><img src="image_28.png" /><figcaption>(a) Envelope of coherence calculated at all measured points</figcaption></figure>

<figure><img src="image_29.png" /><figcaption>(b) FRF for each measured point along 0° at 1190.625 Hz</figcaption></figure>

<figure><img src="image_30.png" /><figcaption>(a) Dispersion curve at 0° direction</figcaption></figure>

<figure><img src="image_31.png" /><figcaption>(b) k-space</figcaption></figure>

<caption>**Fig. 20 Comparison of dispersion curve at 0° direction and k-pace at 534.375 Hz obtained by IWC and AKSI for isotropic steel plate using experimental data.**</caption>
Considering an approximate unit-cell length (9 cm), the corre-
sponding dispersion curve obtained by BWI is shown in
Fig. 16. The dispersion curve shows significant scatter, result-
ing in the incorrect identification of band gaps.

4.3.2. *Structural periodicity is unknown*

In some cases, the periodicity of periodic structures may be
completely unknown, such as when the periodic structures
are covered by other materials. To investigate this condition,
the sampling interval is set to 1 cm instead of 10 cm (the struc-
tural periodicity).

Fig. 17 shows the comparison of dispersion curves extracted by IWC, BWI, and AKSI. Two phenomena can be observed from Fig. 17: (A) The location and width of the band gaps can be identified by the dispersion curves extracted by AKSI and IWC under this stochastic condition. This is because these two methods have no sampling requirement, and their performance is not related to the periodicity of periodic structures. Comparatively, the BWI requires that the periodicity of the sampling interval must be an integer multiple of the structural periodicity. Therefore, it is not able to provide the accurate wavenumber when structural periodicity is unknown. (B) The dispersion curve identified by AKSI is limited to the first Brillouin zone, where the wavelengths are longer than one period of the structure. That means that the real part of the wavenumber is always within [0, 31.4 rad/m] (31.4 rad/m is calculated by $\pi/\Delta$, where $\Delta$ is the length of the unit cell which is 10 cm in this case). This is because the signal model in AKSI is defined as a sum of exponential functions, as shown in Eq. (8), which is also the displacement model of the Bloch wave theorem<sup>43</sup> for periodic structures. Therefore, the wavenumber obtained from AKSI can agree well with that of WFEM, having the characteristics of the Bloch wavenumber. Comparatively, IWC considers the signal model as an inhomogeneous wave propagating in the whole structure without the characteristics of the Bloch wavenumber, resulting in the dispersion curve that can be bigger than 31.4 rad/m. It is worth noting that both formats of dispersion

<figure><img src="image_27.png" /><figcaption>**Fig. 18** Experimental set-up of a thin isotropic steel plate.</figcaption></figure>

<figure><img src="image_28.png" /><figcaption>(a) Envelope of coherence calculated at all measured points</figcaption></figure>

<figure><img src="image_29.png" /><figcaption>(b) FRF for each measured point along 0° at 1190.625 Hz</figcaption></figure>

<figure><img src="image_30.png" /><figcaption>(a) Dispersion curve at 0° direction</figcaption></figure>

<figure><img src="image_31.png" /><figcaption>(b) k-space</figcaption></figure>

<caption>**Fig. 20** Comparison of dispersion curve at 0° direction and k-pace at 534.375 Hz obtained by IWC and AKSI for isotropic steel plate using experimental data.</caption>
curves can be used to describe dispersion characteristics of periodic structures and identify the locations of band gaps.

## 5. Experimental study

Apply AKSI and IWC to extract the dispersion curve of a thin isotropic steel plate based on experimental data. The plate has dimensions of 0.95 m × 0.60 m, and its thickness is 2 mm. The values of mechanical properties are $E = 210$ GPa, $v = 0.3$, $ρ = 7850$ kg/m³. Fig. 18 presents the experimental setup where the free-free boundary condition was used by suspending the plate at a fixed frame, and the plate was excited by a point force at the red point using an electrodynamic shaker. AKSI and IWC are performed on the out-of-plane displacements along the white dotted lines in Fig. 18. In each line, 50 periodic FRFs with a sampling interval of 0.44 cm were measured via PSV-400 scanning vibrometer. The measured data are located far away from the excitation and boundary conditions, avoiding the influence of evanescent waves. With a pseudo-random excitation, measurements were conducted in the frequency domain (resolution of 3.125 Hz, 50 averages).

Fig. 19(a) shows the envelope of coherence calculated at all measurement points. The smaller the value of the coherence function, the more the signal is affected by signal noise. For example, the FRFs curve in the 0° direction at 1190.625 Hz, as shown in Fig. 19(b), presents the samples are disturbed by relatively large signal noise. Figs. 20(a) and (b) present the comparison of the dispersion curves extracted by AKSI and IWC, and the k-space at 534.375 Hz along seven considered directions, respectively. As a widely used nonlinear inverse method, IWC has been validated to be robust to signal noise in many previous works.¹⁷,¹⁸ From Fig. 20, the dispersion curve and k-space extracted by AKSI agree well with the analytical solutions, indicating that AKSI can be comparable with IWC in the dispersion curve and k-space identification under noisy environments. This is mainly because the introduction of integrals of AKSI can improve its robustness to signal noise. It is worth noting that AKSI only requires solving several linear problems, which can reduce the computational cost compared to IWC.

## 6. Conclusions

The Algebraic K-Space Identification (AKSI) technique is proposed in the algebraic estimator framework. The main contribution of the proposed method is to achieve dispersion curve and k-space identification of 2D structures under a series of stochastic conditions, making it possible to study the dynamic behavior of 2D structures in realistic conditions. The main conclusions are summarized as follows:

(1) The proposed method draws inspiration from the algebraic estimator to overcome some main issues of other popular inverse methods. It can extract the wavenumbers from a more stochastic signal model, which considers the influence of external signal noise and the small perturbation on the structural responses. The small perturbation is reflected in the geometric variability of measuring points' coordinates, which can be caused by various uncertainties, such as grid distortion, structural uncertainty, and operational error.

(2) This paper benchmarks the proposed method and three other representative inverse methods, including IWC, INCOME, BWI, on three numerical applications under different stochastic conditions. These numerical cases involve the extraction of the k-space and dispersion curve of an isotropic plate under signal noise conditions, the extraction of the k-space and dispersion curve extraction of an orthotropic plate under small perturbation and non-periodic sampling conditions, and the band gap identification of an LT-periodic structure under the condition of unknown structural periodicity. These numerical investigations show: (A) Compared to INCOME, AKSI performs better under perturbations conditions such as signal noise and small perturbations. In addition, AKSI is not limited to periodic sampling; (B) Compared to BWI, the performance of AKSI is not related to the structural periodicity, significantly improving its application to periodic structures under stochastic conditions; (C) Compared to IWC, AKSI has a better performance in the low-frequency range. In addition, IWC is a nonlinear iterative method with a high computational cost, whereas AKSI only requires solving several linear problems, reducing the computational cost.

(3) The proposed method was experimentally used to extract the dispersion curve and k-space of a large-scale thin steel plate under a noisy environment. The high level of signal noise can be reflected by the low correlation function values. This experimental case can further validate the robustness of the proposed method to signal noise. In addition, the comprehensive comparison between the popular inverse methods under different stochastic conditions can be a benchmark for developing wavenumber identification methods in the future.

Up to now, AKSI can extract the imaginary part of the wavenumbers of the plate only when the shape of the displacement field is close to a plane wave. This is mainly because the mathematical model of AKSI is a plane wave model, which is defined as a sum of exponential functions. This is also a limitation of inverse methods based on the plane wave model, and a deeper investigation into how to overcome this limitation is desirable in future work.

## 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.

## Acknowledgements

This study was supported by the Lyon Acoustics Center of Lyon University, France. The research of X. Li is funded by the China Scholarship Council (CSC).

## References

1. Honda S, Narita Y. Vibration design of laminated fibrous composite plates with local anisotropy induced by short fibers and curvilinear fibers. *Compos Struct* 2011;93(2):902–10.
curves can be used to describe dispersion characteristics of periodic structures and identify the locations of band gaps.

## 5. Experimental study

Apply AKSI and IWC to extract the dispersion curve of a thin isotropic steel plate based on experimental data. The plate has dimensions of 0.95 m × 0.60 m, and its thickness is 2 mm. The values of mechanical properties are $E = 210 \text{ GPa}, \nu = 0.3, \rho = 7850 \text{ kg/m}^3$. Fig. 18 presents the experimental setup where the free-free boundary condition was used by suspending the plate at a fixed frame, and the plate was excited by a point force at the red point using an electrodynamic shaker. AKSI and IWC are performed on the out-of-plane displacements along the white dotted lines in Fig. 18. In each line, 50 periodic FRFs with a sampling interval of 0.44 cm were measured via PSV-400 scanning vibrometer. The measured data are located far away from the excitation and boundary conditions, avoiding the influence of evanescent waves. With a pseudo-random excitation, measurements were conducted in the frequency domain (resolution of 3.125 Hz, 50 averages).

Fig. 19(a) shows the envelope of coherence calculated at all measurement points. The smaller the value of the coherence function, the more the signal is affected by signal noise. For example, the FRFs curve in the 0° direction at 1190.625 Hz, as shown in Fig. 19(b), presents the samples are disturbed by relatively large signal noise. Figs. 20(a) and (b) present the comparison of the dispersion curves extracted by AKSI and IWC, and the k-space at 534.375 Hz along seven considered directions, respectively. As a widely used nonlinear inverse method, IWC has been validated to be robust to signal noise in many previous works.¹⁷,¹⁸ From Fig. 20, the dispersion curve and k-space extracted by AKSI agree well with the analytical solutions, indicating that AKSI can be comparable with IWC in the dispersion curve and k-space identification under noisy environments. This is mainly because the introduction of integrals of AKSI can improve its robustness to signal noise. It is worth noting that AKSI only requires solving several linear problems, which can reduce the computational cost compared to IWC.

## 6. Conclusions

The Algebraic K-Space Identification (AKSI) technique is proposed in the algebraic estimator framework. The main contribution of the proposed method is to achieve dispersion curve and k-space identification of 2D structures under a series of stochastic conditions, making it possible to study the dynamic behavior of 2D structures in realistic conditions. The main conclusions are summarized as follows:

(1) The proposed method draws inspiration from the algebraic estimator to overcome some main issues of other popular inverse methods. It can extract the wavenumbers from a more stochastic signal model, which considers the influence of external signal noise and the small perturbation on the structural responses. The small perturbation is reflected in the geometric variability of measuring points' coordinates, which can be caused by various uncertainties, such as grid distortion, structural uncertainty, and operational error.

(2) This paper benchmarks the proposed method and three other representative inverse methods, including IWC, INCOME, BWI, on three numerical applications under different stochastic conditions. These numerical cases involve the extraction of the k-space and dispersion curve of an isotropic plate under signal noise conditions, the extraction of the k-space and dispersion curve extraction of an orthotropic plate under small perturbation and non-periodic sampling conditions, and the band gap identification of an LT-periodic structure under the condition of unknown structural periodicity. These numerical investigations show: (A) Compared to INCOME, AKSI performs better under perturbations conditions such as signal noise and small perturbations. In addition, AKSI is not limited to periodic sampling; (B) Compared to BWI, the performance of AKSI is not related to the structural periodicity, significantly improving its application to periodic structures under stochastic conditions; (C) Compared to IWC, AKSI has a better performance in the low-frequency range. In addition, IWC is a nonlinear iterative method with a high computational cost, whereas AKSI only requires solving several linear problems, reducing the computational cost.

(3) The proposed method was experimentally used to extract the dispersion curve and k-space of a large-scale thin steel plate under a noisy environment. The high level of signal noise can be reflected by the low correlation function values. This experimental case can further validate the robustness of the proposed method to signal noise. In addition, the comprehensive comparison between the popular inverse methods under different stochastic conditions can be a benchmark for developing wavenumber identification methods in the future.

Up to now, AKSI can extract the imaginary part of the wavenumbers of the plate only when the shape of the displacement field is close to a plane wave. This is mainly because the mathematical model of AKSI is a plane wave model, which is defined as a sum of exponential functions. This is also a limitation of inverse methods based on the plane wave model, and a deeper investigation into how to overcome this limitation is desirable in future work.

## 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.

## Acknowledgements

This study was supported by the Lyon Acoustics Center of Lyon University, France. The research of X. Li is funded by the China Scholarship Council (CSC).

## References

1. Honda S, Narita Y. Vibration design of laminated fibrous composite plates with local anisotropy induced by short fibers and curvilinear fibers. *Compos Struct* 2011;93(2):902–10.
2. Zhang J, Reynders E, De Roeck G, et al. Model updating of periodic structures based on free wave characteristics. *J Sound Vib* 2019;442:281–307.

3. Zhao HJ, Feng Y, Li W, et al. Numerical study and topology optimization of vibration isolation support structures. *Int J Mech Sci* 2022;228:107507.

4. Wang Z, Luo QT, Li Q, et al. Design optimization of bioinspired helicoidal CFRPP/GFRPP hybrid composites for multiple low-velocity impact loads. *Int J Mech Sci* 2022;219:107064.

5. Tufano G, Errico F, Robin O, et al. K-space analysis of complex large-scale meta-structures using the inhomogeneous wave correlation method. *Mech Syst Signal Process* 2020;135:106407.

6. Sun WJ, Thompson D, Toward M, et al. Modelling of vibration and noise behaviour of embedded tram tracks using a wavenumber domain method. *J Sound Vib* 2020;481:115446.

7. Wang XP, Cai J, Zhou ZQ. A Lamb wave signal reconstruction method for high-resolution damage imaging. *Chin J Aeronaut* 2019;32(5):1087–99.

8. Lin RM, Ng TY. Applications of higher-order frequency response functions to the detection and damage assessment of general structural systems with breathing cracks. *Int J Mech Sci* 2018;148:652–66.

9. Denis V, Mencik JM. A wave-based optimization approach of curved joints for improved defect detection in waveguide assemblies. *J Sound Vib* 2020;465:115003.

10. Tian ZH, Xiao WF, Ma ZY, et al. Dispersion curve regression-assisted wideband local wavenumber analysis for characterizing three-dimensional (3D) profile of hidden corrosion damage. *Mech Syst Signal Process* 2021;150:107347.

11. Hamidpour M, Nami MR, Khosravifard A. An effective crack identification method in viscoelastic media using an inverse meshfree method. *Int J Mech Sci* 2021;212:106834.

12. Ichchou MN, Bareille O, Berthaut J. Identification of effective sandwich structural properties via an inverse wave approach. *Eng Struct* 2008;30(10):2591–604.

13. Santoni A, Schoenwald S, Van Damme B, et al. Determination of the elastic and stiffness characteristics of cross-laminated timber plates from flexural wave velocity measurements. *J Sound Vib* 2017;400:387–401.

14. Droz C, Bareille O, Ichchou MN. A new procedure for the determination of structural characteristics of sandwich plates in medium frequencies. *Compos B* 2017;112:103–11.

15. Berthaut J, Ichchou MN, Jezequel L. K-space identification of apparent structural behaviour. *J Sound Vib* 2005;280(3–5):1125–31.

16. Tufano G. K-space analysis of complex large-scale periodic structures [dissertation]. Ecully: Lyon University and Leuven, KU Leuven; 2020.

17. Lajili R, Bareille O, Bouazizi M, et al. Composite beam identification using a variant of the inhomogeneous wave correlation method in presence of uncertainties. *Eng Comput* 2018;35(6):2126–64.

18. Lajili R, Chikhaoui K, Zergoune Z, et al. Impact of the vibration measurement points geometric coordinates uncertainties on two-dimensional k-space identification: application to a sandwich plate with honeycomb core. *Mech Syst Signal Process* 2022;167:108509.

19. Geslain A, Raetz S, Hiraiwa M, et al. Spatial Laplace transform for complex wavenumber recovery and its application to the analysis of attenuation in acoustic systems. *J Appl Phys* 2016;120(13):135107.

20. Yan GQ, Raetz S, Groby JP, et al. Estimation via laser ultrasonics of the ultrasonic attenuation in a polycrystalline aluminum thin plate using complex wavenumber recovery in the vicinity of a zero-group-velocity lamb mode. *Appl Sci* 2021;11(15):6924.

21. Roozen NB, Labelle L, Leclère Q, et al. Non-contact experimental assessment of apparent dynamic stiffness of constrained-layer damping sandwich plates in a broad frequency range using a Nd: YAG pump laser and a laser Doppler vibrometer. *J Sound Vib* 2017;395:90–101.

22. Morandi F, Robin O, Barbaresi L, et al. Benchmarking of methods for the identification of flexural wavenumbers in wooden plates. *Proceedings of the 23rd international congress on acoustics*; 2019. p. 1294–301.

23. Margerit P, Lebée A, Caron JF, et al. The high-resolution wavevector analysis for the characterization of the dynamic response of composite plates. *J Sound Vib* 2019;458:177–96.

24. Boukadia RF, Claeys C, Droz C, et al. An inverse convolution method for wavenumber extraction (INCOME): formulations and applications. *J Sound Vib* 2022;520:116586.

25. Ribeiro LHMS, Dal Poggetto VF, Huallpa BN, et al. Bloch wavenumber identification of periodic structures using Prony's method. *Mech Syst Signal Process* 2022;178:109242.

26. Li XF, Ichchou M, Zine A, et al. An algebraic wavenumber identification (AWI) technique under stochastic conditions. *Mech Syst Signal Process* 2023;188:109983.

27. Li XF, Ichchou M, Zine A, et al. Wavenumber identification of 1D complex structures using algebraic wavenumber identification (AWI) technique under complex conditions. *J Sound Vib* 2023;548:117524.

28. Fliess M, Sira-Ramírez H. An algebraic framework for linear identification. *ESAIM Control Optim Calc Var* 2003;9:151–68.

29. Allemang RJ. The modal assurance criterion—Twenty years of use and abuse. *Sound Vib* 2003;37(8):14–23.

30. Fliess M, Mboup M, Mounier H, et al. Questioning some paradigms of signal processing via concrete examples. *Algebr Meth Flatness Signal Process State Estim*. 2003. p. 1–21.

31. Beltran-Carbajal F, Silva-Navarro G. On the algebraic parameter identification of vibrating mechanical systems. *Int J Mech Sci* 2015;92:178–86.

32. Stoica P, Selen Y. Model-order selection: a review of information criterion rules. *IEEE Signal Process Mag* 2004;21(4):36–47.

33. Mariani A, Giorgetti A, Chiani M. Model order selection based on information theoretic criteria: design of the penalty. *IEEE Trans Signal Process* 2015;63(11):2779–89.

34. Badeau R, David B, Richard G. Selecting the modeling order for the ESPRIT high resolution method: An alternative approach. In: *2004 IEEE international conference on acoustics, speech, and signal processing*. New York: IEEE; 2004. p. ii–1025.

35. Papy JM, De Lathauwer L, Van Huffel S. A shift invariance-based order-selection technique for exponential data modelling. *IEEE Signal Process Lett* 2007;14(7):473–6.

36. Albert R, Galarza CG. Model order selection for sum of complex exponentials. *2021 IEEE URUCON*. New York: IEEE; 2021. p. 561–5.

37. Margerit P, Lebée A, Caron JF, et al. High resolution wavenumber analysis (HRWA) for the mechanical characterisation of viscoelastic beams. *J Sound Vib* 2018;433:198–211.

38. Okumura S, Nguyen VH, Taki H, et al. Rapid high-resolution wavenumber extraction from ultrasonic guided waves using adaptive array signal processing. *Appl Sci* 2018;8(4):652.

39. Wang C, Reddy J, Lee K Shear deformable beams and plates: *Relationships with classical solutions*. Amsterdam: Elsevier Science; 2000. p. 39–54.

40. Caillet J, Carmona JC, Mazzoni D. Estimation of plate elastic moduli through vibration testing. *Appl Acoust* 2007;68(3):334–49.

41. Muthurajan KG, Sanakaranarayanasamy K, Rao BN. Evaluation of elastic constants of specially orthotropic plates through vibration testing. *J Sound Vib* 2004;272(1–2):413–24.

42. Au FTK, Wang MF. Sound radiation from forced vibration of rectangular orthotropic plates under moving loads. *J Sound Vib* 2005;281(3–5):1057–75.

43. Andreassen E, Jensen JS. Analysis of phononic bandgap structures with dissipation. *J Vib Acoust* 2013;135(4):041015.
2. Zhang J, Reynders E, De Roeck G, et al. Model updating of periodic structures based on free wave characteristics. *J Sound Vib* 2019;442:281–307.

3. Zhao HJ, Feng Y, Li W, et al. Numerical study and topology optimization of vibration isolation support structures. *Int J Mech Sci* 2022;228:107507.

4. Wang Z, Luo QT, Li Q, et al. Design optimization of bioinspired helicoidal CFRPP/GFRPP hybrid composites for multiple low-velocity impact loads. *Int J Mech Sci* 2022;219:107064.

5. Tufano G, Errico F, Robin O, et al. K-space analysis of complex large-scale meta-structures using the inhomogeneous wave correlation method. *Mech Syst Signal Process* 2020;135:106407.

6. Sun WJ, Thompson D, Toward M, et al. Modelling of vibration and noise behaviour of embedded tram tracks using a wavenumber domain method. *J Sound Vib* 2020;481:115446.

7. Wang XP, Cai J, Zhou ZQ. A Lamb wave signal reconstruction method for high-resolution damage imaging. *Chin J Aeronaut* 2019;32(5):1087–99.

8. Lin RM, Ng TY. Applications of higher-order frequency response functions to the detection and damage assessment of general structural systems with breathing cracks. *Int J Mech Sci* 2018;148:652–66.

9. Denis V, Mencik JM. A wave-based optimization approach of curved joints for improved defect detection in waveguide assemblies. *J Sound Vib* 2020;465:115003.

10. Tian ZH, Xiao WF, Ma ZY, et al. Dispersion curve regression-assisted wideband local wavenumber analysis for characterizing three-dimensional (3D) profile of hidden corrosion damage. *Mech Syst Signal Process* 2021;150:107347.

11. Hamidpour M, Nami MR, Khosravifard A. An effective crack identification method in viscoelastic media using an inverse meshfree method. *Int J Mech Sci* 2021;212:106834.

12. Ichchou MN, Bareille O, Berthaut J. Identification of effective sandwich structural properties via an inverse wave approach. *Eng Struct* 2008;30(10):2591–604.

13. Santoni A, Schoenwald S, Van Damme B, et al. Determination of the elastic and stiffness characteristics of cross-laminated timber plates from flexural wave velocity measurements. *J Sound Vib* 2017;400:387–401.

14. Droz C, Bareille O, Ichchou MN. A new procedure for the determination of structural characteristics of sandwich plates in medium frequencies. *Compos B* 2017;112:103–11.

15. Berthaut J, Ichchou MN, Jezequel L. K-space identification of apparent structural behaviour. *J Sound Vib* 2005;280(3–5):1125–31.

16. Tufano G. K-space analysis of complex large-scale periodic structures [dissertation]. Ecully: Lyon University and Leuven, KU Leuven; 2020.

17. Lajili R, Bareille O, Bouazizi M, et al. Composite beam identification using a variant of the inhomogeneous wave correlation method in presence of uncertainties. *Eng Comput* 2018;35(6):2126–64.

18. Lajili R, Chikhaoui K, Zergoune Z, et al. Impact of the vibration measurement points geometric coordinates uncertainties on two-dimensional k-space identification: application to a sandwich plate with honeycomb core. *Mech Syst Signal Process* 2022;167:108509.

19. Geslain A, Raetz S, Hiraiwa M, et al. Spatial Laplace transform for complex wavenumber recovery and its application to the analysis of attenuation in acoustic systems. *J Appl Phys* 2016;120(13):135107.

20. Yan GQ, Raetz S, Groby JP, et al. Estimation via laser ultrasonics of the ultrasonic attenuation in a polycrystalline aluminum thin plate using complex wavenumber recovery in the vicinity of a zero-group-velocity lamb mode. *Appl Sci* 2021;11(15):6924.

21. Roozen NB, Labelle L, Leclère Q, et al. Non-contact experimental assessment of apparent dynamic stiffness of constrained-layer damping sandwich plates in a broad frequency range using a Nd: YAG pump laser and a laser Doppler vibrometer. *J Sound Vib* 2017;395:90–101.

22. Morandi F, Robin O, Barbaresi L, et al. Benchmarking of methods for the identification of flexural wavenumbers in wooden plates. *Proceedings of the 23rd international congress on acoustics*; 2019. p. 1294–301.

23. Margerit P, Lebée A, Caron JF, et al. The high-resolution wavevector analysis for the characterization of the dynamic response of composite plates. *J Sound Vib* 2019;458:177–96.

24. Boukadia RF, Claeys C, Droz C, et al. An inverse convolution method for wavenumber extraction (INCOME): formulations and applications. *J Sound Vib* 2022;520:116586.

25. Ribeiro LHMS, Dal Poggetto VF, Huallpa BN, et al. Bloch wavenumber identification of periodic structures using Prony's method. *Mech Syst Signal Process* 2022;178:109242.

26. Li XF, Ichchou M, Zine A, et al. An algebraic wavenumber identification (AWI) technique under stochastic conditions. *Mech Syst Signal Process* 2023;188:109983.

27. Li XF, Ichchou M, Zine A, et al. Wavenumber identification of 1D complex structures using algebraic wavenumber identification (AWI) technique under complex conditions. *J Sound Vib* 2023;548:117524.

28. Fliess M, Sira-Ramírez H. An algebraic framework for linear identification. *ESAIM Control Optim Calc Var* 2003;9:151–68.

29. Allemang RJ. The modal assurance criterion—Twenty years of use and abuse. *Sound Vib* 2003;37(8):14–23.

30. Fliess M, Mboup M, Mounier H, et al. Questioning some paradigms of signal processing via concrete examples. *Algebr Meth Flatness Signal Process State Estim*. 2003. p. 1–21.

31. Beltran-Carbajal F, Silva-Navarro G. On the algebraic parameter identification of vibrating mechanical systems. *Int J Mech Sci* 2015;92:178–86.

32. Stoica P, Selen Y. Model-order selection: a review of information criterion rules. *IEEE Signal Process Mag* 2004;21(4):36–47.

33. Mariani A, Giorgetti A, Chiani M. Model order selection based on information theoretic criteria: design of the penalty. *IEEE Trans Signal Process* 2015;63(11):2779–89.

34. Badeau R, David B, Richard G. Selecting the modeling order for the ESPRIT high resolution method: An alternative approach. In: *2004 IEEE international conference on acoustics, speech, and signal processing*. New York: IEEE; 2004. p. ii–1025.

35. Papy JM, De Lathauwer L, Van Huffel S. A shift invariance-based order-selection technique for exponential data modelling. *IEEE Signal Process Lett* 2007;14(7):473–6.

36. Albert R, Galarza CG. Model order selection for sum of complex exponentials. *2021 IEEE URUCON*. New York: IEEE; 2021. p. 561–5.

37. Margerit P, Lebée A, Caron JF, et al. High resolution wavenumber analysis (HRWA) for the mechanical characterisation of viscoelastic beams. *J Sound Vib* 2018;433:198–211.

38. Okumura S, Nguyen VH, Taki H, et al. Rapid high-resolution wavenumber extraction from ultrasonic guided waves using adaptive array signal processing. *Appl Sci* 2018;8(4):652.

39. Wang C, Reddy J, Lee K Shear deformable beams and plates: *Relationships with classical solutions*. Amsterdam: Elsevier Science; 2000. p. 39–54.

40. Caillet J, Carmona JC, Mazzoni D. Estimation of plate elastic moduli through vibration testing. *Appl Acoust* 2007;68(3):334–49.

41. Muthurajan KG, Sanakaranarayanasamy K, Rao BN. Evaluation of elastic constants of specially orthotropic plates through vibration testing. *J Sound Vib* 2004;272(1–2):413–24.

42. Au FTK, Wang MF. Sound radiation from forced vibration of rectangular orthotropic plates under moving loads. *J Sound Vib* 2005;281(3–5):1057–75.

43. Andreassen E, Jensen JS. Analysis of phononic bandgap structures with dissipation. *J Vib Acoust* 2013;135(4):041015.