# Metamodeling elastic wave propagation using a mixed factorized Fourier encoder-decoder for online laser-ultrasound testing in additive manufacturing

Frédéric Allaire, Filippo Gatti, Roberto Miorelli, Vincent Dorval, Jérôme Laurent

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

## To cite this version:

Frédéric Allaire, Filippo Gatti, Roberto Miorelli, Vincent Dorval, Jérôme Laurent. Metamodeling elastic wave propagation using a mixed factorized Fourier encoder-decoder for online laser-ultrasound testing in additive manufacturing. Engineering Applications of Artificial Intelligence, 2025, 160 (Part C), pp.111893. 10.1016/j.engappai.2025.111893. hal-05011076v2

HAL Id: hal-05011076

https://hal.science/hal-05011076v2

Submitted on 26 Aug 2025

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

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

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

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

Research paper

# Metamodeling elastic wave propagation using a mixed factorized Fourier encoder-decoder for online laser-ultrasound testing in additive manufacturing

Frédéric Allaire<sup>a, ID,*</sup>, Filippo Gatti<sup>b, ID,*</sup>, Roberto Miorelli<sup>a, *</sup>, Vincent Dorval<sup>a, ID</sup>, Jérôme Laurent<sup>a</sup>

## ARTICLE INFO

**Keywords:**

Surrogate modeling

Deep learning

Additive manufacturing

Polycrystalline media

Laser-ultrasound

Surface acoustic waves

## ABSTRACT

Laser-ultrasound (LU) testing has emerged as a promising technique for characterizing the polycrystalline microstructure of metal components produced by wire-laser additive manufacturing (WLAM), with potential for real-time online application. Numerical models simulating elastic waves propagation provide valuable insights into the relationship between microstructural properties and laser-induced displacements, but their computational cost renders them impractical for automated high-throughput characterization. To overcome this limitation, we build a metamodel that maps a wide variety of two-dimensional anisotropic polycrystalline microstructures — simplified but representative of features commonly observed in WLAM — to simulated surface displacements. Addressing this challenging high-dimensional regression problem, several neural network surrogates relying on a novel combination of layers are investigated. Their architectures include usual convolutional encoder-decoder elements and spectral layers inspired from the Fourier neural operator (FNO) framework. These layers are adapted and some variants are proposed to provide versatility in network design. All metamodels can run both a forward and backward pass at least 100 times faster than a single forward call of the original model. The best architecture implies a trade-off between computational cost and accuracy. Notably, the architecture involving the channel-wise factorized variant of the spectral layers, which is characterized by a relatively small number of parameters, achieved the lowest approximation error. The metamodel successfully captures the primary effects of anisotropy on wave propagation, even for low-anisotropy inputs not included in the training data. These findings represent a promising initial step towards addressing inversion problems and facilitating the development of online LU testing protocols in additive manufacturing.

## 1. Introduction

Nondestructive testing (NDT) techniques are widely employed across various industrial sectors to control the quality of manufactured components without causing damage. These techniques are crucial for assessing the safety and durability of parts, especially in critical applications such as nuclear, aerospace, and biomedical. Additive manufacturing (AM) processes, commonly referred to as 3D printing, enables the production of components with complex geometries and minimizes material removal during post-processing. The present study focuses on WLAM, a metal additive manufacturing (MAM) technique in which a metallic wire is melted using a laser source (see for instance (Ghanadi and Pasebani, 2024) for a review). Other MAM processes may involve metal powder feedstock instead of wire and/or employ an energy source other than laser.

### 1.1. Potential and challenges of online laser-ultrasound for AM

Recent reviews on WLAM have emphasized the critical role of monitoring systems and control techniques in optimizing both process efficiency and part quality (Ghanadi and Pasebani, 2024; Abuabiah et al., 2023). Common sensors used in MAM systems include pyrometers, optical cameras, and infrared cameras. However, several NDT techniques can provide valuable supplementary information (Segovia Ramírez et al., 2023). In this context, Shaloo et al. (2022) presented a review of NDT techniques in Wire-Arc Additive Manufacturing (WAAM), discussed their principles and potential for defect detection, which are largely applicable to other MAM processes such as WLAM. Similarly, Honarvar and Varvani-Farahani (2020) reviewed and highlighted the potential of ultrasonic testing (UT) applications in AM. Since all AM processes involve the deposition of material in a layer-by-layer fashion,

<footnote><sup>a</sup> Université Paris-Saclay, CEA, List, F-91120, Palaiseau, France</footnote>

<footnote><sup>b</sup> Université Paris-Saclay, Laboratoire de Mécanique Paris-Saclay, CentraleSupélec, ENS Paris-Saclay, CNRS, Gif-sur-Yvette, 91190, France</footnote>

<footnote>* Corresponding authors.
E-mail addresses: frederic.allaire@thalesgroup.com (F. Allaire), filippo.gatti@centralesupelec.fr (F. Gatti), roberto.miorelli@cea.fr (R. Miorelli).</footnote>
## Nomenclature

Several acronyms can be used to refer to the WLAM process:

<table><tr><td>WLAM</td><td>Wire-laser additive manufacturing<br>Wire-based laser additive manufacturing<br>Wire-feed laser additive manufacturing</td></tr><tr><td>LWAM</td><td>Laser wire-feed additive manufacturing</td></tr><tr><td>LMWD</td><td>Laser metal wire deposition</td></tr><tr><td>W-LDED</td><td>Wire-laser directed energy deposition</td></tr><tr><td>LMD</td><td>Laser metal deposition<br>(NB: may involve powder or wire feedstock, unless 'w' or 'W' is included in the acronym)</td></tr><tr><td>LDED</td><td>Laser directed energy deposition<br>(NB: may involve powder or wire feedstock, unless 'w' or 'W' is included in the acronym)</td></tr></table>

The following abbreviations are used in the present manuscript:

<table><tr><td>LU</td><td>Laser-ultrasound</td><td>ML</td><td>Machine learning</td></tr><tr><td>UT</td><td>Ultrasonic testing</td><td>DL</td><td>Deep learning</td></tr><tr><td>NDT</td><td>Nondestructive testing</td><td>AI</td><td>Artificial intelligence</td></tr><tr><td>WLAM</td><td>Wire-laser additive manufacturing</td><td>NN</td><td>Neural network</td></tr><tr><td>AM</td><td>Additive manufacturing</td><td>CNN</td><td>Convolutional neural network</td></tr><tr><td>MAM</td><td>Metal additive manufacturing</td><td>DNN</td><td>Deep neural network</td></tr><tr><td>WAAM</td><td>Wire-arc additive manufacturing</td><td>NO</td><td>Neural operator</td></tr><tr><td>L-PBF</td><td>Laser-powder bed fusion</td><td>FNO</td><td>Fourier neural operator</td></tr><tr><td>IN718</td><td>Inconel 718</td><td>F-FNO</td><td>Factorized-Fourier neural operator</td></tr><tr><td>SRAS</td><td>Spatially resolved acoustic spectroscopy</td><td>UNO</td><td>U-shaped neural operator</td></tr><tr><td>EBSD</td><td>Electron backscatter diffraction</td><td>PINN</td><td>Physics-informed neural network</td></tr><tr><td>PDE</td><td>Partial differential equation</td><td>Conv2D</td><td>2D convolution</td></tr><tr><td>FEM</td><td>Finite element method</td><td>PReLU</td><td>Parametric rectified linear unit</td></tr><tr><td>(I)FFT</td><td>(Inverse) fast Fourier transform</td><td>SSIM</td><td>Structural similarity index</td></tr><tr><td>DOE</td><td>Design of experiments</td><td>CWT</td><td>Continuous Wavelet Transform</td></tr><tr><td>LHS</td><td>Latin hypercube sample</td><td>EM</td><td>Envelope misfit</td></tr><tr><td>GT</td><td>Ground truth</td><td>PM</td><td>Phase misfit</td></tr></table>

## Elastodynamics

<table><thead><tr><th>Notations</th><th>Meaning</th><th>Specific notation for d = 2</th></tr></thead><tbody><tr><td>d ∈ {1, 2, 3}</td><td>Space dimension</td><td></td></tr><tr><td>Ω ⊂ ℝ<sup>d</sup></td><td>Material domain</td><td></td></tr><tr><td>T = 3 μs</td><td>Final time of wave propagation</td><td></td></tr><tr><td>t ∈ [0, T]</td><td>Time variable</td><td></td></tr><tr><td>x = (x<sub>i</sub>)<sub>1≤i≤d</sub> ∈ Ω</td><td>Space vector and coordinates</td><td>x = (x<sub>i</sub>)<sub>1≤i≤d</sub></td></tr><tr><td>u = (u<sub>i</sub>)<sub>1≤i≤d</sub></td><td>Displacement vector and coordinates</td><td>u = (u<sub>x</sub>)<sub>1≤i≤d</sub></td></tr><tr><td>σ = (σ<sub>ij</sub>)<sub>1≤i,j≤d</sub></td><td>Cauchy stress tensor</td><td>σ = (σ<sub>xx</sub> σ<sub>xy</sub>)<sub>1≤i,j≤d</sub></td></tr><tr><td>ε = (ε<sub>ij</sub>)<sub>1≤i,j≤d</sub></td><td>Strain tensor</td><td></td></tr><tr><td>c = (c<sub>ij</sub>)<sub>1≤i,j,k,j≤d</sub></td><td>Elasticity tensor (4th order)</td><td></td></tr><tr><td>f = (f<sub>i</sub>)<sub>1≤i≤d</sub></td><td>Body force</td><td></td></tr><tr><td>δ<sub>ij</sub></td><td>Kronecker's delta</td><td></td></tr><tr><td>n ∈ N</td><td>Number of nodes in FEM</td><td></td></tr><tr><td>M<sub>n</sub> ∈ ℝ<sup>n×n</sup>, K<sub>n</sub> ∈ ℝ<sup>n×n</sup></td><td>Mass matrix and stiffness matrix in FEM</td><td></td></tr><tr><td>f<sub>n</sub> ∈ ℝ<sup>n</sup></td><td>Nodal vector of external forces in FEM</td><td></td></tr><tr><td>u<sub>n</sub> ∈ ℝ<sup>n</sup></td><td>Nodal displacement vector in FEM</td><td></td></tr><tr><td>ρ = 8200 kg/m<sup>3</sup></td><td>Material density</td><td></td></tr><tr><td>C<sub>11</sub>, C<sub>12</sub>, C<sub>44</sub></td><td>Monocrystal elasticity constants (cubic anisotropy)</td><td></td></tr><tr><td>A = (2C<sub>44</sub>) / (C<sub>11</sub> - C<sub>12</sub>)</td><td>Zener's cubic anisotropy ratio</td><td></td></tr><tr><td>v<sub>QP</sub>, v<sub>P</sub> = 5800 m s<sup>-1</sup></td><td>(Quasi-)longitudinal wave speed</td><td>x<sub>s</sub> = 12.5 mm<br>σ<sub>s</sub> ≈ 0.255 mm</td></tr><tr><td>v<sub>QS</sub>, v<sub>S</sub> = 3100 m s<sup>-1</sup></td><td>(Quasi-)transverse wave speed</td><td>k<sub>el</sub> ∈ [1, 4]</td></tr><tr><td>τ ≈ 30.3 ns</td><td>Temporal width parameter of laser source<br>Horizontal coordinate of laser spot center<br>Gaussian standard deviation for laser spot<br>Elongation factor<br>Elongation angle (relative to vertical axis y)<br>Crystal orientation angle<br>Angle of wave vector direction of qP/qS wave in (x, y) plane</td><td>α<sub>el</sub> ∈ [-π/4, π/4]<br>φ ∈ [-π/4, π/4]<br>θ ∈ [-π, π]</td></tr></tbody></table>
## Metamodel

<table><thead><tr><th>Notation</th><th>Meaning</th><th>Metamodel counterpart</th></tr></thead><tbody><tr><td><i>M</i></td><td>Numerical model</td><td><i>M̃</i></td></tr><tr><td><i>N</i><sub>x</sub> = 2000, <i>N</i><sub>y</sub> = 800</td><td>Size of discretized crystal orientation field <i>φ</i>(<i>x</i>, <i>y</i>)</td><td><i>Ñ</i><sub>x</sub> = 500, <i>Ñ</i><sub>y</sub> = 200</td></tr><tr><td><i>a</i> ∈ ℝ<sup><i>N</i><sub>x</sub>×<i>N</i><sub>y</sub></sup> × ℝ<sup>3</sup></td><td>Set of variable model inputs</td><td><i>ã</i> ∈ ℝ<sup><i>Ñ</i><sub>x</sub>×<i>Ñ</i><sub>y</sub></sup> × ℝ<sup>3</sup></td></tr><tr><td><i>a</i><sub>VT</sub></td><td>Continuous representation common to <i>a</i> and <i>ã</i></td><td></td></tr><tr><td>(<i>x</i>[l])<sub>1≤<i>j</i>≤<i>N</i><sub>r</sub></sub>=81</td><td>Receiver positions on top surface</td><td></td></tr><tr><td>(<i>t</i>[<i>k</i>])<sub>1≤<i>k</i>≤<i>N</i><sub>r</sub></sub>=1000</td><td>Regular discretization of signal time series</td><td></td></tr><tr><td><i>b</i> ∈ ℝ<sup><i>N</i><sub>r</sub>×<i>N</i><sub>r</sub></sup></td><td>Output B-scan: signal amplitude at receiver locations</td><td><i>b̃</i> ∈ ℝ<sup><i>N</i><sub>r</sub>×<i>N</i><sub>r</sub></sup></td></tr><tr><td><i>C</i>' ∈ ℝ<sup>4×<i>N</i><sub>y</sub>×<i>N</i><sub>x</sub></sup></td><td>Spatially discretized elasticity components in Eq. (5)</td><td><i>C̃</i>' ∈ ℝ<sup>4×<i>Ñ</i><sub>y</sub>×<i>Ñ</i><sub>x</sub></sup></td></tr></tbody></table>

## Spectral layers

<table><thead><tr><th>Notation</th><th>Meaning</th></tr></thead><tbody><tr><td><i>F̃</i>, <i>F̃</i><sup>-1</sup></td><td>Fourier transform and its inverse</td></tr><tr><td><i>q</i> ∈ ℝ<sup>*</sup></td><td>Input dimension of functions <i>r</i> and <i>v</i></td></tr><tr><td><i>D</i> ⊂ ℝ<sup><i>q</i></sup></td><td>Input space of functions <i>r</i> and <i>v</i></td></tr><tr><td><b>x</b> = (<i>x</i><sub>i</sub>)<sub>1≤<i>i</i>≤<i>q</i></sub> ∈ <i>D</i></td><td>Input vector space and coordinates</td></tr><tr><td><b>k</b> = (<i>k</i><sub>i</sub>)<sub>1≤<i>i</i>≤<i>q</i></sub> ∈ ℝ<sup><i>q</i></sup></td><td>Mode vector in Fourier space</td></tr><tr><td><i>N</i><sub>in</sub> ∈ ℝ<sup>*</sup></td><td>Number of input channels in a spectral layer</td></tr><tr><td><i>N</i><sub>out</sub> ∈ ℝ<sup>*</sup></td><td>Number of output channels in a spectral layer</td></tr><tr><td><i>N</i> ∈ ℝ<sup>*</sup></td><td>Number of channels in a spectral layer (if constant)</td></tr><tr><td><i>K</i> : <b>v</b> ↦ <i>K</i>[<b>v</b>]</td><td>Kernel integral operator defined in Eq. (12)</td></tr><tr><td><b>v</b> : <i>D</i> → ℝ<sup><i>N</i><sub>in</sub></sup></td><td>Input function of <i>K</i> and spectral layer</td></tr><tr><td><b>r</b> : <i>D</i> → ℝ<sup><i>N</i><sub>out</sub>×<i>N</i><sub>in</sub></sup></td><td>Implicit parameter function of <i>K</i></td></tr><tr><td><b>R</b> ∈ ℝ<sup><i>S</i></sup></td><td>Parameter tensor of spectral layer (shape <i>S</i> depends on layer variant)</td></tr><tr><td><b>Γ</b> ∈ ℝ<sup><i>N</i><sub>out</sub>×<i>N</i><sub>in</sub></sup></td><td>Parameter matrix in channel-wise spectral layer</td></tr><tr><td><b>B</b> ∈ ℝ<sup><i>N</i><sub>out</sub></sup></td><td>Bias parameter vector in a spectral layer</td></tr><tr><td><i>N</i><sub>1</sub>, ..., <i>N</i><sub><i>q</i></sub> ∈ ℝ</td><td>Size of regular discretization of <i>D</i> per dimension</td></tr><tr><td><i>K</i><sub>1</sub>, ..., <i>K</i><sub><i>q</i></sub> ∈ ℝ</td><td>Number of modes involved in spectral layer per dimension</td></tr><tr><td><i>m</i><sub>1</sub>, ..., <i>m</i><sub><i>q</i></sub> ∈ ℝ</td><td>Number of non-negative modes required per dimension</td></tr><tr><td colspan="3"><strong>Others</strong></td></tr><tr><td><i>z</i></td><td>Generic variable/vector / function in the present table</td></tr><tr><td><i>V</i>(<i>z</i><sub>1</sub>, <i>z</i><sub>2</sub>)</td><td>Uniform distribution on interval (<i>z</i><sub>1</sub>, <i>z</i><sub>2</sub>)</td></tr><tr><td><i>b</i><sub>avg</sub></td><td>Average amplitude (over time and space) in B-scan <i>b</i></td></tr><tr><td><i>b̂</i><sub>train</sub></td><td>Average B-scan based on train dataset</td></tr><tr><td><i>n</i><sub>train</sub> = 10 000</td><td>Size of train dataset</td></tr><tr><td><i>n</i><sub>test</sub> = 1000</td><td>Size of test dataset</td></tr><tr><td>|<i>z</i>|</td><td>(Element-wise) absolute value/modulus of <i>z</i></td></tr><tr><td>|<i>z</i>|, [<i>z</i>] ∈ ℝ</td><td>Floor and ceil values of <i>z</i> ∈ ℝ<sub>+</sub></td></tr><tr><td>||<i>b</i>||<sub><i>p</i></sub></td><td><i>L</i><sub><i>p</i></sub>-norm of B-scan <i>b</i> treated as a vector of length <i>N</i><sub>r</sub> × <i>N</i><sub>t</sub></td></tr><tr><td><i>ℓ</i><sub><i>p</i></sub>, <i>L</i><sub><i>p</i></sub>, <i>L</i><sub><i>p</i>,rel</sub></td><td>Error metrics between B-scans, defined in Eqs. (8)–(10)</td></tr><tr><td><i>V</i><sub><i>a</i></sub><i>z</i> ∈ ℝ<sup><i>N</i><sub>x</sub>×<i>N</i><sub>y</sub></sup> × ℝ<sup>3</sup></td><td>Gradient of <i>z</i> : <i>a</i> ↦ <i>z</i>(<i>a</i>) ∈ ℝ (<i>a</i> is treated as a vector of length <i>N</i><sub>x</sub> × <i>N</i><sub>y</sub> + 3)</td></tr><tr><td><i>V</i><sub><i>ã</i></sub><i>z̃</i> ∈ ℝ<sup><i>Ñ</i><sub>x</sub>×<i>Ñ</i><sub>y</sub></sup> × ℝ<sup>3</sup></td><td>Gradient of <i>z̃</i> : <i>ã</i> ↦ <i>z̃</i>(<i>ã</i>) ∈ ℝ (<i>ã</i> is treated as a vector of length <i>Ñ</i><sub>x</sub> × <i>Ñ</i><sub>y</sub> + 3)</td></tr><tr><td><i>E</i>[<i>z</i>]</td><td>Expected value of <i>z</i> with regard to input distribution described in Section 2.2</td></tr><tr><td><i>f</i><sub>c</sub> = 60 MHz</td><td>Cut frequency of low-pass filter</td></tr><tr><td><i>S</i><sub>0</sub> = 0.005 G Pa<sup>-1</sup></td><td>Fixed scaling parameter in neural networks</td></tr></tbody></table>

The above notations are used in the body of the present manuscript, additional notations are defined in the appendix. The numerical values of variables set in the present study are reported. Notation **x** intentionally appears both in “Elastodynamics” and “Spectral layers”: while it would coincide for straightforward application of FNO, it is not the case for the neural networks considered in the present study.

they present an opportunity for carrying out NDT not only offline but also online. In this regard, LU techniques show considerable promise, as they are non-contact methods (all-optical) that can be applied to high-temperature materials, making them well-suited for online monitoring in MAM systems.

LU techniques involve generating a short laser pulse that briefly heats and induces localized thermal expansion in the material. This process initiates the propagation of elastic waves through the medium at ultrasonic frequencies. The resulting material displacements can be

measured and processed to characterize the materials' properties (Lian et al., 2024), such as anisotropy, elasticity, and microstructure. LU techniques are usually non-destructive; however, they can cause material ablation at high laser spot energy densities. For instance, Malmström et al. (2022) demonstrated the use of online LU testing to estimate grain size in high-temperature austenitic steel during hot strip mill processes. While LU is often performed offline, it can be used to measure elasticity tensors of anisotropic alloys, such as through transient grating spectroscopy at the monocrystal scale (Stoklasová et al., 2021), or using
<figure><img src="image_4.png" /><figcaption>**Fig. 1.** The schematic configuration of a LU system for online monitoring in WLAM involves the following setup: after material deposition, a pulsed source laser heats the surface, inducing localized thermal expansion and generating elastic waves that propagate through the material. A receiver interferometer laser (on the left) is used to measure the vertical surface displacements resulting from wave propagation.</figcaption></figure>

spatially resolved acoustic spectroscopy (SRAS) over a polycrystalline surface (Dryburgh et al., 2022). Additionally, SRAS may also be used for rapid estimation of crystal orientation mappings, which are more commonly obtained using electron backscatter diffraction (EBSD) techniques. In the context of MAM, Harke et al. (2022) employed offline LU combined with numerical simulation to detect defects in a component built by Laser Powder Bed Fusion (L-PBF). The study identified surface and subsurface defects by analyzing ultrasound scattering caused by these irregularities, demonstrating the potential of LU for defect characterization in L-PBF processes.

LU-based online characterization of WLAM beads is schematically represented in Fig. 1. The LU system relies on a laser source (emitter) and a laser interferometry device (receiver), with both laser beams focalized on the top surface of the deposited bead. In this configuration, the source generates a broadband ultrasonic signal, while the receiver measures the out-of-plane absolute displacements resulting from elastic wave propagation through the polycrystalline medium. These measured signals contain valuable surface and sub-surface information about the material properties between the source and receiver at the bead scale (roughly 1 mm in height). This configuration is technically feasible for online monitoring, as demonstrated by Bakre et al. (2022) where online LU measurements of powder-based MAM were shown to have quality comparable to that of offline measurements. In the current configuration, the received signals would ideally be processed automatically to detect anomalies in expected microstructural properties, entailing adjustments to the WLAM process parameters. For instance, occurrence of defects that are difficult to identify with other sensors, such as cracks and porosities, may warrant a change in the deposition strategy or an early stop of the process, depending on defect severity. Some studies showed the possibility of using LU with a broadband source signal to estimate elasticity constants (Gasteau et al., 2016) or detecting defects in L-PBF components (Harke et al., 2022) in offline settings. However, the application of LU for online diagnostic remains challenging, particularly in rapidly and automatically extracting microstructural information from a limited set of measurements. A common approach to facilitate inversion in LU-based characterization techniques involves the excitation of narrow frequency ranges (e.g., via grating, heterodyne source reception) instead of using a raw broadband source signal. However, this technique is typically applicable to smaller scales (monocrystal) and would be too time-consuming for the online characterization of an entire bead.

## 1.2. Modeling ultrasound propagation in WLAM microstructures

A complex interplay of phenomena related to elastic wave propagation must be considered when applying LU techniques to WLAM

components. Notably, grain morphology and polycrystal texture, which are affected by the strong thermal gradients typical of MAM, influence wave propagation. This influence is manifested in varying phase wave velocity depending on the direction of propagation and wave scattering at grain boundaries. These effects are particularly pronounced in highly anisotropic materials. Additionally, other influential factors, such as surface roughness (Zhang et al., 2024) and sensor operating conditions can contribute to noise in measurements, which may subsequently impact the accuracy of the resulting data.

Analytical and numerical models based on classical elastodynamics equations are complementary approaches for relating the anisotropic properties of polycrystalline materials to ultrasonic data. On the one hand, analytical models typically rely on some physical approximations and assumptions regarding grain morphology and other material properties to estimate physical quantities such as average wave velocities and scattering attenuation as a function of frequency, grain size, and elastic constants (Stanke and Kino, 1984; Weaver, 1990; Kube and Turner, 2015; Arguelles and Turner, 2017). On the other hand, numerical models, often employing finite element methods (FEM), can simulate elastic wave propagation either in 2D or 3D media. Numerous studies have demonstrated strong agreement between analytical and FEM-based approaches (Van Pamel et al., 2015; Ryzy et al., 2018; Sha et al., 2020; Huang et al., 2022). In theory, FEM fully captures the physics of multiple scattering of elastic waves, with the primary source of approximation arising from the numerical scheme used to solve the system of partial differential equations (PDE). Another advantage of FEM over analytical methods is their ability to determine displacements at any point and time within an explicitly defined microstructure, thus enabling straightforward comparison with experimental ultrasonic signals. However, this capability often requires fine spatiotemporal resolution, consequently leading to significantly higher computational costs, particularly when modeling 3D media.

Regarding LU testing, some studies (Every et al., 2013; Zhang et al., 2024; Tu et al., 2024) have adopted a multi-physics modeling approach to include thermal, optical, and electromagnetic effects induced by the laser source, thereby simulating wave propagation with higher fidelity, often using FEM. However, to the authors' knowledge, most of these studies consider propagation through isotropic media, with only occasional consideration of anisotropy (Tu et al., 2024), and none explicitly address polycrystalline media. Purely elastodynamics models can partially account for the thermal and electromagnetic effects through a boundary condition representative of the LU-source (Scruby and Drain, 1990); this approach has been applied to the numerical simulation of elastic wave propagation in polycrystals (Grabec et al., 2022).
Regardless of the model used, it is essential to account for spe-
cific material properties. For instance, grain elongation and size can
have a strong effect on elastic wave scattering in polycrystals (see
for instance Arguelles et al., 2016; Bai et al., 2020; Huang et al.,
2021; Victoria Giraldo and Tie, 2024). In WLAM, microstructures often
exhibit large, elongated grains and textured polycrystals (i.e., pre-
ferred crystal orientations) due to the strong thermal gradients present
during the manufacturing process. These microstructural features are
influenced by process parameters and may vary between deposited
layers. Using idealized microstructures as input to generate ultrasound
time-series as output, numerical models provide a forward mapping
of wave propagation. However, these models are too slow for real-
time applications, even in 2D. The inversion problem which consists in
inferring microstructural properties from online LU data is even more
challenging. Inversion is inherently ill-posed, more computationally
expensive (often requiring multiple forward model executions), and
highly sensitive to discrepancies between the models and observed data
due to model assumptions and simplifications, or noisy measurements.

1.3. *Applications of machine learning to MAM and UT*

Machine learning (ML) has recently received significant attention
from the scientific community due to its potential applications to indus-
trial MAM (Safdar et al., 2024). Although, to the authors' knowledge,
no review specifically addressing ML in WLAM exists, the interested
reader may find valuable insights in reviews covering applications
of ML and/or artificial intelligence (AI) to MAM (Guo et al., 2022),
laser-based MAM (Cai et al., 2023), and WAAM (He et al., 2023).
The following major contributions of ML, applicable to WLAM, can be
identified.

1. *Offline process design* primarily focuses on defining an operating window for MAM process parameters. ML models facilitate this by providing a mapping between variables describing input process parameters such as laser power and scanning speed to system behavior including characteristics like bead profile.

2. *Classical online monitoring and control* aims to detect and act upon process anomalies observed at the external macroscopic scale. ML models are trained to process data from cameras and/or pyrometers acquired online as input to estimate melt pool size, part geometry and/or temperature. Due to the relatively straightforward relation between these output quantities and input process parameters, such models can be integrated into closed-loop control algorithms based on differences between expected and estimated values. Online application of LU testing to WLAM relates to this category, although the relationship between microstructure and MAM process conditions for control purposes is not straightforward.

3. For *asynchronous control*, the ML-based mapping involves input data which are acquirable online and output data that are typically measured offline. These outputs are usually related to input data and MAM process parameters in a complex manner. For instance, *Jamnikar et al.* (2022) describes the measurement of microstructural properties, namely alpha lath thickness and beta grain sizes of Ti-6Al-4V alloy, using scanning electron microscopy of the cross-sections of WLAM beads. Then, they trained (offline) a convolutional neural network (CNN) model mapping melt pool images and temperatures previously acquired online to the aforementioned microstructural data. This approach facilitates real-time online predictions of the microstructural characteristics.

4. *Path planning* represents another domain where ML (more gener-
ally AI) techniques are applied in an offline manner. The primary
objective of these techniques is to find admissible and/or optimal
trajectories for metal deposition, tailored to achieve the desired
part geometry.

5. *Virtualization* is dedicated to the development of human-machine interfaces, which may be used offline for training purposes or online as support tools. While it predominantly falls within the realm of AI, it may incorporate ML techniques.

Regarding the first three areas, most studies rely on experimental data. However, numerical models can also be used to generate a database of synthetic data to train ML models and/or evaluate the feasibility of these approaches. As a matter of fact, metamodeling has been applied to MAM process models to provide near-real-time predictions of 3D temperature field (Hemmasian et al., 2023). Such methods are mostly used for offline design studies, although they show potential for innovative online control strategies.

ML techniques, most particularly deep neural networks (DNN), have often been combined with ultrasonic sensing measurements, notably for automation purposes (Cantero-Chinchilla et al., 2022). In LU testing, Lv et al. (2023) used broadband signals and DNN to automatically predict the presence of defects. Common LU applications such as grain size estimation benefit from ML techniques (Yu et al., 2023). Using either time series (1D) or image-like (2D) UT input data, CNN have demonstrated strong performance for classification and regression problems. These include defect detection with 1D input (Cai et al., 2021), regression of corrosion height with 1D input followed by multiple model executions to generate a 1D corrosion profile (Cantero-Chinchilla et al., 2023), and defect localization and sizing (i.e., regression with three output scalar quantities) using 2D input (Miorelli et al., 2021). A key advantage of DNN is their ability to quickly process high-dimensional UT data that are sometimes difficult to interpret even after dedicated processing. Common hurdles in the applicability of DNN are their generalizability and interpretability. Furthermore, in UT applications, experimental data acquisition can be costly and may only apply to very specific use cases due to the high dependence on part geometry, material properties, defect type, defect location, defect size, UT device, etc. Also, detailed material information requires time-consuming measurement techniques such as computed tomography (CT-X) and EBSD. This lack of information and data poses a significant limitation, particularly in MAM. Even though a large quantity of LU-induced signals can be acquired under various process conditions, microstructural information will only be known partially. An intuitive approach for compensating limited data availability is to make use of synthetic UT datasets by means of simulation and/or data augmentation as highlighted in a recent review (Uhlig et al., 2023). In a similar context, although not limited to UT techniques, Wu et al. (2024) reviewed the prediction of residual stresses in laser-based AM components and pointed out the potential of ML based on experimental and/or synthetic data.

1.4. *Contributions and related works*

In the present study, we propose to rely on numerical simulation of LU-induced elastic wave propagation in 2D polycrystals based on an elastodynamics model solved by FEM. While considering 3D media would yield more physically-accurate data (especially with a multi-physics approach), the time and computer memory required to generate a sufficiently large ensemble of such simulations is prohibitive. This modeling approach was deemed as the best compromise between physical representativity and computational cost, but is still too time consuming for online LU-testing of WLAM microstructures. We propose to address this problem by means of metamodeling, which consists in providing an approximation of the output of interest at a significantly lower computational cost compared to the original model. Metamodels are commonly employed when the original model is computationally expensive and/or in applications that require many forward model executions, such as uncertainty quantification, optimization, and inverse problems (Zhao et al., 2024). Most metamodels involve a statistical regression or interpolation technique whose parameters are fitted using a dataset of input-output pairs generated from the original model.
In most metamodeling studies, only a handful of model parameters are considered as input variables and a few scalar output quantities of interest are selected. In the present study, however, both input microstructure and output surface displacements are high-dimensional data for which traditional methods are generally ill-suited and suffer from the curse of dimensionality, even after application of dimension reduction techniques such as principal component analysis. On a positive side, DNN are well-suited to handle the complexity of high-dimensional data, although they often require a large training database to achieve satisfactory accuracy.

In the field of metamodeling of elastic wave propagation, there are a number of objectives that can be pursued. Some studies aim to emulate the full propagation over time in a variety of configurations. For instance, starting from the initial state, the displacement field at a given time can be used as input together with other variables to predict the subsequent displacement field. In Gantala and Balasubramaniam (2022), such a DNN is trained given single and multi-point excitation sources with variable positions. Alternatively, time can be considered as one of the variable input parameters, which is the strategy proposed in Hemmasian et al. (2023). Nevertheless, the reproduction of the whole displacement fields with these methods requires several metamodel calls, which may prove inefficient and unnecessary for control purposes given that the in-domain displacements are not quantifiable. Furthermore, to reduce computational time, it may be preferable to consider the displacement field at all required times and locations as the output, so that the metamodel is only called once with a given model configuration. Regarding inverse problems, Singh et al. (2022) used a DNN and a generative adversarial network (GAN) to map simulated output UT data with perturbations to corresponding input microstructures.

To handle the potentially high-dimensional input and output data involved in PDE problems, the concept of neural operators (NO) has recently emerged as a promising approach. A core idea of NO is to emulate the behavior of PDE operators which are seen as a mapping between input and output functions of time and/or space regardless of the discretization used to numerically solve the PDE. Two prominent NO model architectures are Fourier Neural Operators (FNO Li et al., 2021) and DeepONet (Lu et al., 2021). For instance, a DeepONet architecture was used to train a surrogate of a high-cost computational code to predict compressible air flow (i.e., a multi-variable field) around an airfoil given variable geometry for shape optimization (Shukla et al., 2024). Notable variants of FNO are factorized Fourier neural operators (F-FNO Tran et al., 2023) and U-shaped neural operators (UNO Rahman et al., 2023), which both showed good performance in metamodeling earthquake-induced surface displacements based on variable heterogeneous 3D geology (Lehmann et al., 2024; Zou et al., 2024; Lehmann et al., 2025) – a context analogous to that of the present case study. We follow a similar approach, while making necessary adaptations due to material anisotropy, among other features. Another novelty proposed in the present work is the inclusion of spectral layers – a key element of FNO models – in a more standard CNN architecture, for which several variants are detailed and compared.

The remainder of the paper is organized into four main sections. In Section 2, the generation of the synthetic database of numerical simulations and the metamodeling problem are explained and then deep-learning (DL) surrogates and the training strategy are proposed. In Section 3, the main results of the proposed approach are shown and analyzed. These results are then discussed in Section 4. Finally, a conclusion is provided in Section 5.

## 2. Methods

The present study is a numerical investigation of the application of DL surrogates based on synthetic datasets. Several model simplifications allow for easier data representation while accounting for the

influence of anisotropy on wave propagation and addressing computational challenges. Some information regarding the specific WLAM system that motivated the present study can be found in Roch et al. (2023).

The present section is composed of six parts. The mathematical and numerical model utilized for the execution of a single elastodynamics simulation is described in Section 2.1. Then, in Section 2.2, the associated design of experiments (DOE) to obtain synthetic datasets is presented, and the metamodeling problem is formulated. In Section 2.3, a number of accuracy evaluation metrics are introduced. Section 2.4 provides detailed information on layers related to FNO models that are central to the present work. The main DNN surrogate is described in Section 2.5 and some notable variants are mentioned. Finally, Section 2.6 provides details on aspects related to DNN training and implementation.

### 2.1. Modeling and simulation of elastic waves propagation induced by laser in a polycrystalline medium

We will consider an anisotropic material over a spatial domain $\Omega \subset \mathbb{R}^d$, $d \in \{1, 2, 3\}$ with Cartesian coordinates $x_i$, $i \in [1, d]$. Regarding material properties, mass density is denoted by $\rho$ and it is assumed that the material exhibits linear anisotropic elasticity, which is described by Hooke's law: $\sigma = c : \epsilon$, where $\sigma = (\sigma_{ij})_{i,j}$ is the Cauchy stress tensor, $c = (c_{ijkl})_{i,j,k,l}$ is the material's 4th order elasticity tensor, and $\epsilon = (\epsilon_{ij})_{i,j}$ is the strain tensor. Considering that the material is subject to a body force expressed via source term $f = (f_i)_i$, the classical elastodynamics equation for displacement vector $u = (u_i)_i$ in $\Omega$ for time $t \in [0, T]$, $T > 0$ writes:

$$ \forall i \in [1, d], \left\{ \rho \delta_{il} \frac{\partial^2}{\partial t^2} - \frac{\partial}{\partial x_k} c_{kilj} \frac{\partial}{\partial x_l} \right\} u_i = f_i, \quad (1) $$

where $\delta_{il}$ stands for Kronecker's delta and using Einstein's notation. For further theoretical background, the interested reader may refer to Nayfeh (1995) (notably chapters 2 and 3) or Royer and Valier-Brasier (2022). The problem combining Eq. (1) with given initial condition and boundary conditions admits a unique solution for unknown displacement $u$. In this context, the simulation of elastic wave propagation entails the provision of an approximate solution to the aforementioned problem using a numerical model. FEM involves a variational formulation of the problem, followed by a discrete representation of vector fields $u$ and $f$. This generally leads to the introduction of numerical vectors $u_n$ and $f_n$ both in $\mathbb{R}^n$, where $n$ is the number of nodes, with corresponding mass and stiffness matrices $M_n$ and $K_n$, which are related as follows:

$$ M_n \frac{d^2}{dt^2} u_n + K_n u_n = f_n. \quad (2) $$

More details regarding the FEM used in the present study and its implementation can be found in Dorval et al. (2025).

The general simulation configuration in the present study is represented in Fig. 2. A 2D model is considered ($d = 2$) where spatial domain $\Omega$ is a rectangle with dimensions of 25 mm in length and 10 mm in height. Wave propagation is simulated up to final time $T = 3$ µs. First-order square elements (Q1) of side 12.5 µm are used for FEM. Null displacements are assumed at $t = 0$. The body source term $f$ is also null, but Neumann boundary conditions are used instead, in which case FEM will also lead to Eq. (2) with an adjusted $f_n$ representative of surface external forces. Notwithstanding the absence of thermal and electromagnetic considerations in the model, theoretical elastodynamics models showed that a vertical point force (resp. radial force dipole) with a Gaussian pulse temporal profile is a suitable model of the laser source in ablation (resp. thermoelastic) regime (Scruby and Drain, 1990). The boundary condition for the top surface is designed to provide a fair representation of the LU pulse in terms of force excitation, whereas the other three surfaces are displacement-free. Consequently, the following top boundary condition for normal stress
<figure><img src="image_5.png" /><figcaption>Fig. 2. Representation of a typical 2D medium in the simulations with source at the middle of the top surface and receivers located up to 5 mm from the source. Crystal orientation angles and grain elongation angle $\alpha_{el}$ are defined counterclockwise relative to the upward vertical $e_y$.</figcaption></figure>

component $\sigma_{yy}$, which is similar to the source laser's irradiance, is used:

$$ \sigma_{yy}(x, y = 0, t) = -\exp\left\{-\frac{(x - x_s)^2}{2\sigma_s^2}\right\} \frac{t}{\tau} \exp\left\{-\left(\frac{t}{\tau}\right)^2\right\}, \quad (3) $$

where the spatial profile is Gaussian with a spatial standard deviation $\sigma_s$ centered on the laser spot's center at coordinate $x_s = 12.5$ mm, whereas the temporal profile is a smooth positive pulse parameterized via time constant $\tau > 0$. It should be noted that in the present study, at the top surface ($y = 0$), shear stress $\sigma_{xy}$ is null, unlike normal stress $\sigma_{yy}$, which is representative of ablation regime. Here, the spatial extent is assumed to match that of the laser spot. On the one hand, nanosecond lasers typically exhibit shorter pulse durations with a full width at half maximum (FWHM) around 5~10 ns for which spot diameters down to approximately 0.2 mm can be achieved. Conversely, the spatial resolution must be significantly refined given such small values to prevent numerical dispersion effects in the FEM. A trade-off between realism and computational burden was reached by selecting higher values $\sigma_s \approx 0.255$ mm and $\tau \approx 30.3$ ns, which correspond to a 1-mm laser spot diameter and a temporal pulse with a FWHM of 35 ns.

The material's physical properties are selected to closely align with the known values of Inconel 718 (IN718) alloy at room temperature. Mass density is set to $\rho = 8200 \text{ kg/m}^3$, and reference values of compression and shear wave speeds for an equivalent homogeneous isotropic medium are respectively set to $v_p = 5800 \text{ m/s}$ and $v_s = 3100 \text{ m/s}$. At the grain scale, IN718 is commonly assumed to exhibit cubic anisotropy, but several sets of values for its cubic elasticity constants $C_{11}$, $C_{12}$, and $C_{44}$ are available in the literature. Zener's cubic anisotropy ratio $A = 2C_{44}/(C_{11} - C_{12})$ is around 2.4 for IN718, indicative of a moderate degree of positive anisotropy ($A > 1$). The propagation medium is assumed to be a polycrystal with large elongated grains and average grain surface area around $4 \text{ mm}^2$. The grain shapes are generated using a Voronoi tessellation, following the method described in Dorval et al. (2025). First, the coordinates of the so-called Voronoi nuclei (a.k.a. seeds, or centroids) are independently sampled according to a uniform distribution in $\Omega$. Then, a grain is defined as the set of points in $\Omega$ belonging to the closest nucleus among the set of nuclei according to the Euclidean distance. Due to the restriction to $\Omega$, grains located close to the domain boundaries can be much smaller than average. This procedure can be adapted to apply overall grain elongation based on given aspect ratio $k_{el} \ge 1$ and direction defined via angle $\alpha_{el}$. Finally, crystal orientations are constant grain-wise and sampled independently. In the present study, one of the crystal axes is forced along the $x_3$-axis (perpendicular to the propagation plane) for the sake of simplicity. Therefore, crystal orientation is defined by a single rotation angle $\phi$ which is assigned independently from one grain to another according

to uniform distribution $\mathcal{U}(-\pi/4, \pi/4)$. This polycrystal shows what we refer to as a transverse texture, additional details can be found in Appendix A.1.

The influence of anisotropy is illustrated by simulation results at $t \approx 1.3 \mu\text{s}$ for six different media represented in Fig. 3. Four homogeneous (monocrystal) media are considered: isotropic and anisotropic with variable values of A and/or $\phi$. Phase wave speeds for quasi-compression waves $v_{qp}$ and quasi-shear waves $v_{qs}$ both depend on the propagation direction. In this plane where propagation direction is indicated via angle $\theta$, the extrema of $v_{qp}$ are obtained for $\theta \equiv \phi \mod \pi/2$ (lowest) and $\theta \equiv \phi + \pi/4 \mod \pi/2$ (highest), hence the apparent deformation of the circular outer wavefront compared to the isotropic case. Extrema of $v_{qs}$ are obtained at the same values of $\theta$ but maximum and minimum are swapped, hence the $(\pi/4)$-tilted alterations to the inner wavefronts. The higher A, the higher the variations of wave speeds with $\theta$. For polycrystals, differences in crystal orientations for adjacent grains leads to reflection and refraction when waves reach the grains' interface. This scattering phenomenon is more pronounced when A increases as shown for the two polycrystalline media represented in Fig. 3, which have the same microstructure but different values of A. Scattering influences the surface displacements in the form of seemingly noisier signals compared to homogeneous media.

In practice, the vertical component of displacement at the top surface can be quantified through the use of a laser interferometry device. Therefore, the output of interest is $u_y(x, y = 0, t)$. We focus on time series for $N_r = 81$ receivers located up to 5 mm next to the laser source on both sides with a 0.125-mm spacing between receivers. In accordance with the simulation settings, reflections occurring at the left, right, and bottom domain boundaries will not reach these receivers. In consideration of the stability condition of the FEM, the time step is about 1.7 ns in the simulations. Output time series are interpolated to a reference time vector with a time step of 3 ns, resulting in a length of $N_t = 1000$. A time series at a given receiver is called an A-scan and their concatenation over the $N_r$ consecutive receivers is a B-scan, which is the output of interest in the present study. The units of the output B-scan are arbitrary due to the linearity of the wave equations. They are scaled via a fixed factor so that most amplitude values are contained in $[-1, 1]$, but lower values are obtained at early times for receivers closest to the source.

## 2.2. Design of experiments and problem formulation

In order to construct a metamodel, an ensemble of simulations is generated. In all simulations, the spatial and temporal extents of the simulation are identical, as are the source terms. The values of
<figure><img src="image_6.png" /><figcaption>Fig. 3. Snapshots at $t \approx 1.3 \mu s$ of the Euclidean norm $||u||_2$ of the displacement vector at each location in different media from 2D FEM simulations are presented. Top left: $A = 1$ (isotropic); top right: $A = 1.5$, $\phi = 0^{\circ}$; middle left: $A = 1.5$, $\phi = -20^{\circ}$; middle right: $A = 2.5$, $\phi = -20^{\circ}$; bottom: same polycrystal for $A = 1.5$ (left) and $A = 2.5$ (right).</figcaption></figure>

<caption>**Table 1**</caption>

Input distributions used to generate ensembles of simulations. In the overall input distribution, these four inputs are independent.

<table><thead><tr><th>Input</th><th>Distribution</th></tr></thead><tbody><tr><td>A</td><td>V(1.5, 2.5)</td></tr><tr><td>log k<sub>el</sub></td><td>V(0, log 4)</td></tr><tr><td>α<sub>el</sub></td><td>V(-π/4, π/4)</td></tr><tr><td>Microstructure (prior to elongation)</td><td>Voronoi & transverse texture</td></tr></tbody></table>

monocrystal elasticity constants are defined in accordance with the anisotropy ratio $A$ as specified in the following Eq. (4) with set values of $\rho$, $v_p$, and $v_s$. This ensures that, according to the Hill average, the equivalent homogeneous medium of polycrystals with transverse texture is the same for all $A$ (see Eq. (A.16) and further details in Appendix A.1).

$$ \left\{ \begin{aligned} C_{11}(A) &= \rho v_p^2 - \frac{(A-1)(A+3)}{A^2+6A+1} \rho v_s^2 \\ C_{12}(A) &= \rho v_p^2 - \frac{A^2+10A+5}{A^2+6A+1} \rho v_s^2 \\ C_{44}(A) &= \frac{4A(A+1)}{A^2+6A+1} \rho v_s^2. \end{aligned} \right. \qquad (4) $$

For each simulation, a unique microstructure is randomly generated as described in Section 2.1. The global grain elongation parameters $k_{el}$ and $\alpha_{el}$ also vary between simulations. This variability reflects the significant differences in local WLAM microstructures, which can fluctuate across locations and between parts, with elongation primarily influenced by thermal gradients and, consequently, process parameters. The overall input distribution is based on independent distributions reported in Table 1. It should be noted that $k_{el}$ ranges in [1, 4], but the uniform distribution is applied to its logarithm.

In this DOE, the only variable input quantities are the components of the elasticity tensor in the coordinate system $(x, y)$, which vary from grain to grain along with orientation angle $\phi$. In general, there are

nine components in a 2D model and some of them are equal. Setting the redundant terms aside, six components remain. Since a transverse texture is used in the present study, it is possible to further reduce the elasticity tensor to four individual components that only depend on the crystal orientation angle $\phi \in [-\pi/4, \pi/4]$ and the monocrystal elasticity constants using the following formula (see Eq. (A.6) and further details in Appendix A.1):

$$ \left\{ \begin{aligned} C'_{11} &= \left( (3 + \cos(4\phi))C_{11} + (1 - \cos(4\phi))C_{12} + 2(1 - \cos(4\phi))C_{44} \right)/4 \\ C'_{12} &= \left( (1 - \cos(4\phi))C_{11} + (3 + \cos(4\phi))C_{12} - 2(1 - \cos(4\phi))C_{44} \right)/4 \\ C'_{66} &= \left( (1 - \cos(4\phi))C_{11} - (1 - \cos(4\phi))C_{12} + 2(1 + \cos(4\phi))C_{44} \right)/4 \\ C'_{16} &= -\sin(4\phi)(C_{11} - C_{12} - 2C_{44})/4. \end{aligned} \right. \qquad (5) $$

Considering a regular discretization of sizes $N_x$ and $N_y$ along dimensions $x$ and $y$, respectively, let us denote by $C' \in \mathbb{R}^{4 \times N_x \times N_y}$ the representation of the four components defined in Eq. (5). The variable part of the input is completely represented by $C'$ for $N_x$ and $N_y$ matching the FEM resolution, thus $(N_x, N_y) = (2000, 800)$ in the present study.

Furthermore, the monocrystal elasticity constants are assumed to be global parameters, so the input can be further reduced to the spatial discretization of $\phi$ and the value of $(C_{11}, C_{12}, C_{44})$. Let us denote this simplified representation via $a \in \mathbb{R}^{N_x \times N_y \times \mathbb{R}^3}$. Further reduction is technically possible by considering variable input $A$ defined in Eq. (4), although it provides a less explicit representation of material parameters. The microstructure is generated from a Voronoi tessellation from which any discretization can theoretically be applied, with the corresponding input being denoted as $a_{VT}$. The output B-scan $b \in \mathbb{R}^{N_x \times N_y}$ can be seen as a deterministic function $M$ of input $a$:

$$ \begin{array}{l} M: \quad \mathbb{R}^{N_x \times N_y} \times \mathbb{R}^3 \rightarrow \mathbb{R}^{N_x \times N_y} \\ a \quad \mapsto \quad b \end{array} \qquad (6) $$
<figure><img src="image_7.png" /><figcaption>Fig. 4. Examples of input-output pairs in the test dataset are presented. Top of each quadrant: Input $\tilde{a}$ corresponds to anisotropy factor $A$ and field of crystal orientations $\phi(x, y)$ at 50-µm resolution, which is coarser than that of the simulations. Bottom of each quadrant: output $b$ is the simulated B-scan with dimensions $(N_r, N_t) = (1000, 81)$. Receiver positions range from $x = 7.5$ mm to $x = 17.5$ mm, hence the blank spaces for $x$ outside of this range.</figcaption></figure>

that consists in (numerically) solving the wave equation and extracting $b$. The objective of the present study is to build a metamodel of $\mathcal{M}$ defined in Eq. (6), i.e., a model that is capable of providing a good approximation of output B-scan $b$ for any input $a$ generated according to Table 1, at a much faster rate than the original forward model $\mathcal{M}$. While a fine resolution is necessary in the FEM, using a coarser resolution of $a_{VT}$ is acceptable to represent the microstructure given the large average grain size, and serves as a straightforward dimension reduction strategy. Therefore, the metamodeling problem is defined as finding a computationally-effective surrogate $\tilde{\mathcal{M}}$ expressed as follows:

$$ \tilde{\mathcal{M}}: \begin{array}{rcl} \mathbb{R}^{\tilde{N}_x \times \tilde{N}_y} & \times & \mathbb{R}^3 \\ \tilde{a} & \mapsto & \tilde{b} \end{array}, \qquad (7) $$

where $\tilde{N}_x \le N_x$, $\tilde{N}_y \le N_y$, $\tilde{a}$ and $a$ are both discretizations of the same input microstructure $a_{VT}$. It is also required that for any $a_{VT}$ drawn from Table 1, $\tilde{\mathcal{M}}(\tilde{a}) \approx \mathcal{M}(a)$, i.e. $\tilde{b} \approx b \in \mathbb{R}^{N_r \times N_t}$. In the following, unless specified otherwise, coarser input resolution $(\tilde{N}_x, \tilde{N}_y) = (500, 200)$ is used for lower computational burden, which corresponds to 50 µm in both dimensions.

Four examples of input-output pairs are shown in Fig. 4. In the first two examples, grain elongation is rather high, yet in different overall directions. In the other two examples, grains have slight elongation.

Despite these clearly visible differences in microstructural properties, there is no obvious correlation with the output B-scans. In these four examples, as in all simulations in the DOE, wave propagation from the source results in a somewhat symmetric B-scan relative to axis $x = 12.5$ mm (source location). In an isotropic medium, one would expect a perfectly symmetric B-scan indicating the propagation of two wavefronts (shown in Fig. A.3) with a faster small-amplitude positive pulse characteristic of surface-skimming waves, followed by a high-amplitude oscillation (first positive, then negative) about twice slower, characteristic of Rayleigh waves. Additional small wavefronts that are mostly parallel to the two main ones appear distinctly; however, these are expected to be less pronounced in 3D simulations as more directions are available for propagation and, therefore, scattering. In either case, compared to isotropic propagation, anisotropy and scattering are likely to induce asymmetry and additional features in the B-scan compared to the isotropic case. The higher the anisotropy ratio $A$, the more pronounced these effects become, as somewhat visible in Fig. 4. Still, identifying the influence of crystal orientations from the B-scan alone is challenging. On average, this influence decreases with depth and distance from both source and receiver.

In ML terminology, the design of metamodel $\tilde{\mathcal{M}}$ defined in (7) is a regression problem where both input and output are high-dimensional. Most classical regression and interpolation techniques are ill-suited to
<figure><img src="image_8.png" /><figcaption>Fig. 5. Comparison of simulated B-scan (GT: ground truth) in test dataset against two baseline predictors (Pred): $\hat{b}_{train}$ (top row) and low-pass filtered B-scan (bottom row). Rightmost column: average in the test dataset of point-wise absolute error $\ell_1(x, t)$.</figcaption></figure>

solve such problems. Another challenging aspect of the problem is the simulation budget: even though computational time is fairly low for each simulation, a synthetic database is quite small compared to applications such as computer vision and natural language processing. A numerical DOE usually involves the generation of two ensembles of simulations, which are respectively used to train and evaluate the model. The corresponding input-output pairs $(\tilde{a}, b = M(a))$ are referred to as *train* and *test* datasets in the following. For the test dataset of size $n_{test}$, the input samples are generated independently following input distribution in Table 1. For the train dataset of size $n_{train}$, the 3 scalar inputs $A$, $\log k_{el}$, and $\alpha_{el}$ are obtained from a randomly generated Latin hypercube sample (LHS) in $[0, 1]^3$, then they are scaled to the desired intervals via affine transform. To obtain a more space-filling LHS, a centered $L_2$-discrepancy criterion is minimized by means of random permutations of the points' coordinates (for further information regarding space-filling designs, the interested reader may refer to Damblin et al., 2013). This approach does not suit the microstructure prior to elongation, so $n_{train}$ samples are generated independently according to the last line of Table 1.

## 2.3. Evaluation metrics and metamodel baselines

Several error metrics can be used to quantitatively evaluate the approximation performance of metamodel $\tilde{M}$. The output B-scan $b$ corresponds to a vertical displacement as a function of receiver position $x$ and time $t$. For specific values of $x$ and $t$, one may consider a point-wise error defined as follows for $p \in \mathbb{N}^*$:

$$ \ell_p\{\tilde{b}, b\}(x, t) = |\tilde{b}(x, t) - b(x, t)|^p. \qquad (8) $$

A global extension of this point-wise error related to the $L_p$-norm $\|\cdot\|_p$ where B-scans are treated as vectors of length $J \times K$ can be defined for regular discretization $(x^{[J]}, t^{[k]})_{1 \le i \le J; 1 \le k \le K}$ (here, $J = N_r$ and $K = N_t$) as the average of point-wise errors:

$$ \mathcal{L}_p(\tilde{b}, b) = \left( \frac{1}{JK} \sum_{j,k} \ell_p\{\tilde{b}, b\}(x^{[J]}, t^{[k]}) \right)^{1/p} = \|\tilde{b} - b\|_p. \qquad (9) $$

Common values of $p$ are 1 and 2, which relate to mean absolute error (a.k.a. $L_1$-loss) and root mean square error, respectively. Several relative variants of the above metrics using division by some positive

value can be considered. In the present work, we use the following definition of global relative error for non-constant B-scans:

$$ \mathcal{L}_{p, \text{rel}}(\tilde{b}, b) = \frac{\mathcal{L}_p(\tilde{b}, b)}{\mathcal{L}_p(b_{\text{avg}}, b)} = \frac{\|\tilde{b} - b\|_p}{\|b_{\text{avg}} - b\|_p}, \qquad (10) $$

where $b_{avg} = \frac{1}{JK} \sum_{j,k} b(x^{[J]}, t^{[k]})$ is the average signal amplitude (over receiver position and time) of B-scan $b$. For a dataset of size $I$ comprised of simulated B-scans $(b_i)_{1 \le i \le I}$ with corresponding approximations $(\tilde{b}_i)_{1 \le i \le I}$, the (statistical) averages of the above error metrics can be used to evaluate metamodel accuracy; according to the central limit theorem, the test dataset provides a statistical estimation of their expected values over the input distribution. Some of these metrics may also serve during training as loss functions. Other evaluation metrics commonly used in computer vision (SSIM: structural similarity index) and seismology (EM: envelope misfit and PM: phase misfit, both based on continuous wavelet transform) are calculated; their specifications are provided in Appendix A.2.

To better assess approximation quality, simple predictors can be defined to serve as baselines. One of them is the average B-scan based on the set of B-scans $(b_i)_{1 \le i \le n_{train}}$ in the train dataset, defined as follows:

$$ \forall(x, t), \hat{b}_{\text{train}}(x, t) = \frac{1}{n_{\text{train}}} \sum_{i=1}^{n_{\text{train}}} b_i(x, t). \qquad (11) $$

Using $\hat{b}_{train}$ regardless of the input microstructure is simplistic yet provides fair accuracy. Another simple way to obtain a baseline with more “physical” meaning is to consider a low-pass filter. Here, a fourth-order low-pass Butterworth filter with cut frequency $f_c = 60$ MHz is used and is assumed to represent limitations of the measurement system for very high frequencies. The application of said filter to all A-scans in a given B-scan is also quite fast and provides a very good approximation.

An illustration of these two predictors is provided in Fig. 5 for one of the B-scans in the test dataset; the corresponding average $\ell_1(x, t)$ for the test dataset are also represented. The average B-scan $\hat{b}_{train}$ is much closer to that of an isotropic medium than an anisotropic polycrystalline one, and shows two notable features of surface wave propagation: surface-skimming waves (faster front) and Rayleigh waves (slower front). Using the filtered B-scans leads to predictions much closer to the ground truth. Both baselines exhibit errors that are higher in amplitude around the Rayleigh wavefronts compared to the scattered wavefronts. Filtering results in the highest difference at early times for
<table><caption>**Table 2**
Average error metrics evaluated on the test dataset for two baseline models.</caption><thead><tr><th></th><th>Low-pass filter<br>f<sub>c</sub> = 60 MHz</th><th>b&#x0302;<sub>train</sub></th></tr></thead><tbody><tr><td>L<sub>1</sub> (a.u.) &times;10<sup>3</sup></td><td>2.52</td><td>20.62</td></tr><tr><td>L<sub>1,rel</sub> (%)</td><td>3.42</td><td>28.00</td></tr><tr><td>L<sub>2,rel</sub> (%)</td><td>5.82</td><td>29.05</td></tr><tr><td>L<sub>3,rel</sub> (%)</td><td>9.03</td><td>26.72</td></tr><tr><td>1 - SSIM (%)</td><td>0.17</td><td>4.48</td></tr><tr><td>EM</td><td>0.034</td><td>0.355</td></tr><tr><td>PM</td><td>0.024</td><td>0.146</td></tr></tbody></table>

the centermost receivers, which is not the case using $\hat{b}_{train}$ as predictor. Quantitatively, several error metrics of both baselines are reported in Table 2, confirming the better accuracy of the low-pass filtering of very high frequencies.

## 2.4. Spectral layer in FNO and variants

A core element of FNO (Li et al., 2021) is the introduction of a Fourier integral operator $\mathcal{K}_{\phi}$, which can be defined for vector function $v: x \in D \subset \mathbb{R}^q \mapsto v(x) \in \mathbb{R}^N$ as follows:

$$ \forall x \in D, \mathcal{K}\{v\}(x) = \mathcal{F}^{-1}\{\mathcal{F}\{r\} \cdot \mathcal{F}\{v\}\}(x), \quad (12) $$

where $\mathcal{F}$ indicates the Fourier transform operator, $\mathcal{F}^{-1}$ is its inverse, $\mathcal{F}\{r\}$ is the Fourier transform of a periodic tensor function $r: D \to \mathbb{R}^{N \times N}$ which is implicit in practice, and $\cdot$ acts as an element-wise product with regard to the input space of its operands (here, $\mathbb{C}^q$) and as a matrix-vector product with regard to their output spaces. Operator $\mathcal{K}$ defined in Eq. (12) maps $v$ to output function $\mathcal{K}\{v\}: D \to \mathbb{R}^N$, i.e., to the same functional space. According to the convolution theorem, calculating any of the $N$ channels of $\mathcal{K}$ is equivalent to carrying out a convolution product with each channel of $v$ (each time with a different real-valued function of $D$), then summing over the $N$ convolution products. Implementations of such operators often involve Fast Fourier Transform (FFT) which assumes a regular discretization of $D$. Such implementation is referred to as a spectral layer in the following. Due to the discrete nature of data in practice, a finite amount of Fourier coefficients expressed via complex-valued tensor $R$ related to $\mathcal{F}\{r\}$ are instantiated and compose the trainable parameters of the spectral layer. To handle finer discretizations of $D$, a spectral layer usually ignores the coefficients of high-order modes not accounted for via $R$, thus behaving as an operator to some extent. This also allows to filter out high-order modes that may not provide useful information.

The functioning of a spectral layer in dimension $q=2$ is represented in Fig. 6. In this representation, the dot product involving $R$ in Eq. (12) is decomposed into two steps: $N$ element-wise tensor multiplications, then summation over the $N$ blocks. A spectral layer is often parameterized by the number of non-negative modes kept in each dimension, $m_i = \lceil K_i/2 \rceil$ where $K_i$ is the total number of low-order modes involved in the $i$-th dimension. For real-valued data, the Fourier transform is Hermitian, which can be enforced by using “real” variants of the FFT that ignore the negative modes for one of the mode coordinates $k_i$, thus dealing with data and parameters in Fourier space with about half the represented size.

In the present study, four notable adaptations of the aforementioned spectral layer are proposed and investigated.

1. The number of input and output channels can technically differ. This represents a straightforward extension where function $r$ maps to $\mathbb{R}^{N_{in} \times N_{out}}$ instead of $\mathbb{R}^{N \times N}$.

2. Use of factorized Fourier transform operator, which relies on $q$ separate 1D Fourier transforms instead of a $q$D transform. The $q$ separate branches are summed at the end of the layer. This approach is a core element of F-FNO proposed in Tran et al. (2023) and involves much fewer parameters than its not-factorized counterpart.

3. A channel-wise variant where $r$ maps to $\mathbb{R}^{N_{in}}$ instead of $\mathbb{R}^{N_{in} \times N_{out}}$ is proposed. The tensor dot product in Eq. (12) is replaced by a single element-wise tensor multiplication followed by $N_{out}$ linear combinations of the $N_{in}$ channels via parameter matrix $\Gamma \in \mathbb{C}^{N_{out} \times N_{in}}$. This reduces the number of parameters by almost a factor $N_{out}$ in most cases.

4. Addition of a bias term per channel via vector $\mathbf{B} \in \mathbb{R}^{N_{out}}$ applied to the Fourier coefficients before application of the IFFT except those set to 0.

<table><caption>**Table 3**
Number of parameters for spectral layers whether Fourier factorization and/or channel-wise layers are implemented. $K_i$ is the total number of low-order modes requested by the user. In practical implementations, the actual amount of Fourier coefficients is about half the indicated value due to the use of real-valued input/output data, hence the 1/2 factors. Addition of a bias term per channel would involve $N_{out}$ supplementary real parameters.</caption><thead><tr><th>Neural network layer</th><th>Number of complex parameters</th></tr></thead><tbody><tr><td>Spectral</td><td>$\frac{1}{2} N_{in} N_{out} \prod_{i=1}^{q} K_i$</td></tr><tr><td>Spectral channel-wise</td><td>$N_{in} (N_{out} + \frac{1}{2} \prod_{i=1}^{q} K_i)$</td></tr><tr><td>Factorized spectral</td><td>$\frac{1}{2} N_{in} N_{out} \sum_{i=1}^{q} K_i$</td></tr><tr><td>Factorized spectral channel-wise</td><td>$N_{in} (q N_{out} + \frac{1}{2} \sum_{i=1}^{q} K_i)$</td></tr></tbody></table>

Combining these adaptations leads to a modified spectral layer represented in Fig. 7. The amount of parameters in a spectral layer whether any of these adaptations is included or not is reported in Table 3. The idea of using channel-wise operations to reduce the number of network parameters is not new in DL. For instance, Ji et al. (2021) used channel-wise fully-connected (a.k.a. dense) layers for inversion of ground processing radar data. To the authors' knowledge, similar variants have rarely been proposed for FNO-related models; akin to channel-wise operations is the use of group convolution for FNO proposed in Kim et al. (2024). As for the use of combinations of 1D convolutions in DNN instead of multi-dimensional convolutions, it is one of the many improvements proposed for the architecture of DNN Inception (Szegedy et al., 2016) to reduce the number of network parameters.

Contrary to most applications of FNO, our input and output fields have differently defined domains. Here, the input is a function of spatial coordinates $x$ and $y$, whereas the output B-scan is a function of $x$ and time $t$. Also, the range of spatial coordinate $x$ is not the same between input and output. For the latter, $x$ corresponds to specific receiver positions. Raw application of FNO would conserve the input and output dimensions, which makes little sense in this context. Several strategies can be proposed to handle the data size difference between input and output via neural network layers. For instance, interpolations can either increase or reduce data size, while pooling layers can reduce it. Among layers with trainable parameters, convolutions with stride and transpose convolutions can respectively decrease and increase data size. More specific to FNO is the strategy proposed for UNO (Rahman et al., 2023): size is increased via zero-padding of higher-order modes in Fourier space, whereas discarding null high-order coefficients reduces data size. For F-FNO variants, this strategy is well-suited for size increase but ill-suited for size decrease. Here, nearest neighbor interpolations are used to increase size.

In the present work, the spectral layers are implemented with fixed output size. If the input size is different from the output's, the amount of zero-coefficients is adapted before application of the IFFT to return the desired output size. This is a notable difference from the usual approach represented in Figs. (6, 7) where Fourier coefficients are handled so that the output size matches that of the input. Due to the use of interpolations, this subtle difference does not come into play in the present study.

The number of non-negative modes involved in a spectral layer is set to the maximum value dimension-wise: $m_i = \lceil N_i/2 \rceil$ where $N_i$ is the corresponding output size in dimension $i$. The so-called Fourier layers of
<figure><img src="image_9.png" /><figcaption>Fig. 6. Typical spectral layer for 2D data with usual “image-like” format (C, H, W). Here, C = N for both input and output. The number of modes involved in dimensions relative to x₁ and x₂ are K₁ and K₂, respectively. A blank cell indicates that the (complex) coefficient in Fourier space is set to 0. Usually, this allows output to have dimensions matching those of the input, i.e., (H, W) = (N₂, N₁).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 7. Adaptation of the spectral layer represented in Fig. 6 including different number of channels for input (N<sub>in</sub>) and output (N<sub>out</sub>), factorized Fourier transform, channel-wise variant, and bias term. The number of modes involved in dimensions relative to x<sub>1</sub> and x<sub>2</sub> are K<sub>1</sub> and K<sub>2</sub>, respectively. A blank cell indicates that the (complex) coefficient in Fourier space is set to 0.</figcaption></figure>

FNO and its variants often involve residual connection, linear operators,
and an activation layer in addition to a spectral layer. Unless specified
otherwise, a spectral layer is simply followed by an activation function
in the present work, which is another notable difference from these
Fourier layers.

2.5. Deep neural network architectures

The first layer of our neural network surrogates is a pre-processing
that applies Eq. (5) to input microstructure ã, thus yielding tensor C'
of shape (4, Ñ_y, Ñ_x). The elasticity constants are expressed in GPa,
whereas the angles are expressed in radians (NB: shown in degrees in
the figures). A scaling factor S₀ is then applied. Here, parameter S₀ =
0.005 GPa⁻<fnref n="1" /> is fixed (not trainable), which leads to scaled coefficients
ranging roughly from -0.2 to 1.6. This can be seen as a form of input
normalization.

The main DNN architecture used in the present work is represented
in Fig. 8 which emphasizes the data shape along the network layers.
The DNN is a feed forward network with the aforementioned pre-
processing layer followed by three groups of layers. It is referred to as
a mixed factorized-Fourier encoder-decoder (MFFED) in the following be-
cause it involves typical convolutional encoder and decoder groups sep-
arated by an intermediate group mostly composed of factorized spectral
layers (as represented in Fig. 7). Overall, the architecture is similar to
that of a U-Net without skip-connections, involving 1 195 133 trainable
parameters.

* The first group acts as an encoder increasing the number of layers from 4 to 128 and decreasing the size along the other two dimensions from 200 × 500 to 25 × 63. It involves successions of 2D convolutions (Conv2D) with 5 × 5 kernels, parametric rectified linear unit (PReLU) activations, and 2 × 2 average pooling. There are three such blocks followed by two {Conv2D, PReLU} blocks. This group is composed of 680 357 trainable parameters, which includes the convolution kernels and five negative slope coefficients

(a single one per PReLU layer). A deterministic algorithm is used
for all Conv2D operations to ensure reproducibility.<fnref n="1" />

* The presence of an intermediate group is the main difference com-
pared to usual U-Net-like architectures. There are three factorized
spectral layers as represented in Fig. 7 (i.e., channel-wise and with
bias) using all the available modes in both dimensions. Each of
these layers is followed by a PReLU activation and two interpola-
tion layers are used in-between to progressively increase the data
size, except for the number of channels that remains fixed. Near-
est neighbor interpolation (2D) is used, which duplicates certain
rows and columns but not all of them since the size increases by
a factor lower than 2. This interpolation method was preferred to
bilinear interpolation for which reproducibility is not guaranteed
during gradient backpropagation using CUDA software, unless a
fixed algorithm is used<fnref n="2" /> (resulting in higher computational time).
The spectral layers involve about 50 000 complex parameters.
Counting the number of trainable float parameters (both real
and imaginary parts) plus the three PReLU coefficients add up
to 269 451 trainable parameters for this group.

* The third and final group acts as a decoder modifying data shape from 128 × 63 × 81 to that of the B-scans. It consists of four blocks each composed of nearest neighbor interpolation, 5 × 5 Conv2D, and PReLU activation. These four blocks progressively reduce the number of channels while increasing size along the height, leading to data with shape 8 × 1000 × 81. Finally, a 5 × 5 Conv2D kernel is applied to yield a single channel, and no activation is used thereafter. This group is composed of 272 341 trainable parameters that include the convolution kernels and four PReLU coefficients.

Six alternatives to MFFED with a different intermediate group rep-
resented in Fig. A.1 are considered in the present study.

<footnote n="1">¹ see for instance the PyTorch documentation regarding deterministic algorithms, last consulted on 10.23.2024.</footnote>

<footnote n="2">² see footnote 1.</footnote>
<figure><img src="image_11.png" /><figcaption>Fig. 8. Schematic representation of the MFFED architecture, which is the main neural network surrogate model used in the present study. Numbers indicate the data shape (excluding batch size) as they are processed by the network. For a given data block, the adjacent numbers indicate channels, height and width from top to bottom, with bold numbers highlighting shape changes.</figcaption></figure>

* MFFED_2: same as MFFED but without the channel-wise variant and with intermediate number of channels equal to 32 instead of 128 (except for input and output of the group that remain at 128 channels).

* MFED: same as MFFED_2, but where spectral layers are channel-wise and not factorized.

* MFED_2: same as MFED, but without the channel-wise variant and an intermediate number of channels equal to 4.

* CED: compared to the previous models, each spectral layer is replaced by a succession of three separable Conv2D layers with kernel sizes 1 × 21, 21 × 1, and 1 × 1 where the first-two are channel-wise while the third is not. Zero-padding is used to preserve data shape.

* MDED: compared to the previous models, each spectral layer is replaced by a channel-wise dense layer with PReLU activation. The data are transformed (for all 128 channels) from shape 25 × 63 to 2 × 3, then 2 × 3, and finally 63 × 81.

* ED: a single 2D nearest neighbor interpolation layer modifies the shape to 128 × 125 × 81 (no intermediate group).

Some models have a lower number of either channels or width/height compared to MFFED in order to reduce the number of trainable parameters and limit overfitting due to overparameterization. The main performance evaluation criteria to compare these metamodels are approximation error (cf. Section 2.3) and computational time, other factors such as memory usage are also of interest.
## 2.6. Training strategy and implementation details

The FEM elastodynamics simulations rely on a C++ implementation. Aside from that, most of the methods described previously are implemented in Python with common scientific computing and ML libraries (numpy, sciPy, pandas, etc.), notably PyTorch for DNN design, training, and evaluation.

Regarding the DOE, train and test datasets have respective sizes $n_{train} = 10000$ and $n_{test} = 1000$. During training, the train dataset is randomly split in (sub-)train and validation datasets using a 9:1 ratio. The validation dataset is involved in learning rate scheduling and helps to assess convergence during training, whereas the test dataset is used to evaluate a model only once training is completed. The loss function used to fit the DNN parameters is $L_1$ defined in Eq. (9) that is minimized using an Adam optimizer (Kingma and Ba, 2015) with a batch size of 5, no regularization, and usual parameters: notably $10^{-3}$ for learning rate and (0.9, 0.999) for betas. Additionally, a "reduce on plateau" scheduler is applied to the learning rate, which is decreased by a factor of 2 if the loss does not decrease by at least 1% over 5 consecutive epochs. The optimizer is set to run up to a maximum of 100 epochs and may stop before then if the learning rate gets under $5 \times 10^{-5}$ (i.e., after five reductions); this is expected to reduce the time required for training with little effect on metamodel accuracy. The DNN weight values saved after training are those obtained at the end of the epoch with lowest loss on the validation dataset. The same strategy is applied to all the aforementioned DNN architectures.

GPU-acceleration is used to speed up DNN computations, and is only turned off for comparison with simulation runtimes (FEM implementation is CPU-only) on the same hardware. Performance of the trained models are evaluated using an Intel Core i7-12800H CPU that may be accelerated via a NVIDIA RTX A1000 (4GB) GPU. Training relies on higher-performance hardware, namely a NVIDIA RTX A5000 (24GB) GPU with an Intel Xeon Gold 5220R CPU.

## 3. Results

In Section 3.1, the training and performance of the neural network surrogates is presented and examined. Then, Section 3.2 focuses on metamodel MFFED, more specifically on its predictions and improvements over baseline $\hat{b}_{train}$. Influence of model architecture and training procedure on performance is analyzed in Section 3.3.

### 3.1. General comparison between the different neural network surrogates

The evolution of the loss function $L_1$ during training is represented in Fig. 9 for six out of the seven neural networks mentioned in Section 2.5. Convergence failed for ED ($L_1$ remained around 0.06 for both datasets), so the training curves are not represented. In contrast, the other models converged and training was completed in 40 to 50 epochs, which required about 1 to 3.5 h. Several error metrics other than the loss were calculated during training, they show similar evolution. It is likely that ED is not sufficiently expressive to provide a good approximation of the simulated B-scans, this issue is overcome by the addition of a few intermediate layers in the other models. Additional information regarding training of these models is reported in Table 4 together with various performance metrics. For all models except ED, $L_1$ in the validation dataset at the end of the first epoch is lower than that of baseline $\hat{b}_{train}$ in the test dataset. In the subsequent epochs, these models further outperform $\hat{b}_{train}$ but do not reach the accuracy of the low-pass filter baseline.

Metamodels in Table 4 are sorted by increasing $L_{1,rel}$ but the rankings would be the same using any other error metric, except for CED and MFED_2 whose values are very close. In terms of approximation error, the best model is MFFED, although it is also the one with highest computational time. The models with factorized spectral layers (MFFED and MFED_2) obtained better approximation than their non-factorized counterparts.

In this regard, the channel-wise variant seems beneficial only when factorized layers are involved.

Depending on the amount of trainable parameters in the intermediate layers, two groups of models can be identified. Models MFFED, MFED, and CED (resp. MFFED_2, MFED_2, and MDED) involved in the left (resp. right) half of Fig. 9 have an additional number of trainable parameters lower than 500 000 (resp. higher than 1 000 000). Focusing on the latter group, it seems that the higher the number of parameters, the higher the gap between train and validation errors, and the worse the approximation error. However, this does not apply to the other group, suggesting that the quality of approximation is influenced by both layer type and size.

In metamodeling, one is typically interested in the reduction in computational time of a forward call compared to the original model. However, with neural networks, the gradient of any real-valued function derived from the output can be calculated via backpropagation at low additional cost. For instance, the partial derivatives of the $L_p$-error between predicted and observed B-scans can be a relevant quantity for gradient-based inversion and optimization techniques, but their computational cost would be prohibitive without such surrogates. Therefore, the runtimes for computing $\nabla_g \|\hat{b}\|_1$ (the gradient of the $L_1$-norm of the surrogate B-scans relative to the metamodel inputs) are provided in the last two rows of Table 4.

The models with fewer parameters are not necessarily faster to compute: CED and MFFED have intermediate layers involving larger data sizes, hence the higher runtimes; whereas runtime is at most twice that of ED for the others. A forward (multi-CPU) call of $M$ takes approximately 30 s, which amounts to 3~4 days to generate the whole training database using one machine. Using the same hardware, a forward and backward call of MFFED is about a hundred times shorter, even though the original model is quite fast. The gain is even higher for the other networks, reaching up to about 400 (excluding ED). GPU-acceleration is straightforward for these models and allows to further reduce computational time by a factor of 4~5 using a moderately powerful GPU, amounting to about 65 ms for MFFED.

In terms of accuracy and computational time, models CED and MFED are outperformed on both criteria by models MFFED_2 and MDED, respectively. The choice between models MFFED, MFFED_2, MFED_2, and MDED amounts to a trade-off between these two criteria, where MFFED is the most accurate of the four models, whereas MDED is the fastest.

To investigate the influence of inputs on metamodel outputs, $|\nabla_a \|\hat{b}\|_1|$ averaged over all the instances in the test dataset is represented in Fig. 10 for all metamodels except ED. For each scalar input, the corresponding coordinate is an estimate of its average first-order contribution to the output, referred to as influence in the following, for short. Strictly speaking, it is neither a global sensitivity nor a global feature importance measure, only a proxy thereof. Calculation of the partial derivative relative to $A$ is straightforward using the derivatives relative to each elastic constant and the chain rule, which involves the derivatives of Eq. (4). All six models yield a similar value for the influence of $A$, which is higher than that of any point-wise input crystal orientation.

For all models, crystal orientation has higher influence close to the surface, especially under the receiver positions, then tends to decrease with the distance from the closest receiver. Also, crystal orientation influence is mostly symmetric relative to axis $x = 12.5$ mm. Most likely, these properties also apply to the original model. However, non-null values are obtained at locations far from both source and receivers (notably in the bottom corners) whereas they have no influence on the simulated B-scans in a 3-µs time frame. Also, some non-smooth patterns such as spots and horizontal or vertical bands are visible, mostly close to the domain boundaries, although they are unlikely to apply to $M$. Interestingly, MFFED is the metamodel with best accuracy while featuring a high-influence area of crystal orientations that is arguably more concentrated around the receivers and somewhat smoother compared to the other metamodels.
<figure><img src="image_12.png" /><figcaption>**Fig. 9.** Training curves of the neural networks. For each model, the train dataset is divided in subsets “train” and “validation” of size 9000 and 1000, respectively. For a given model/color, the vertical segment indicates the best epoch for which the loss on the validation dataset is lowest. The circular markers in the curves indicate epochs after which the learning rate was reduced by a factor of 2.</figcaption></figure>

<caption>**Table 4**</caption>

Properties, training summary, and performance of neural networks surrogates. Performance is evaluated using the parameter values obtained at the end of the “best” epoch, i.e., with lowest validation loss.

<table><thead><tr><th></th><th>MFFED</th><th>MFFED_2</th><th>CED</th><th>MFED_2</th><th>MDED</th><th>MFED</th><th>ED</th></tr></thead><tbody><tr><td>L<sub>1</sub> (a.u.) ×10<sup>3</sup></td><td>6.69</td><td>7.44</td><td>7.86</td><td>7.92</td><td>8.01</td><td>8.37</td><td>54.63</td></tr><tr><td>L<sub>1,rel</sub> (%)</td><td>9.09</td><td>10.10</td><td>10.68</td><td>10.76</td><td>10.88</td><td>11.65</td><td>74.14</td></tr><tr><td>L<sub>2,rel</sub> (%)</td><td>8.24</td><td>9.13</td><td>9.71</td><td>9.70</td><td>9.87</td><td>10.26</td><td>99.20</td></tr><tr><td>L<sub>3,rel</sub> (%)</td><td>6.86</td><td>7.59</td><td>8.12</td><td>8.07</td><td>8.27</td><td>8.54</td><td>106.91</td></tr><tr><td>1 − SSIM (%)</td><td>0.37</td><td>0.45</td><td>0.51</td><td>0.51</td><td>0.53</td><td>0.57</td><td>79.16</td></tr><tr><td>EM</td><td>0.125</td><td>0.135</td><td>0.144</td><td>0.142</td><td>0.143</td><td>0.149</td><td>86.25</td></tr><tr><td>PM</td><td>0.050</td><td>0.055</td><td>0.059</td><td>0.059</td><td>0.059</td><td>0.061</td><td>54.00</td></tr><tr><td></td><td>MFFED</td><td>MFFED_2</td><td>CED</td><td>MFED_2</td><td>MDED</td><td>MFED</td><td>ED</td></tr><tr><td>Additional parameters<sup>a</sup></td><td>242 435</td><td>1 085 827</td><td>65 667</td><td>3 511 691</td><td>5 788 035</td><td>483 651</td><td>0</td></tr><tr><td>Epochs (best/last)</td><td>44/47</td><td>42/46</td><td>42/46</td><td>27/41</td><td>20/40</td><td>31/45</td><td>3/27</td></tr><tr><td>Epoch runtime (s)</td><td>258</td><td>131</td><td>159</td><td>111</td><td>101</td><td>119</td><td>96</td></tr><tr><td>Training runtime (hours h minutes)</td><td>3 h 20</td><td>1 h 40</td><td>2 h 00</td><td>1 h 15</td><td>1 h 10</td><td>1 h 30</td><td>0 h 45</td></tr><tr><td>V<sub>a</sub>||b||<sub>1</sub> CPU runtime<sup>b</sup> (ms)</td><td>304</td><td>137</td><td>172</td><td>83</td><td>77</td><td>103</td><td>66</td></tr><tr><td>V<sub>a</sub>||b||<sub>1</sub> GPU runtime<sup>b</sup> (ms)</td><td>65</td><td>26</td><td>37</td><td>20</td><td>15</td><td>23</td><td>13</td></tr></tbody></table>

## 3.2. Detailed analysis of MFFED predictions

To represent the predictive capabilities of MFFED, three members in the test dataset tagged as low, avg, and high are identified. As can be seen in Fig. 11, the tag names are consistent with prediction quality of both MFFED and $\hat{b}_{train}$; among the test dataset, low is one of the predictions with lowest error, high has one of the highest, and avg is in between. Naming is also indicative of the anisotropy ratio A. Predicted B-scans using MFFED for the three tagged members are represented in Fig. 12 and compared to their simulated counterparts. The two main wavefronts are well predicted, although the prediction error tends to be higher in areas corresponding to the Rayleigh wavefronts. Fine details in the simulated B-scans are not retained: the error patterns suggest that part of the contributions of the scattered wavefronts are not accounted for. The effects of anisotropy are mostly captured by the metamodel but

the marks of the scattered wavefronts appear blurred with the small-amplitude ones missing. For both models, the error tends to increase with anisotropy. This is expected for $\hat{b}_{train}$ because scattered waves tend to have higher amplitudes in the B-scan as A increases, thus increasing the difference from the smooth baseline B-scan. Error increase with A applies not only to B-scan error metrics but also to A-scan error metrics such as EM and PM when they are calculated individually for each receiver (see Fig. A.2 in Appendix A.3).

To examine the error distribution in the B-scan, the average in the test dataset of point-wise error $\ell_1(x, t)$ of both models is shown in Fig. 13 together with the corresponding increase factor (for instance, a value of 0.7 at a given point means that the average error of MFFED at said point is 30% lower than that of the baseline). This value is lower for MFFED at most points and is otherwise negligible. The latter occurs before the arrival of the first waves: signal amplitude is null for all simulations and for $\hat{b}_{train}$ so error is null for the baseline. A small

<footnote><sup>a</sup> Trainable parameters in the intermediate group of network layers. Contribution of the encoder and decoder groups amounts to 952 698 parameters.</footnote>

<footnote><sup>b</sup> Gradient is calculated for all instances in the test dataset sequentially, the average runtime is reported first using CPU only, then using GPU-acceleration.</footnote>
<figure><img src="image_13.png" /><figcaption>**Fig. 10.** Average absolute value of partial derivative relative to all scalar inputs for different neural networks surrogates. In the legend, · refers to the function mapping the input to the L₁-norm of the estimated B-scan.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Fig. 11.** *L*₁,*rel* error on all instances of the test dataset for predictors *ŷ*<sub>train</sub> and MFFED as a function of input anisotropy ratio *A*. Tags *low*, *avg*, and *high* refer to specific dataset members.</figcaption></figure>

regularization coefficient (10⁻⁸) is added to the denominator of the increase factor to avoid division by zero, but its value remains very high in that case. Elsewhere, the error of MFFED is higher around the arrival of the Rayleigh waves but is also reduced the most relative to the baseline (up to ~85% lower). For higher values of *t*, the average error is

reduced by about 30% to 60%; this improvement tends to decrease with time and is smaller for the centermost receivers. This latter observation is tough to explain, no outlier prediction with considerably high error for the centermost receivers was found in the test dataset. The scattered waves are predicted fairly well using MFFED but, compared to the
<figure><img src="image_15.png" /><figcaption>**Fig. 12.** Examples of predicted B-scans using metamodel MFFED (Pred) for given input and comparison with corresponding simulated B-scans (GT: ground truth). Tags low, avg, and high refer to specific dataset members. Moreover, both input and output of instance avg correspond to example (1) shown in Fig. 4.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Fig. 13.** Average $\ell_1(x, t)$ error on test dataset for predictors $\hat{b}_{\text{train}}$ (left) and MFFED (middle). Right: ratio of point-wise average error of MFFED by that of $\hat{b}_{\text{train}}$. NB: Error ratio is much higher than 2 in the red areas due to negligible denominator.</figcaption></figure>

two main wavefronts that are highly influenced by subsurface crystal orientations, scattering may occur further from the surface and is arguably more challenging to emulate.

To further study the influence of input $A$ on the predictions, this parameter has been modified while keeping fixed crystal orientations; the corresponding simulations/predictions applied to instance avg are shown in Fig. 14. In the simulations, increasing $A$ leads to B-scans that increasingly deviate from the isotropic case ($A = 1$) with more pronounced effects of scattering. This is also the case for the metamodel predictions up to $A = 2.5$, even though the B-scans are blurry, and the error between the two increases with $A$. This is expected for $A \in$

[1.5, 2.5] since the input is still in the training domain, which is not the case outside of this interval. For $A \in$ [1, 1.5], the predictions remain consistent with their simulated counterpart, and the error is low overall. For $A > 2.5$, however, the predictions are not consistent with the physics anymore, resulting in extremely poor predictions. The same behavior was observed for other instances in the test dataset: acceptable B-scans can still be predicted up to some threshold in [2.5, 3.5]. This behavior was not checked for all instances in the test dataset as this would require a considerable amount of additional simulations. Also, a given neural network architecture trained following the same procedure (but
<figure><img src="image_17.png" /><figcaption>**Fig. 14.** Influence of varying anisotropy ratio *A* at fixed input crystal orientations of instance avg in the test dataset (Pred) against its simulated counterpart (GT). For each line, the four fields in gray scale represent the value of the first layer in the network at given *A*, obtained from Eq. (5). The $C'_i$ have different extrema each and range from black (lowest) to white (highest). The rest of the legend, omitted for readability, is consistent with previous figures such as Fig. 12.</figcaption></figure>

with different RNG) can also lead to a different threshold value. For $A \in [1, 1.5]$, the predictions are always very good.

### 3.3. Influence of model architecture and training procedure on performance

Several modifications to metamodel architectures, training strategy, and hyperparameters have been attempted. Rather than detailing the results of these experiments, the following provides a summary of these results.

First, the spatial resolution of the metamodel input data could be higher: it was 12.5 µm in the FEM simulations, but only 50 µm for the neural networks. Networks with architectures similar to those proposed in Section 2.5 were tailored to the FEM-resolution, typically with more parameters per layer and/or a few more layers in the encoder and intermediate groups of layers. This led to metamodels that were

more computationally expensive, but training did not result in better approximation error in the test dataset. The latter was almost equal to that of the best models or lower. A possible explanation is that a coarser resolution is sufficient to capture the main features of the polycrystalline microstructures: there are about 100 pixels in most grains and the boundaries where scattering occurs are still estimated fairly well.

Reducing the number of parameters in the intermediate group of layers either through number of channels (MDED) or height and width (networks other than ED and MDED) can cause failure during optimization (increase in the loss function at the start or during early epochs) or result in higher validation error. Conversely, increasing the number of parameters and/or the number of layers (mostly attempted in the intermediate group) tends to increase the gap between train and validation errors on both sides, leading to worse metamodels as
they are more computationally expensive and overfitting, with lower approximation quality. Compared to the other metamodels, MFFED has higher computational time and the memory burden is not negligible due to the larger size of data running through the intermediate layers of the network. Reducing overfitting for all networks except ED could be beneficial, but common strategies such as dropout layers, batch normalization layers, and additional regularization terms in the loss function did not result in better approximation error: in some cases, the train loss increased but the test error did not decrease.

Regarding fully convolutional models such as CED, optimization failed to converge when square “full” convolution kernels were used in the intermediate convolution layers; this issue was solved by using large separable kernels instead. This seems consistent with the better performance of models with factorized spectral layers that operate separately on width and height (MFFED and MFFED_2) over their counterparts without factorization (MFED and MFED_2).

The training process seems to be quite robust for the proposed network architectures in terms of overall approximation quality. Runs of the optimization procedure with different RNG resulted in similar global error metrics on the test dataset (up to a difference of a few percents). Different values of initial learning rate and batch size did not lead to better results. Actually, a very small batch size (5, here) seems preferable, possibly due to the relatively small size of the training dataset. This is consistent with the results presented in Lehmann et al. (2024) where an F-FNO model applied to seismic wave propagation surrogate modeling provided good approximation quality using a batch size of 8 or 16 during training.

# 4. Discussion

## 4.1. Interpretation and implications

In the training domain (input distribution in Table 1), the networks are able to approximately capture scattered wavefronts and their dependence with A, so they must be sensitive to differences in the $C'_{ij}$ of adjacent grains, which lead to scattering. Regarding the relevance of calculating the actual $C'_{ij}$ via a pre-processing layer, predictions of MFFED were found to be much closer to the ground truth compared to a slight variant of MFFED lacking said layer (cf. Fig. A.3), more details in Appendix A.3. This noteworthy result suggests that models such as MFFED may be trained on a database restricted to individuals with higher anisotropy and still perform well when anisotropy is lower, especially when they include the pre-processing layer. For instance, a training database restricted to $A \in [2, 2.5]$ with the same size may well lead to better approximation quality in this range with little decrease in approximation quality for lower values of A, however this hypothesis warrants further investigation.

Study of the different variants indicates that factorized spectral layers lead to better accuracy than their not-factorized counterparts. A similar result was observed in Lehmann et al. (2024), where F-FNO models obtained better accuracy compared to FNO and other variants such as UNO. In the present study, a notable advantage of MFFED over the other architectures is that it seems less prone to overfitting (a lower train-validation gap was obtained) yet complex enough to provide good approximation quality despite the relatively low number of parameters. In MFFED_2, the channel-wise variant is not used, which increases the number of parameters, yet the model obtains good accuracy while being faster to compute than MFFED thanks to the use of fewer channels. Also, models with lower approximation error on the train dataset have potential to become more accurate if trained on a larger training database, i.e., the trade-off between accuracy and computational time depends on $n_{train}$.

Spectral layers (more specifically, Fourier layers) are most commonly employed in FNO frameworks where both input and output share the same definition domain. The present study demonstrates that spectral layers can be effective outside of this framework and

can be integrated in convolutional encoder-decoder architectures. This approach is compatible with different definition domains for the input and the output, but forgoes the flexibility in discretizing the definition domain, which is a noteworthy feature of FNO models. To reduce computational burden, we recommend applying spectral layers at a low resolution or after a few pooling layers. Factorized spectral layers, in particular, appear to be promising alternatives to convolutional layers at the bottleneck of the network, potentially improving accuracy. Using channel-wise variants and adjusting the number of channels provides flexibility in terms of the number of parameters. A comparison between models MFFED_2 and CED suggests that reducing the number of channels, in particular, may lead to faster computations while maintaining good accuracy.

## 4.2. Applicability to online LU testing

Regarding the application of MFFED in practice, the computational time using moderate-performance GPU is comparable to the pulse repetition rate of conventional nanosecond lasers used as LU excitation source (20 Hz), although higher rates (~1 kHz) are technically possible depending on the equipment, the desired pulse energy, and acquisition time requirements. Despite this potential, it is difficult to determine whether these metamodels are actually efficient enough for online LU inspection. Inversion techniques may require a high amount of metamodel calls, so the computational time might still be limiting but, unlike the original model, not prohibitive.

Several elements that significantly influence LU-induced displacements in WLAM components were accounted for in the present study, namely: anisotropy, local crystal orientations of specific grains (notably those closest to the source), which can be large and elongated, and scattering in polycrystals. The setting is similar to that of Singh et al. (2022), in that simulations of UT transmission in polycrystals with transverse textures are used to train DNNs. However, it differs in its objective: while Singh et al. (2022) focuses on inversion with UT data from several source locations, our study concentrates on forward metamodeling based on a single source location. Additionally, the present study accounts for uncertainty regarding material elasticity constants through Eq. (4), and considers crystal orientations varying in $[-\pi/4, \pi/4]$, rather than just half of the interval. This allowed for a detailed study of the influence of anisotropy on forward metamodeling, but actual microstructure variability is only partially represented. Notably, a fiber texture with a preferred orientation following the elongation direction is often observed in MAM, which does not align well with a transverse texture. In practice, three Euler angles representing crystal orientation should be considered as input of the metamodel, instead of a single angle. For consistency, Eq. (A.19) should be used instead of Eq. (4) to relate elasticity constants with A. It is worth noting that the magnitude of anisotropy is rather straightforward for crystals with cubic symmetry and can be summarized using a single quantity; however, this is not the case for other systems. If a 2D model were used without further simplifications, six elasticity components would need to be considered in the pre-processing layer, instead of the four provided in Eq. (5), with different formulas consistent with Eq. (A.2). In a 3D model, there would be up to 21 elastic constants to consider.

Actual IN718-based WLAM microstructures have variable grain size with potentially higher aspect ratios, higher anisotropy at high temperatures, rough surface, etc., these elements are not accounted for in the present study. Microstructures with smaller grains may require smaller input resolutions to ensure that enough pixels per grain are present, which will increase computational time. Considering the relatively large grain size implied by our input distribution in Table 1, an input resolution of 50 µm seems sufficient, and the approximation quality is arguably mostly limited by the relatively small size of the training database (10 000). However, microstructures with smaller grains are technically out-of-distribution and our study of the extrapolation capabilities of MFFED showed that the prediction can be either excellent
(for $A \in [1, 1.5]$) or extremely poor (e.g. for $A = 3.0$). The influence of varying grain size is out of the scope of the present study, and later studies should consider input distribution that cover a range of grain sizes as large as possible, as good generalization capabilities are not guaranteed out-of-distribution.

Regarding LU techniques, a thermoelastic (non-destructive) regime is generally used in practice, if possible. Such regime should be driven by a boundary condition on shear stress $\sigma_{xy}$ instead of normal stress in Eq. (3). In practice, it is difficult to provide a detailed equivalent source term because it might be relevant to use both normal and shear components, with their proportions dependent on the energy of the laser pulse. For further discussion on this topic, please refer to Section 5.10 in *Scruby and Drain* (1990).

The proxy for input influence shown in Fig. 10 suggests that it may be possible to estimate the crystal orientations a few millimeters below the source and receivers. While some limitations apply to the proxy for input influence shown in Fig. 10, the choice of more appropriate techniques is not trivial. Indeed, input independence is a common assumption in global sensitivity analysis (Saltelli et al., 2008), which enables a unique decomposition of variance, allowing for the estimation of Sobol' indices (Sobol', 2001) based on a high number of metamodel calls. However, the family of input microstructures implies that crystal orientations at different locations are not independent. While Sobol' indices could technically be estimated using random “white noise” crystal orientation fields within reasonable computational time, this would not yield valuable information. Other sensitivity measures that do not require input independence could be considered, but measures such as mutual information, which focus on dissimilarity between joint and conditional distributions, would necessitate conditional and/or joint density estimation of the output for all inputs, making it too time-consuming. Additionally, a Hoeffding-like decomposition would still be lacking. Derivative-based sensitivity measures (Sobol' and Kucherenko, 2009) allow for ranking input influence and are computationally efficient in this context, especially since DNN are used. Simpler measures, such as correlation coefficients, could also be considered. Regardless of the sensitivity measure, selecting a meaningful real-valued function based on an output B-scan is challenging. The $L_1$-norm was chosen due to its simplicity, but other functions could prove to be more relevant.

Despite the simplifications, one may attempt to use the forward metamodel for inversion. However, during gradient descent, the input crystal orientations are unlikely to remain piece-wise constant unless a strong prior on the microstructure is made. Furthermore, first-order derivative information may not be sufficient since scattering at grain boundaries depends on crystal orientations on both sides. In a nutshell, the topology and high dimensionality of the input contribute to the ill-posedness of the inverse problem, making it unlikely to recover a full “actual” microstructure based on a surrogate of the forward model. Still, the performance of the metamodels is acceptable and the proposed architectures are promising for inversion tasks with lower target dimension such as defect/anomaly detection, localization, and sizing. Same for online LU testing, but actual application requires further efforts, most notably experimental data.

## 5. Conclusion

The present study investigates the use of simulated data in the perspective of online LU testing of additive manufacturing components. Simulation of elastic wave propagation in 2D polycrystals with a LU excitation source relies on a FEM model. A distribution of input microstructures is assigned to this forward model, accounting for variable crystal orientations, uncertain anisotropy and large elongated grains typical in WLAM. The output of interest is a B-scan expressing surface out-of-plane displacements around the source at early times. Despite the short computational time of this forward 2D model (30 s per forward call), it is prohibitive in high-frequency online inspection. We propose to address this bottleneck using a metamodel, whose design

is a challenging regression problem with high-dimensional input and output. We find that deep learning models thoughtfully designed and trained for this purpose allow for both fast and accurate approximation of anisotropy and scattering effects.

Several neural networks are trained based on a fairly small database (10 000 training samples). The architectures include fully convolutional encoder–decoder elements, their novelty lies in the inclusion of spectral layers, strongly inspired from the (factorized) Fourier layers at the core of FNO and F-FNO models. Several variants of said layers are proposed and investigated. Among the trained models, the lowest approximation error is obtained using model MFFED that relies on a channel-wise factorized variant of the spectral layers, characterized by a considerably lower number of parameters compared to other variants. Despite this apparent simplification, said layers still imply a high number of operations and large data running through the network, making MFFED the slowest of all trained models. Still, MFFED is quite fast: it is able to run both a forward and backward call about 100 times faster than one execution of the forward model on the same hardware (CPU-only). Using modest-performance GPU, speed-up can be higher and runtime of one MFFED call (forward and backward) takes about 65 ms.

MFFED is also able to capture some of the effects of anisotropy on wave propagation, not only in-distribution but also in out-of-domain cases where anisotropy is weaker. In the surrogate B-scans, the two main surface wavefronts are predicted quite well (the average error is reduced by 85% compared to a simple baseline) and the features of waves scattered at grain boundaries are also accounted for, but the improvement is not as high (average error is reduced by 30 ~ 60% compared to the baseline). Approximation error of MFFED generally increases with anisotropy, in which case metamodeling is arguably more challenging. The extrapolation capacity of MFFED for weaker anisotropy suggests that it should be trained focusing on cases with higher anisotropy, but further investigation is required on this matter.

A novel aspect of this work lies in the inclusion of spectral layers at the bottleneck of a CNN encoder–decoder architecture. The different variants of these layers allow for versatility in neural network design, notably with regards to the number of trainable parameters and the size of the data passing through the network. The present case study shows that the trade-off between computational cost and metamodel accuracy is not straightforward. Still, a few insights were gained from this case study for neural network design in a similar context. First, using the factorized variants of the spectral layers appears to lead to reduced overfitting and better approximation. Also, to reduce computational burden, we recommend using the spectral layers with low resolution data, which also applies after a few pooling layers if high-resolution input data are involved.

Even though the proposed metamodels show potential for online LU inspection of WLAM components, further work is required for practical application. In practice, the surface of WLAM beads is likely to be rough and measurement noise will lead to discrepancies compared to either simulated or predicted B-scans. In this regard, other discrepancies may be caused by the limits of the forward model itself, considering a 2D domain involves simplifications that neglect 3D effects. Likewise, electromagnetic, optical, and thermal phenomena are only accounted for via a simplified elastic source term, which should be adapted to a thermoelastic LU regime. Moreover, temperature will be quite influential on mechanical properties in an online context, but this aspect was neglected in the present study. There is also room for improvement in terms of modeling the WLAM microstructures by using more realistic crystal orientations (correlation between elongation and one of the crystal axes instead of the unrealistic transverse texture), variable grain size and aspect ratios, which are usually lower for the bottom layers and higher for the top-most layers during MAM processes. To that end, experimental EBSD data can be useful and could be completed by synthetic data generated using numerical models (see for instance Camus et al., 2023). The method proposed in the present study is a valuable
starting point for future work on forward metamodeling and inversion focusing on more realistic data.

Regarding inversion problems, which are arguably of higher interest, yet significantly more challenging, promising DNN-based approaches in related fields such as electromagnetic and seismic wave propagation have been proposed (e.g., Ji et al., 2021; Zhu et al., 2023). The physics-embedded framework proposed in Guo et al. (2023) regroups several strategies and techniques that can make use of simulated and/or observed data. An interesting perspective is to include physical loss terms related to the model's equations, commonly used in physics-informed neural networks (PINN), which would require all components of the displacement vector and first-order derivatives (boundary conditions) or second-order derivatives (elastic wave equations) given our forward model. Perhaps most important for actual ultrasound testing is to draw the best from simulated and experimental data. Other approaches based on DNN such as domain adaptation (Pyle et al., 2022), multi-fidelity training (Granados et al., 2023; Nerlikar et al., 2024), and generative neural networks (Gatti and Clouteau, 2020; Granados et al., 2024) recently showcased potential in this regard.

### CRediT authorship contribution statement

**Frédéric Allaire:** Writing – review & editing, Writing – original draft, Visualization, Validation, Software, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. **Filippo Gatti:** Writing – review & editing, Supervision, Methodology, Conceptualization. **Roberto Miorelli:** Writing – review & editing, Supervision, Methodology, Conceptualization. **Vincent Dorval:** Writing – review & editing, Software, Methodology, Conceptualization. **Jérôme Laurent:** Writing – review & editing, Supervision, Project administration, Methodology, Funding acquisition, Conceptualization.

### Declaration of competing interest

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

### Acknowledgments

This work is part of the COLUMBO project (grant ANR-21-CE08-0026) funded by the French National Agency for Research (ANR). The authors would like to thank Alexandre Imperiale for the implementation of the FEM code used in the present study and Fanny Lehmann for fruitful discussion on F-FNO models and available code on GitHub.

## Appendix

### A.1. Anisotropic elastodynamics equations, cubic symmetry, and texture

For polycrystals, let us assume that a grain is characterized by homogeneous elasticity tensor $c$. In this case, Eq. (1) restricted to said grain can be simplified as follows:

$$ \forall i \in [1, d], \rho \frac{\partial^2 u_i}{\partial t^2} - \sum_{j,k,l} c_{kil} \frac{\partial^2 u_j}{\partial x_k \partial x_j} = f_i. \quad (A.1) $$

For another grain, crystal orientation may be different, same for $c$. Continuity equations for displacement and stress should be accounted for to solve Eq. (1) in a polycrystal. In the present study, the FEM scheme ensures continuity at element boundaries. Voigt notation is used to relate the 4th-order elasticity tensor components to those of its 2nd-order counterpart $C = (C_{ij})_{1 \le i,j \le 6}$ with the following correspondence between indices: $1 \leftrightarrow 11, 2 \leftrightarrow 22, 3 \leftrightarrow 33, 4 \leftrightarrow 23, 5 \leftrightarrow 13, 6 \leftrightarrow 12$. In 2D, index 3 (resp. indices 3, 4, and 5) is omitted in the 4th-order (resp. 2nd-order) representation. For further theoretical background regarding crystallography, the interested reader may refer to Kocks et al. (2000).

Let $P = (P_{ij})_{i,j}$ denote an orthogonal transformation matrix from so-called local coordinate system to so-called reference coordinate system, where $c^{\text{loc}} = (c_{ijkl}^{\text{loc}})_{i,j,k,l}$ and $c^{\text{ref}} = (c_{ijkl}^{\text{ref}})_{i,j,k,l}$ are the respective 4th-order elasticity tensor in said coordinate systems. The following relationship applies between these two tensors:

$$ c_{mnop}^{\text{ref}} = \sum_{i,j,k,l} P_{m,l} P_{n,j} P_{o,k} P_{p,l} c_{ijkl}^{\text{loc}}. \quad (A.2) $$

Let us consider a material with cubic symmetry and 2nd-order elasticity tensor $C^{\text{loc}} = (C_{ij}^{\text{loc}})_{1 \le i,j \le 6}$ in local coordinate system matching crystal directions $\langle 100 \rangle$ (i.e., normal to a cube face). $C^{\text{loc}}$ can be written using only three elasticity constants $C_{11}, C_{12}$, and $C_{44}$ as follows:

$$ C^{\text{loc}} = \begin{pmatrix} C_{11} & C_{12} & C_{12} & 0 & 0 & 0 \\ C_{12} & C_{11} & C_{12} & 0 & 0 & 0 \\ C_{12} & C_{12} & C_{11} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{44} & 0 \\ 0 & 0 & 0 & 0 & 0 & C_{44} \end{pmatrix}. \quad (A.3) $$

Considering that the crystal coordinate system is rotated from the reference by rotation by an angle $\phi$ about the $x_3^{\text{loc}}$-axis, the corresponding transformation matrix is the inverse rotation matrix, i.e.,

$$ P = \begin{pmatrix} \cos(\phi) & \sin(\phi) & 0 \\ -\sin(\phi) & \cos(\phi) & 0 \\ 0 & 0 & 1 \end{pmatrix} \quad (A.4) $$

for which the resulting 2nd-order elasticity tensor $C^{\text{ref}} = (C_{ij}^{\text{ref}})_{1 \le i,j \le 6}$ in reference coordinate system reads

$$ C^{\text{ref}} = \begin{pmatrix} C_{11}^{\text{ref}} & C_{12}^{\text{ref}} & C_{12} & 0 & 0 & C_{16}^{\text{ref}} \\ C_{12}^{\text{ref}} & C_{11}^{\text{ref}} & C_{12} & 0 & 0 & -C_{16}^{\text{ref}} \\ C_{12} & C_{12} & C_{11} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{44} & 0 \\ C_{16}^{\text{ref}} & -C_{16}^{\text{ref}} & 0 & 0 & 0 & C_{66}^{\text{ref}} \end{pmatrix}. \quad (A.5) $$

$$ \text{where} \quad \begin{cases} C_{11}^{\text{ref}} = \left( \frac{(3 + \cos(4\phi))C_{11} + (1 - \cos(4\phi))C_{12} + 2(1 - \cos(4\phi))C_{44}}{4} \right) \\ C_{12}^{\text{ref}} = \left( \frac{(1 - \cos(4\phi))C_{11} + (3 + \cos(4\phi))C_{12} - 2(1 - \cos(4\phi))C_{44}}{4} \right) \\ C_{66}^{\text{ref}} = \left( \frac{(1 - \cos(4\phi))C_{11} - (1 - \cos(4\phi))C_{12} + 2(1 + \cos(4\phi))C_{44}}{4} \right) \\ C_{16}^{\text{ref}} = -\sin(4\phi)(C_{11} - C_{12} - 2C_{44})/4. \end{cases} \quad (A.6) $$

In the present study, polycrystalline materials with monocrystal cubic symmetry are considered and all crystals have one $\langle 100 \rangle$ direction imposed along the $x_3$-axis. In this specific case, a counterclockwise rotation through an angle $\phi$ about the $x_3$-axis is sufficient to determine the orientation of the remaining cube faces grain-wise. The general term of the elasticity tensor in global coordinate system is therefore given by Eqs. (A.5) and (A.6). We use the term transverse texture to describe such a material with no preferred orientation other than $x_3$, i.e. such that $\phi \sim U(-\pi/4, \pi/4)$.

For polycrystals, one may use an averaging scheme to calculate an equivalent homogeneous material whose elasticity tensor is denoted by $C^X = (C_{ij}^X)_{1 \le i,j \le 6}$. For a polycrystal with transverse texture, common averaging schemes result in the following:

$$ C^X = \begin{pmatrix} C_{11}^X & C_{12}^X & C_{12} & 0 & 0 & C_{16}^X \\ C_{12}^X & C_{11}^X & C_{12} & 0 & 0 & -C_{16}^X \\ C_{12} & C_{12} & C_{11} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{44} & 0 \\ C_{16}^X & -C_{16}^X & 0 & 0 & 0 & C_{66}^X \end{pmatrix}. \quad (A.7) $$
Using the Voigt average (Voigt, 1928) (X=V) yields

$$ egin{cases} C_{11}^V = \frac{C_{11}^3 + C_{12}^2 + 2C_{44}}{4} \\ C_{12}^V = \frac{C_{11}^3 + 3C_{12}^2 - 2C_{44}}{4} \\ C_{66}^V = \frac{C_{11}^3 - C_{12}^2 + 2C_{44}}{4} \\ C_{16}^V = 0. \end{cases} \quad (A.8) $$

The Reuss average (Reuss, 1929) (X=R) focuses on the compliance tensor and results in the following for a transverse texture:

$$ egin{cases} C_{11}^R = \frac{C_{11}^2 - C_{12}^2 + 6C_{11}C_{44} - 2C_{12}C_{44}}{2(C_{11} - C_{12} + 2C_{44})} \\ C_{12}^R = \frac{C_{11}^2 - C_{12}^2 + 2C_{11}C_{44} + 2C_{12}C_{44}}{2(C_{11} - C_{12} + 2C_{44})} \\ C_{66}^R = \frac{2C_{11}C_{44} - 2C_{12}C_{44}}{C_{11} - C_{12} + 2C_{44}} \\ C_{16}^R = 0. \end{cases} \quad (A.9) $$

One may notice that $C_{12}^X = C_{11}^X - 2C_{66}^X$ for both Voigt and Reuss averages, thus the equivalent homogeneous material exhibits transverse isotropy. More precisely, in the plane normal to the $x_3$-axis, wave propagation is isotropic. The Hill average (Hill, 1952) (X=H) is a commonly-used homogenization scheme that consists in using the arithmetic average of Voigt's and Reuss's schemes and also leads to a transverse isotropic material in the present case.

Note that these results are different when applied to untextured polycrystals (i.e., all crystal orientations are possible and have equal probability density), which is a common hypothesis when applying Voigt, Reuss, or Hill averages. For an untextured polycrystal (still assuming monocrystal cubic symmetry), all three schemes result in average elasticity tensor denoted by $C^{0,X} := (C_{ij}^{0,X})_{1 \le i,j \le 6}$ with the following form:

$$ C^{0,X} = \begin{pmatrix} C_{11}^{0,X} & C_{12}^{0,X} & C_{12}^{0,X} & 0 & 0 & 0 \\ C_{12}^{0,X} & C_{11}^{0,X} & C_{12}^{0,X} & 0 & 0 & 0 \\ C_{12}^{0,X} & C_{12}^{0,X} & C_{11}^{0,X} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{44}^{0,X} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{44}^{0,X} & 0 \\ 0 & 0 & 0 & 0 & 0 & C_{44}^{0,X} \end{pmatrix}, \quad (A.10) $$

$$ \text{where } C_{12}^{0,X} = C_{11}^{0,X} - 2C_{44}^{0,X}, \quad (A.11) $$

i.e., the equivalent material is isotropic. For untextured polycrystals, the Voigt average yields

$$ egin{cases} C_{11}^{0,V} = \frac{3C_{11} + 2C_{12} + 4C_{44}}{5} \\ C_{44}^{0,V} = \frac{C_{11} - C_{12} + 3C_{44}}{5} \end{cases}, \quad (A.12) $$

whereas the Reuss average yields

$$ egin{cases} C_{11}^{0,R} = \frac{C_{11}^3 - 2C_{12}^2 + 2C_{11}C_{12} + 8C_{11}C_{44} - 4C_{12}C_{44}}{3(C_{11} - C_{12}) + 4C_{44}} \\ C_{44}^{0,R} = \frac{5C_{11}C_{44} - 5C_{12}C_{44}}{3(C_{11} - C_{12}) + 4C_{44}} \end{cases} \quad (A.13) $$

A common measure of anisotropy for a cubic monocrystal is Zener's cubic anisotropy ratio $A := 2C_{44}/(C_{11} - C_{12})$. For polycrystalline alloys, even though the monocrystal is anisotropic, the alloy is often defined via reference compression and shear wave speeds $v_p$ and $v_s$, assuming isotropic behavior at large scales. Consistent with the assumption of an overall equivalent isotropic material defined via elasticity constants $C_{11}^0$ and $C_{44}^0$, one may write $v_p = \sqrt{C_{11}^0}/\rho$ and $v_s = \sqrt{C_{44}^0}/\rho$. Using Eqs. (A.8) and (A.9) related to transverse texture according to $C_{11}^0 = C_{11}^X$ and $C_{44}^0 = C_{66}^X$, it is possible to express monocrystal (cubic) elasticity coefficients as a function of $A$, $v_p$, and $v_s$ consistently with the Voigt average:

$$ egin{cases} C_{11}^{A,V} = \rho v_p^2 - \frac{A-1}{A+1} \rho v_s^2 \\ C_{12}^{A,V} = \rho v_p^2 - \frac{A+3}{A+1} \rho v_s^2 \\ C_{44}^{A,V} = \frac{2A}{A+1} \rho v_s^2, \end{cases} \quad (A.14) $$

consistently with the Reuss average:

$$ egin{cases} C_{11}^{A,R} = \rho v_p^2 - \frac{A-1}{2A} \rho v_s^2 \\ C_{12}^{A,R} = \rho v_p^2 - \frac{3A+1}{2A} \rho v_s^2 \\ C_{44}^{A,R} = \frac{A+1}{2} \rho v_s^2, \end{cases} \quad (A.15) $$

and consistently with the Hill average:

$$ egin{cases} C_{11}^{A,H} = \rho v_p^2 - \frac{(A-1)(A+3)}{A^2+6A+1} \rho v_s^2 \\ C_{12}^{A,H} = \rho v_p^2 - \frac{A^2+10A+5}{A^2+6A+1} \rho v_s^2 \\ C_{44}^{A,H} = \frac{4A(A+1)}{A^2+6A+1} \rho v_s^2. \end{cases} \quad (A.16) $$

Monocrystal elasticity coefficients with values $C_{11}^{A,X}$, $C_{12}^{A,X}$, and $C_{44}^{A,X}$ defined in Eqs. (A.14)–(A.16) can be used to study the influence of the anisotropy degree (as expressed by $A$) for a polycrystal with transverse texture while keeping the same equivalent homogeneous equivalent medium (transverse isotropic, here) characterized by $\rho$, $v_p$, and $v_s$. Another advantage of this approach is for uncertainty quantification of cubic elasticity constants of polycrystalline alloys. Indeed, these can be somewhat uncertain (e.g., several sets of values available in the scientific literature) even though $\rho$, $v_p$, and $v_s$ are often known with relatively high certainty. In that case, one may attribute most of the uncertainty to $A$.

Lastly, in practice, it is often assumed that polycrystalline media are untextured. In that case, for reference, one may apply the same approach to Eqs. (A.12) and (A.13) according to $C_{11}^0 = C_{11}^{0,X}$ and $C_{44}^0 = C_{44}^{0,X}$, which yields instead

$$ egin{cases} C_{11}^{A,V} = \rho v_p^2 - \frac{4(A-1)}{3A+2} \rho v_s^2 \\ C_{12}^{A,V} = \rho v_p^2 - \frac{4A+6}{3A+2} \rho v_s^2 \\ C_{44}^{A,V} = \frac{5A}{3A+2} \rho v_s^2, \end{cases} \quad (A.17) $$

$$ egin{cases} C_{11}^{A,R} = \rho v_p^2 - \frac{4(A-1)}{5A} \rho v_s^2 \\ C_{12}^{A,R} = \rho v_p^2 - \frac{8A+2}{5A} \rho v_s^2 \\ C_{44}^{A,R} = \frac{2A+3}{5} \rho v_s^2, \end{cases} \quad (A.18) $$

$$ \text{and} \quad \left\{ \begin{aligned} C_{11}^{A,H} &= \rho v_p^2 - \frac{4(A-1)(A+4)}{3A^2+19A+3} \rho v_s^2 \\ C_{12}^{A,H} &= \rho v_p^2 - \frac{4A^2+32A+14}{3A^2+19A+3} \rho v_s^2 \\ C_{44}^{A,H} &= \frac{5A(2A+3)}{3A^2+19A+3} \rho v_s^2. \end{aligned} \right. \quad (A.19) $$

## A.2. Details of some evaluation metrics employed

Accuracy metrics defined in Section 2.3 are defined based on the $L_p$-norm, and may be complemented by other (dis)similarity metrics commonly used in other fields.

Due to the two-dimensional structure of B-scans, the structural similarity index (SSIM) commonly used in computer vision can be computed using the following formula:

$$ SSIM(\tilde{b}, b) = \frac{(2\tilde{b}_{avg}b_{avg} + c_1)(2Cov + c_2)}{(\tilde{b}_{avg}^2 + b_{avg}^2 + c_1)(\tilde{V} + V + c_2)}, \quad (A.20) $$

where $c_1 = 1 \times 10^{-4}$, $c_2 = 9 \times 10^{-4}$,

$$ b_{avg} = \frac{1}{JK} \sum_{j,k} b(x^{[j]}, t^{[k]}) \quad (A.21) $$

is the average signal amplitude of B-scan $b$,

$$ V = \frac{1}{JK} \sum_{j,k} (b(x^{[j]}, t^{[k]}) - b_{avg})^2 \quad (A.22) $$
<figure><img src="image_18.png" /><figcaption>Fig. A.1. Intermediate group of layers for the different neural network architectures. All networks have the same architecture for the encoder and decoder groups of layers.</figcaption></figure>

is the variance of amplitude in B-scan b, and

$$ \text{Cov} = \frac{1}{JK} \sum_{j,k} (b(x^{[j]}, t^{[k]}) - b_{\text{avg}})(\tilde{b}(x^{[j]}, t^{[k]}) - \tilde{b}_{\text{avg}}) \quad (\text{A.23}) $$

is the covariance between amplitudes of B-scans b and $\tilde{b}$. For convenience, its complement to 1 (i.e., 1 − SSIM) is provided in the paper so that only “error” metrics are used.

To compare two sets of time series, several comparison metrics based on continuous wavelet transform (CWT) were introduced, notably in seismology which involves displacement-related time series (Kristeková et al., 2009). Given a so-called mother wavelet function $\psi_0 : \eta \in \mathbb{R} \mapsto \psi_0(\eta) \in \mathbb{C}$ and time series $s : t \mapsto s(t) \in \mathbb{R}$, one may replace the scale and location coefficients in the CWT to yield a wavelet-based time-frequency transform. The coefficient for time series $s$ at time $t$ and frequency $f \in \mathbb{R}$ can be generally expressed as follows:

$$ W(s)(t, f) = \sqrt{\beta|f|} \int_{-\infty}^{+\infty} s(\zeta)\psi_0^*((\zeta - t)\beta f) d\zeta, \quad (\text{A.24}) $$

where $\beta > 0$ can be adapted depending on the mother wavelet, and $\psi^*$ denotes the complex conjugate of $\psi \in \mathbb{C}$. Following guidelines in Torrence and Compo (1998), wavelets similar to the studied signal are recommended for CWT-based wavelet analysis, which commonly involves Morlet wavelets. Here, 2nd-order Paul wavelets are similar to the source signal. Their mother wavelet function can be defined as

follows:

$$ \psi_0(\eta) = -\sqrt{\frac{8}{3\pi}(1 - i\eta)^{-3}}, \quad (\text{A.25}) $$

where $i^2 = -1$. For this mother wavelet, $\beta = 5/(4\pi)$ is a suitable value in Eq. (A.24). Then, a set of wavelet coefficients are calculated at L specific frequencies $\{f^{[1]}, \dots, f^{[L]}\}$, usually with logarithmic spacing. In the present study, 32 frequencies ranging from $f^{[1]} = 5$ MHz to $f^{[32]} = 50$ MHz are used. Misfit coefficients focusing on the envelope and phase of the signals are defined based on a weighted average of their respective coefficients defined in Eq. (A.24). For time-series $s$ and $\tilde{s}$ with notations $W(t, f) = W(s)(t, f)$ and $\tilde{W}(t, f) = W(\tilde{s})(t, f)$, the normalized envelope misfit (EM) and phase misfit (PM) defined in Kristeková et al. (2009) respectively write:

$$ \text{EM}\{\tilde{s}, s\} = \sqrt{\frac{\sum_{k,l} (|\tilde{W}(t^{[k]}, f^{[l]})| - |W(t^{[k]}, f^{[l]})|)^2}{\sum_{k,l} |W(t^{[k]}, f^{[l]})|^2}}, \quad (\text{A.26}) $$

$$ \text{and PM}\{\tilde{s}, s\} = \frac{1}{\pi} \sqrt{\frac{\sum_{k,l} (|W(t^{[k]}, f^{[l]})| \Delta \arg(t^{[k]}, f^{[l]}))^2}{\sum_{k,l} |W(t^{[k]}, f^{[l]})|^2}}, \quad (\text{A.27}) $$

where $\Delta \arg(t, f)$ is the difference (in radians) between the arguments of $W(t, f)$ and $\tilde{W}(t, f)$. In the present work, the formula for envelope
<figure><img src="image_19.png" /><figcaption>**Fig. A.2.** (EM, PM) calculated for simulated A-scans against predicted counterparts for all 81 receivers and all 1000 instances in the test dataset using two predictors: baseline $\hat{b}_{train}$ (top) and metamodel MFFED (bottom). Right plots are zoomed-in representations of the left ones. Each point is colored according to the value of input A. Due to the high number of points, kernel density estimation is carried out. Five iso-density contours are shown (except in the bottom-left plot) encompassing 25, 50, 75, 90, and 95 % of the 81000 points.</figcaption></figure>

misfit in Eq. (A.26) can be applied to compare two A-scans related to receiver at position x by evaluating $EM\{\tilde{b}, b\}(x) = EM\{\tilde{b}(x, .), b(x, .)\}$. Spatial averaging is used to evaluate $EM\{\tilde{b}, b\} = \frac{1}{J} \sum_j EM\{\tilde{b}, b\}(x^{[j]})$ as a B-scan misfit score. The same applies to PM defined in Eq. (A.27).

### A.3. Complements regarding neural network surrogates

Specific parts of the neural network architectures described in Section 2.5 are represented in Fig. A.1. The representation of the intermediate group of model MFFED is the same as that of Fig. 8. The other networks have the same layers in the encoder and decoder groups, so only the layers in the intermediate group (green data blocks) are represented.

Envelope misfit (EM) and phase misfit (PM) defined in Eqs. (A.26) and (A.27) are two wavelet-based scores used to compare time series, so they can be calculated to compare predicted and simulated A-scans of all receivers individually, rather than being averaged to provide an error metric between two B-scans. Application to baseline predictor $\hat{b}_{train}$ and metamodel MFFED leads to the plots in Fig. A.2. Same as error metrics on B-scans, EM and PM on A-scans are lower using MFFED compared to the baseline, roughly by a factor of 3. For both predictors, the smallest misfits are obtained for low values of A, tend to increase with A, which was also the case for other metrics. Both scores are strongly influenced by the norm of the simulated A-scan: the higher this norm, the lower the scores, which is expected given the definition of EM and PM. The time series tend to fluctuate increasingly with A, so the norm of the simulated A-scans is positively correlated with A, hence the influence of A on EM and PM.

The extrapolation capabilities of MFFED for $A \in [1, 1.5]$ are discussed in Section 3.2. The isotropic case $A = 1$ is further investigated by comparing the quality of predictions using MFFED, baseline $\hat{b}_{train}$,

and a variant of MFFED. The latter involves a modification of the first pre-processing layer of the network: instead of computing the $C'_{ij}$ (4 channels) from crystal orientations and scalars $C_{11}, C_{12}, C_{44}$ via Eq. (5), the first channel of “MFFED modified” is the crystal orientation angles (in radians) and the other three channels are constant fields each equal to one of the $C'_{ij}$ in GPa divided by pre-processing factor $S_0$ to obtain values scaled around $\pm 1$. Both original and modified networks are trained following the same procedure. In the test dataset, the approximation error of both models is roughly the same. The main difference occurs out of the training domain for $A \in [1, 1.5]$, in which case MFFED performs better than its variant, as shown in Fig. A.3. The B-scan predicted by MFFED is very close to the simulated B-scan with only a slight error regarding the Rayleigh wavefront. For MFFED modified, when the crystal orientations of instance avg are used as input, the error is higher in this area and elsewhere because some traces of scattering remain in the prediction. Prediction by MFFED is also much better than that of baseline $\hat{b}_{train}$ for which the surface-skimming and Rayleigh wavefronts are more spread out than in the isotropic B-scan due to the variations of wave speeds with anisotropy and crystal orientations in the train dataset. To conclude, using Eq. (5) has two advantages: 1) it results in better predictions in the out-of-domain cases with lower anisotropy and guarantees the first network layer to be unique in the isotropic case, whereas the output undesirably depends on the input crystal orientation for MFFED when $A = 1$.

## Data availability

Data will be made available on request.
<figure><img src="image_20.png" /><figcaption>Fig. A.3. Simulated output B-scan with isotropic (A = 1) input microstructure (GT: “ground truth”, top right) against predictions using several metamodels (Pred, middle column) and difference (Pred–GT to the right). Corresponding input fields of the first NN layer are shown for predictors MFFED and its *modified* variant that involves stacking the scalar inputs $C_{ij}$ as three additional channels to crystal orientations (here, those of instance avg in the test dataset). The $C_{ij}$ and $C'_{ij}$ each range from black (lowest) to white (highest). The rest of the legend, omitted for readability, is consistent with previous figures such as Fig. 12 and Fig. 14.</figcaption></figure>

## References

Abuabiah, M., Mbodj, N.G., Shaqour, B., Herzallah, L., Juaidi, A., Abdallah, R., Plapper, P., 2023. Advancements in laser wire-feed metal additive manufacturing: A brief review. Materials 16 (5), 2030, https://doi.org/10.3390/ma16052030.

Arguelles, A.P., Kube, C.M., Hu, P., Turner, J.A., 2016. Mode-converted ultrasonic scattering in polycrystals with elongated grains. J. Acoust. Soc. Am. 140 (3), 1570–1580, https://doi.org/10.1121/1.4962161.

Arguelles, A.P., Turner, J.A., 2017. Ultrasonic attenuation of polycrystalline materials with a distribution of grain sizes. J. Acoust. Soc. Am. 141 (6), 4347–4353, https://doi.org/10.1121/1.4984290.

Bai, X., Tie, B., Schmitt, J.-H., Aubry, D., 2020. Comparison of ultrasonic attenuation within two- and three-dimensional polycrystalline media. Ultrasonics 100, 105980, https://doi.org/10.1016/j.ultras.2019.105980.

Bakre, C., Afzalimir, S.H., Jamieson, C., Nassar, A., Reutzel, E.W., Lissenden, C.J., 2022. Laser generated broadband Rayleigh waveform evolution for metal additive manufacturing process monitoring. Appl. Sci. 12 (23), 12208, https://doi.org/10.3390/app122312208.

Cai, Y., Song, Y., Ni, P., Liu, X., Li, X., 2021. Subwavelength ultrasonic imaging using a deep convolutional neural network trained on structural noise. Ultrasonics 117, 106552, https://doi.org/10.1016/j.ultras.2021.106552.

Cai, Y., Xiong, J., Chen, H., Zhang, G., 2023. A review of in-situ monitoring and process control system in metal-based laser additive manufacturing. J. Manuf. Syst. 70, 309–326, https://doi.org/10.1016/j.jmsy.2023.07.018.

Camus, T., Maisonnette, D., Baulin, O., Senninger, O., Guillemot, G., Gandin, C.-A., 2023. Three-dimensional modeling of solidification grain structures generated by laser powder bed fusion. Materialia 30, 101804, https://doi.org/10.1016/j.mtla.2023.101804.

Cantero-Chinchilla, S., Simpson, C.A., Ballisat, A., Croxford, A.J., Wilcox, P.D., 2023. Convolutional neural networks for ultrasound corrosion profile time series regression. NDT E Int. 133, 102756, https://doi.org/10.1016/j.ndteint.2022.102756.

Cantero-Chinchilla, S., Wilcox, P.D., Croxford, A.J., 2022. Deep learning in automated ultrasonic NDE - Developments, axioms and opportunities. NDT E Int. 131, 102703, https://doi.org/10.1016/j.ndteint.2022.102703.

Damblin, G., Couplet, M., Iooss, B., 2013. Numerical studies of space-filling designs: optimization of Latin Hypercube Samples and subprojection properties. J. Simul. 7 (4), 276–289, https://doi.org/10.1057/jos.2013.16.

Dorval, V., Leymarie, N., Imperiale, A., Demaldent, E., Lhuillier, P.-É., 2025. Numerical estimation of ultrasonic phase velocity and attenuation for longitudinal and shear waves in polycrystalline materials. Ultrasonics 148, 107517, https://doi.org/10.1016/j.ultras.2024.107517.

Dryburgh, P., Li, W., Pieris, D., Fuentes-Domínguez, R., Patel, R., Smith, R.J., Clark, M., 2022. Measurement of the single crystal elasticity matrix of polycrystalline materials. Acta Mater. 225, 117551, https://doi.org/10.1016/j.actamat.2021.117551.

Every, A.G., Utegulov, Z.N., Veres, I.A., 2013. Laser thermoelastic generation in metals above the melt threshold. J. Appl. Phys. 114 (20), 203508, https://doi.org/10.1063/1.4832483.
Gantala, T., Balasubramaniam, K., 2022. DPAI: A Data-driven simulation-assisted-Physics learned AI model for transient ultrasonic wave propagation. Ultrasonics 121, 106667, https://doi.org/10.1016/j.ultras.2021.106671.

Gasteau, D., Chigarev, N., Ducouso-Ganjehi, L., Gusev, V.E., Jenson, F., Calmon, P., Tournat, V., 2016. Single crystal elastic constants evaluated with surface acoustic waves generated and detected by lasers within polycrystalline steel samples. J. Appl. Phys. 119 (4), 043103, https://doi.org/10.1063/1.4940367.

Gatti, F., Clouteau, D., 2020. Towards blending physics-based numerical simulations and seismic databases using generative adversarial network. Comput. Methods Appl. Mech. Engrg. 372, 113421, https://doi.org/10.1016/j.cma.2020.113421.

Ghanani, N., Pasebani, S., 2024. A review on wire-laser directed energy deposition: Parameter control, process stability, and future research paths. J. Manuf. Mater. Process. 8 (2), 84, https://doi.org/10.3390/jmmp8020084.

Grabec, T., Veres, I.A., Ryzy, M., 2022. Surface acoustic wave attenuation in poly-crystals: Numerical modeling using a statistical digital twin of an actual sample. Ultrasonics 119, 106585, https://doi.org/10.1016/j.ultras.2021.106585.

Granados, G.E., Gatti, F., Miorelli, R., Robert, S., Clouteau, D., 2024. Generative domain-adapted adversarial auto-encoder model for enhanced ultrasonic imaging applications. NDT E Int. 148, 103234, https://doi.org/10.1016/j.ndteint.2024.103234.

Granados, G., Miorelli, R., Gatti, F., Robert, S., Clouteau, D., 2023. Towards a multi-fidelity deep learning framework for a fast and realistic generation of ultrasonic multi-modal Total Focusing Method images in complex geometries. NDT & E Int. 139, 102906, https://doi.org/10.1016/j.ndteint.2023.102906.

Guo, S., Agarwal, M., Cooper, C., Tian, Q., Gao, R.X., Guo, W., Guo, Y., 2022. Machine learning for metal additive manufacturing: Towards a physics-informed data-driven paradigm. J. Manuf. Syst. 62, 145-163, https://doi.org/10.1016/j.jmsy.2021.11.003.

Guo, R., Huang, T., Li, M., Zhang, H., Eldar, Y.C., 2023. Physics-embedded machine learning for electromagnetic data imaging: Examining three types of data-driven imaging methods. IEEE Signal Process. Mag. 40 (2), 18-31, https://doi.org/10.1109/MSP.2022.3198805.

Harke, K.J., Calta, N., Tringe, J., Stobbe, D., 2022. Laser-based ultrasound interrogation of surface and sub-surface features in advanced manufacturing materials. Sci. Rep. 12, 3309, https://doi.org/10.1038/s41598-022-07261-w.

He, F., Yuan, L., Mu, H., Ros, M., Ding, D., Pan, Z., Li, H., 2023. Research and application of artificial intelligence techniques for wire arc additive manufacturing: a state-of-the-art review. Robot. Comput.-Integr. Manuf. 82, 102525, https://doi.org/10.1016/j.rcim.2023.102525.

Hemmasian, A., Ogoke, F., Akbari, P., Malen, J., Beuth, J., Barati Farimani, A., 2023. Surrogate modeling of melt pool temperature field using deep learning. Addit. Manuf. Lett. 5, 100123, https://doi.org/10.1016/j.addlet.2023.100123.

Hill, R., 1952. The elastic behaviour of a crystalline aggregate. Proc. Phys. Soc. Sect. A 65 (5), 349, https://doi.org/10.1088/0370-1298/65/5/307.

Honarvar, F., Varvani-Farahani, A., 2020. A review of ultrasonic testing applications in additive manufacturing: Defect evaluation, material characterization, and process control. Ultrasonics 108, 106227, https://doi.org/10.1016/j.ultras.2020.106227.

Huang, M., Rokhlin, S.I., Lowe, M.J.S., 2022. Appraising scattering theories for polycrystals of any symmetry using finite elements. Phil. Trans. R. Soc. A 380 (20210382), https://doi.org/10.1098/rsta.2021.0382.

Huang, M., Sha, G., Huthwaite, P., Rokhlin, S.I., Lowe, M.J.S., 2021. Longitudinal wave attenuation in polycrystals with elongated grains: 3D numerical and analytical modeling. J. Acoust. Soc. Am. 149 (4), 2377-2394, https://doi.org/10.1121/10.003955.

Jarnnikar, N.D., Liu, S., Brice, C., Zhang, X., 2022. In situ microstructure property prediction by modeling molten pool-quality relations for wire-feed laser additive manufacturing. J. Manuf. Process. 79, 803-814, https://doi.org/10.1016/j.jmapro.2022.05.013.

Ji, Y., Zhang, F., Wang, J., Wang, Z., Jiang, P., Liu, H., Sui, Q., 2021. Deep neural network-based permittivity inversions for ground penetrating radar data. IEEE Sensors J. 21 (6), 8172-8183, https://doi.org/10.1109/JSEN.2021.3050618.

Kim, M., Park, J., Shin, J., 2024. Efficient Fourier neural operators by group convolution and channel shuffling. In: ICLR 2024 Workshop on AI4DifferentialEquations in Science. available at OpenReview: https://openreview.net/forum?id=TyglZCofeE3.

Kingma, D.P., Ba, J., 2015. Adam: A method for stochastic optimization. In: Proceedings of the 3rd International Conference on Learning Representations. https://doi.org/10.48550/arXiv.1412.6980.

Kocks, F., Tomé, C., Wenk, H.-R., 2000. Texture and Anisotropy. Preferred Orientations in Polycrystals and Their Effect on Material Properties. Cambridge University Press.

Kristeková, M., Kristek, J., Moczo, P., 2009. Time-frequency misfit and goodness-of-fit criteria for quantitative comparison of time signals. Geophys. J. Int. 178 (2), 813-825, https://doi.org/10.1111/j.1365-246X.2009.04177.x.

Kube, C.M., Turner, J.A., 2015. Acoustic attenuation coefficients for polycrystalline materials containing crystallites of any symmetry class. J. Acoust. Soc. Am. 137 (6), EL476-EL482, https://doi.org/10.1121/1.4921676.

Lehmann, F., Gatti, F., Bertin, M., Clouteau, D., 2024. 3D elastic wave propagation with a Factorized Fourier Neural Operator (F-FNO). Comput. Methods Appl. Mech. Engrg. 420, 116718, https://doi.org/10.1016/j.cma.2023.116718.

Lehmann, F., Gatti, F., Clouteau, D., 2025. Multiple-input Fourier neural operator (MIFNO) for source-dependent 3D elastodynamics. J. Comput. Phys. 527, 113813, https://doi.org/10.1016/j.jcp.2025.113813.

Li, Z., Kovachki, N.B., Azizzadenesheli, K., Liu, B., Bhattacharya, K., Stuart, A.M., Anandkumar, A., 2021. Fourier neural operator for parametric partial differential equations. In: International Conference on Learning Representations. available at OpenReview: https://openreview.net/forum?id=c8P9NQVtmmnO.

Lian, Y., Du, F., Xie, L., Hu, Q., Jin, P., Wang, Y., Lu, Z., 2024. Application of laser ultrasonic testing technology in the characterization of material properties: A review. Measurement 234, 114855, https://doi.org/10.1016/j.measurement.2024.114855.

Lu, L., Jin, P., Pang, G., Zhang, Z., Karniadakis, G.E., 2021. Learning nonlinear operators via DeepONet based on the universal approximation theorem of operators. Nat. Mach. Intell. 3, 218-229, https://doi.org/10.1038/s42256-021-00302-5.

Lv, G., Guo, S., Chen, D., Feng, H., Zhang, K., Liu, Y., Feng, W., 2023. Laser ultrasonics and machine learning for automatic defect detection in metallic components. NDT E Int. 133, 102752, https://doi.org/10.1016/j.ndteint.2022.102752.

Malmström, M., Jansson, A., Hutchinson, B., Lönnqvist, J., Gillgren, L., Bäcke, L., Sollander, H., Bärwald, M., Hochhard, S., Lundin, P., 2022. Laser-ultrasound-based grain size gauge for the hot strip mill. Appl. Sci. 12 (19), https://doi.org/10.3390/app121910048.

Miorelli, R., Fisher, C., Kulakovskyi, A., Chapuis, B., Mesnil, O., d'Almeida, O., 2021. Defect sizing in guided wave imaging structural health monitoring using convolutional neural networks. NDT E Int. 122, 102480, https://doi.org/10.1016/j.ndteint.2021.102480.

Nayfeh, A.H., 1995. Wave propagation in layered anisotropic media. In: North-Holland Series in Applied Mathematics and Mechanics, Elsevier.

Nerlikar, V., Miorelli, R., Recoquillay, A., d'Almeida, O., 2024. A physics-embedded deep-learning framework for efficient multi-fidelity modeling applied to guided wave based structural health monitoring. Ultrasonics 141, 107325, https://doi.org/10.1016/j.ultras.2024.107325.

Pyle, R.J., Bevan, R.L.T., Hughes, R.R., Ali, A.A.S., Wilcox, P.D., 2022. Domain adapted deep-learning for improved ultrasonic shock characterization using limited experimental data. IEEE Trans. Ultrason. Cerroelectr. Freq. Control 69 (4), 1485-1496, https://doi.org/10.1109/TUFFC.2022.3151397.

Rahman, M.A., Ross, Z.E., Azizzadenesheli, K., 2023. U-NO: U-shaped neural operators. Trans. Mach. Learn. Res. available at OpenReview: https://openreview.net/forum?id=j3oQF9coJd.

Reuss, A., 1929. Berechnung der fließgrenze von mischkristallen auf grund der plastizitätsbedingung für einkristalle. ZAMM - J. Appl. Math. Mech. / Z. Angew. Math. Mech. 9 (1), 49-58, https://doi.org/10.1002/zamm.19290090104.

Roch, C., Tournier, C., Lavernhe, S., 2023. Process based modelling of power density for wire laser additive manufacturing using a coaxial head. Addit. Manuf. 73, 103648, https://doi.org/10.1016/j.addma.2023.103648.

Royer, D., Valier-Brasier, T., 2022. Elastic Waves in Solids 1. John Wiley & Sons, Ltd, https://doi.org/10.1002/9781119902942.

Ryzy, M., Grabec, T., Sedlák, P., Veres, I.A., 2018. Influence of grain morphology on ultrasonic wave attenuation in polycrystalline media with statistically equixed grains. J. Acoust. Soc. Am. 143 (1), 219-229, https://doi.org/10.1121/1.5020785.

Safdar, M., Paul, P.P., Lamouche, G., Wood, G., Zimmermann, M., Hannesen, F., Bescond, C., Wanjara, P., Zhao, Y.F., 2024. Fundamental requirements of a machine learning operations platform for industrial metal additive manufacturing. Comput. Ind. 154, 104037, https://doi.org/10.1016/j.compind.2023.104037.

Saltelli, A., Ratto, M., Andres, T., Campolongo, F., Cariboni, J., Gatelli, D., Saisana, M. Tarantola, S., 2008. Global Sensitivity Analysis: The Primer. Wiley, Chichester (England).

Scruby, C.B., Drain, L.E., 1990. Laser Ultrasonics Techniques and Applications, first ed. Taylor & Francis, https://doi.org/10.1201/9780203749098.

Segovia Ramírez, I., García Márquez, F.P., Papaelias, M., 2023. Review on additive manufacturing and non-destructive testing. J. Manuf. Syst. 66, 260-286, https://doi.org/10.1016/j.jmsy.2022.12.005.

Sha, G., Huang, M., Lowe, M.J.S., Rokhlin, S.I., 2020. Attenuation and velocity of elastic waves in polycrystals with generally anisotropic grains: Analytic and numerical modeling. J. Acoust. Soc. Am. 147 (4), 2442-2465, https://doi.org/10.1121/10.0001087.

Shaloo, M., Schnall, M., Klein, T., Huber, N., Reitinger, B., 2022. A review of non-destructive testing (NDT) techniques for defect detection: Application to fusion welding and future wire arc additive manufacturing processes. Materials 15 (10), https://doi.org/10.3390/ma15103697.

Shukla, K., Oommen, V., Peyvan, A., Penwarden, M., Plewacki, N., Bravo, L., Ghoshal, A., Kirby, R.M., Karniadakis, G.E., 2024. Deep neural operators as accurate surrogates for shape optimization. Eng. Appl. Artif. Intell. 129, 107615, https://doi.org/10.1016/j.engappai.2023.107615.

Singh, J., Tant, K., Curtis, A., Mulholland, A., 2022. Real-time super-resolution mapping of locally anisotropic grain orientations for ultrasonic non-destructive evaluation of crystalline material. Neural Comput. Appl. 34, 4993-5010, https://doi.org/10.1007/s00521-021-06670-8.

Sobol', I.M., 2001. Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates. Math. Comput. Simulation 55 (1), 271-280, https://doi.org/10.1016/S0378-4754(00)00270-6, The Second IMACS Seminar on Monte Carlo Methods.

Sobol', I., Kucherenko, S., 2009. Derivative based global sensitivity measures and their link with global sensitivity indices. Math. Comput. Simulation 79 (10), 3009-3017, https://doi.org/10.1016/j.matcom.2009.01.023.
Stanke, F.E., Kino, G.S., 1984. A unified theory for elastic wave propagation in polycrystalline materials. J. Acoust. Soc. Am. 75 (3), 665-681, https://doi.org/10.1121/1.390577.

Stoklasová, P., Grabec, K., Zoubková, K., Sedlák, P., Krátký, S., Seiner, H., 2021. Laser-ultrasonic characterization of strongly anisotropic materials by transient grating spectroscopy. Exp. Mech. 61, 663-676, https://doi.org/10.1007/s11340-021-00698-6.

Szegedy, C., Vanhoucke, V., Ioffe, S., Shlens, J., Wojna, Z., 2016. Rethinking the Inception architecture for computer vision. In: 2016 IEEE Conference on Computer Vision and Pattern Recognition. CVPR, pp. 2818-2826, https://doi.org/10.1109/CVPR.2016.308.

Torrence, C., Compo, G.P., 1998. A practical guide to wavelet analysis. Bull. Am. Meteorol. Soc. 79 (1), 61-78, https://doi.org/10.1175/1520-0477(1998)079<0061:APGTWA>2.0.CO;2.

Tran, A., Mathews, A., Xie, L., Ong, C.S., 2023. Factorized Fourier neural operators. In: The Eleventh International Conference on Learning Representations. available at OpenReview: https://openreview.net/forum?id=tmiIMPI4IPA.

Tu, X.L., Zhang, J., Gambaruto, A.M., Wilcox, P.D., 2024. Finite element modelling strategy for determining directivity of thermoelastically generated laser ultrasound. Ultrasonics 138, 107252, https://doi.org/10.1016/j.ultras.2024.107252.

Uhlig, S., Alkhasli, I., Schubert, F., Tschöpe, C., Wolff, M., 2023. A review of synthetic and augmented training data for machine learning in ultrasonic non-destructive evaluation. Ultrasonics 134, 107041, https://doi.org/10.1016/j.ultras.2023.107041.

Van Pamel, A., Brett, C.R., Huthwaite, P., Lowe, M.J.S., 2015. Finite element modelling of elastic wave scattering within a polycrystalline material in two and three dimensions. J. Acoust. Soc. Am. 138 (4), 2326-2336, https://doi.org/10.1121/1.4984290.

Victoria Giraldo, J.C., Tie, B., 2024. Theoretical and numerical simulation of bulk and surface ultrasonic waves scattering in microstructures with elongated grains. In: CSMA2024, 16e Colloque national en calcul des structures. CSMA, Giens (Var), France, available at https://hal.science/hal-04490154.

Voigt, W., 1928. Lehrbuch Der Kristallphysik. B.G. Teubner, Leipzig.

Weaver, R., 1990. Diffusivity of ultrasound in polycrystals. J. Mech. Phys. Solids 38 (1), 55-86, https://doi.org/10.1016/0022-5096(90)90021-U.

Wu, S.-H., Tariq, U., Joy, R., Sparks, T., Flood, A., Liou, F., 2024. Experimental, computational, and machine learning methods for prediction of residual stresses in laser additive manufacturing: A critical review. Materials 17 (7), https://doi.org/10.3390/ma17071498.

Yu, H., Yin, A., Xu, Z., Zhang, J., Wu, J., Xu, X., Zhang, Z., 2023. Grain size characterization of TA1 with GA-BP neural network using laser ultrasonics. Optik 275, 170600, https://doi.org/10.1016/j.ijleo.2023.170600.

Zhang, Z., Li, B., Huang, Y., Qin, F., Chen, L., 2024. Surface roughness characterization of additive manufactured Ti-6Al-4 V based on laser ultrasonic signal. Opt. Laser Technol. 177, 111070, https://doi.org/10.1016/j.optlastec.2024.111070.

Zhao, H., Fu, C., Zhang, Y., Zhu, W., Lu, K., Francis, E.M., 2024. Dimensional decomposition-aided metamodels for uncertainty quantification and optimization in engineering: A review. Comput. Methods Appl. Mech. Engrg. 428, 117098, https://doi.org/10.1016/j.cma.2024.117098.

Zhu, M., Feng, S., Lin, Y., Lu, L., 2023. Fourier-DeepONet: Fourier-enhanced deep operator networks for full waveform inversion with improved accuracy, generalizability, and robustness. Comput. Methods Appl. Mech. Engrg. 416, 116300, https://doi.org/10.1016/j.cma.2023.116300.

Zou, C., Azizzadenesheli, K., Ross, Z.E., Clayton, R.W., 2024. Deep neural Helmholtz operators for 3-D elastic wave propagation and inversion. Geophys. J. Int. 239 (3), 1469-1484, https://doi.org/10.1093/gji/ggae342.