NDT&E International 116 (2020) 102327

Contents lists available at ScienceDirect

NDT and E International

journal homepage: http://www.elsevier.com/locate/ndteint

Optimization algorithms for ultrasonic array imaging in homogeneous

anisotropic steel components with unknown properties
Corentin Ménard a, *, Sébastien Robert a, Roberto Miorelli a, Dominique Lesselier b
a
CEA-LIST, Department of Imaging and Simulation for Non-Destructive Testing, Bât. 565, PC 120, 91191, Gif-sur-Yvette Cedex, France
b
Université Paris-Saclay, CNRS, CentraleSupélec, Laboratoire des signaux et systèmes, 91190, Gif-sur-Yvette, France

A R T I C L E I N F O A B S T R A C T

Keywords: In this paper, we investigate the ability of an adaptive ultrasonic method to image point-like reflectors inside
Ultrasonic array anisotropic and homogeneous nuclear steels with unknown properties, such as V-shape welds or cladded com­
Total Focusing Method ponents. The optimization scheme combines the Total Focusing Method (TFM) imaging algorithm with a
Adaptive imaging
customized gradient ascent method to improve both the quality and the reliability of ultrasound images. A
Anisotropic materials
Optimization algorithm
statistical analysis of its robustness is performed with simulated echoes and using a surrogate model to speed up
Surrogate-model the computation times of the TFM images. Then, the optimization procedure is evaluated with several experi­
mental cases and provides highly enhanced images with a 5 MHz array. The positioning of the artificial defects of
2.0 mm diameter is estimated with less than 1 mm error with respect to their actual position, and the signal-to-
noise ratio is increased by up to 10 dB. The elastic properties are also estimated with less than 10% error when
compared to their actual values.

1. Introduction properties are not known during inspection. If the ray-based recon­
struction model significantly deviates from the wave propagation in the
Ultrasonic arrays are becoming a standard technology in many in­ physical medium, the sum of delayed inter-element signals becomes
dustries since they allow the application of high-performance detection incoherent at every focusing point, leading to strong image degradation
and imaging methods. The Total Focusing Method (TFM) [1] is a recent and unreliable diagnosis. We therefore propose an optimization pro­
example of an imaging method that is becoming a reference in nonde­ cedure based on a gradient ascent algorithm to iteratively enhance the
structive testing thanks to its ability to provide high-resolution images. reliability of the TFM images in unknown anisotropic structures. In the
Its basic principle consists in synthetically focusing the inter-element following, the algorithm is referred to as Gradient Ascent Optimization
signals recorded with an array in pulse-echo mode on every point of (GAO). The semi-analytical ray-tracing algorithm [8] of the CIVA soft­
the region of interest using a “delay-and-sum” strategy. As the set of ware is called at each iteration for the TFM reconstruction. The first
recorded signals physically corresponds to the matrix of the impulse image is computed assuming an isotropic propagation (as in most on-site
responses of the medium [2,3], one of the advantages of TFM is that inspections when the material properties are not well known) and, if a
focusing is achieved a posteriori, without applying delay laws at the time low-amplitude echo is visible, the procedure is started. During the
of acquisition. Thus, adaptive imaging algorithms can be applied on a optimization, the properties of the inspected structure are modified until
single set of data to optimize images in media with uncertain or un­ the defect echo amplitude is maximized. Once convergence reached,
known properties. For example, in Ref. [4,5] adaptive TFM under both imaging and material characterization are achieved, as a reliable
complex and unknown surfaces are proposed, and in Ref. [6], the image of the defect and estimates of the material properties in the vi­
principle is extended to multi-mode TFM imaging [7] when the profiles cinity of the flaw are provided. If the material is assumed homogeneous,
of the outer and inner surfaces are unknown. then the material characterization is achieved for the whole structure.
In the present work, TFM is used to image point-like reflectors in This optimization procedure is inspired by the works of Pudovikov et al.
anisotropic butt welds and cladded steel components, whose elastic [9] and Juengert et al. [10] where the image of a well known reflector (e.

* Corresponding author.
E-mail addresses: corentin.menard@cea.fr (C. Ménard), sebastien.robert@cea.fr (S. Robert), roberto.miorelli@cea.fr (R. Miorelli), dominique.lesselier@l2s.
centralesupelec.fr (D. Lesselier).

https://doi.org/10.1016/j.ndteint.2020.102327
Received 2 March 2020; Received in revised form 29 May 2020; Accepted 9 July 2020
Available online 16 August 2020
0963-8695/© 2020 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license
(http://creativecommons.org/licenses/by-nc-nd/4.0/).
C. Ménard et al. NDT and E International 116 (2020) 102327

g., a flat-bottom hole) is iteratively optimized. The difference with the

following study lies in that our procedure maximizes the detection
amplitude of a reflector without prior information about its properties
(location, geometry, size or orientation). The present approach is an
interesting alternative to destructive characterization techniques, which
require a sample to be taken from the structure in order to characterize
its properties using specific various experimental setups, such as ultra­
sonic transmission measurements with two monolithic transducers [11]
or two transducer arrays [12], pitch-catch ultrasonic measurements
with two arrays [13,14], EBSD (Electron Backscatter Diffraction) scans
[15], SAW (Surface Acoustic Wave) scans [16], etc.
The efficiency of the optimization is experimentally exemplified by Fig. 1. 2D geometry describing the aforementioned notations for TFM imaging:
the inspection of a V-shape weld and a cladded steel component. Their propagation times tn (r) and tm (r) for a couple of emitter/receiver (n, m) and a
properties are assumed to be homogeneous and anisotropic, with an pixel P located by r in the ROI, inside an isotropic medium.
orthotropic crystal structure. For the V-shape weld, 5 unknown pa­
rameters govern the 2D wave propagation: 4 stiffness constants and the following. Thus, the image formation needs the computation of N ×
orientation of the crystal structure. Regarding the cladded component, 5 Nx × Nz propagation times for a resolution of Nx × Nz pixels of the region
unknowns are again considered: the same 4 stiffness constants and the of interest (ROI). Nx and Nz stand for the number of pixels along the X
thickness of the anisotropic layer. In order to evaluate the robustness of and Z axes, respectively.
the procedure, a statistical analysis is conducted with simulated echoes. In order to image inside anisotropic structures with TFM, the
For this purpose, a surrogate model using a Kernel Ridge Regression computation of the N propagation times tn (r) is carried out in the CIVA
(KRR) bypasses the forward model in CIVA and shortens the computa­ software by using a semi-analytical forward model, the so-called “pencil
tion times needed for the statistical study. KRR is a fast non-linear method”. It was implemented by Gengembre et al. [8] to simulate wave
regression algorithm that performs well when interpolating a database propagation in isotropic or anisotropic media. Further modeling studies
of reasonable size and dimensions [17], which is why it is suited for the extended this model to continuously variable anisotropic media [21].
current study. Thanks to the surrogate model, a sensitivity analysis is This method consists in tracing collections of rays with a small diver­
also carried out to indicate the relevant parameters to be varied during gence, instead of single rays, in order to account for the intensity
the procedure. Finally, the local GAO is tested against a global opti­ decrease of the wave along its path. The extension to continuously
mizer: the Particle Swarm Optimization (PSO) [18,19]. PSO is a simple variable materials involves dividing the propagation path into elemen­
algorithm to implement, which swiftly produces results for a given tary sub-interfaces across which the wave travels, and applies the pencil
optimization problem and does not require many trials to tune its method at each of these elementary interfaces. Propagation times are
hyper-parameters. computed in a first step for solving the eikonal equation (for the axial
The paper is organized as follows: Section 2 describes the optimi­ quantities) and the amplitude is retrieved from solving the transport
zation procedure for the TFM imaging; Section 3 is devoted to the study equation in a second step (for the paraxial quantities).
of its performances and robustness with numerical data and using a The optimization procedure therefore holds in iteratively modifying
surrogate model; Section 4 discusses the optimization results when the map of propagation times in CIVA, in order to optimize the quality of
dealing with experimental data. the TFM image, step after step.

2. Optimization method
2.2. Optimization procedure
In this section, a few reminders on TFM imaging for anisotropic and
The synoptic diagram of the optimization procedure is shown in
inhomogeneous media in the CIVA software are provided, especially
Fig. 2 and consists of 4 steps. After acquiring the inter-element response
concerning the computation of the inter-element propagation times to
matrix, the theoretical propagation times for TFM are computed with
synthetically focus on every point. The different steps of the optimiza­
the pencil method in CIVA (Step 1), using randomly initialized material
tion procedure are then presented and the GAO, customized for this
properties. Then, based on these propagation times, a TFM image is
specific study, is next detailed.
formed (Step 2) following Eq. (1). From this image, an estimator of
quality is extracted (Step 3) in order to take the optimization decision. If
2.1. TFM imaging in anisotropic structures the image is still not satisfactory, the optimization algorithm indicates
how to modify the unknown material properties (Step 4) and returns
Noting N the number of array elements, the principle of TFM consists back to Step 1. The optimization loop stops when the image is consid­
in coherently summing the N × N signals of the inter-element impulse ered optimized enough. The estimator of the image quality at Step 3 is
response matrix S(t) [20]. A component snm (t) of S(t) corresponds to the the maximum amplitude of the echo in the TFM image (Amax). This
signal received by element m when element n is excited by an electric indicator tends to increase as the theoretical propagation times converge
pulse. As we are interested in complex and unknown materials, the TFM towards those of the waves in the physical medium, provided that the
method is considered in its simplest form, without weighting factors. echo corresponds to direct ray paths and that the wave polarization
The amplitude at a pixel P located by position vector r (Fig. 1) is written (longitudinal or transverse) is known. If the estimation of the propaga­
as: tion times in the medium is not correct, the signals would be summed at
⃒
⃒ ∑N ∑ N
⃒
⃒ wrong pixels, leading to a decrease in the amplitude of the echo.
⃒ ⃒
A(r) = ⃒ s̃nm [tn (r) + tm (r) ] ⃒. (1)
⃒ m=1 n=1 ⃒
2.3. Gradient Ascent Optimization
s̃nm (t) = snm (t) + jH [snm (t)] is the analytic signal associated with the
recorded signal snm (t), where H [ ⋅] denotes the Hilbert transform. tn (r) Step 4 of the optimization procedure uses a local optimizer to pro­
and tm (r) are the propagation times from n to P and from P to m, duce optimized images. It can also give fair estimates of the properties of
respectively, and satisfy tn (r) = tm (r) when n = m, since images will be the inspected material, in the case of a homogeneous anisotropy. This
computed with direct rays paths without mode conversion in the algorithm is a modified Gradient Ascent Optimization (GAO) [22],

2
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 2. Synoptic diagram of the TFM optimization procedure involving 4 main post-processing steps. Step 1 is the computation of the propagation times from a set of
weld parameters. Steps 2 and 3 consist in the computation of the TFM image and in the measurement of the image quality descriptor, respectively. Step 4 determines
how to modify the material properties to increase the image quality, and is achieved thanks to an optimization algorithm.

which locally follows the gradient of the Amax descriptor to retrieve the ( ) ( )
optimum solution corresponding to the global maximum inside the so­ 1− j
J
cos2 π j
2 J
+ ηJ
lution space. The GAO is customized with the two enhancements ηj = η0 , (4)
mentioned in Ref. [22]: a decaying learning rate and a periodic warm 1 + ηJ
restart. The gradient-based optimization may be written as:
where η0 and ηJ are the values of the learning rate at the beginning and
( )
xj+1 = xj + ηj ∇f xj (2) the end of the decaying phase, respectively. The warm restart and the re-
initialization of the decaying phase happen after J iterations. The new
ηj is the learning rate at iteration step j ∈ N of the gradient ascent, restart position x0 is the last position xJ reached during the previous
which follows a decaying law as j increases. A position, designated by decaying phase.
the material property vector xj (in our case, xj ∈ R5 ), is associated to an
Amax value f(xj ). The warm restart forces the decaying learning rate to 3. Numerical study of the robustness of the method using a
reset after a fixed number of iterations J, or after converging, which surrogate model
allows the algorithm to escape from local maxima and saddle points.
Since the convergence of the algorithm is not appraised in terms of a loss In this section, the optimization procedure is evaluated with a noise-
function, a few warm restarts are needed to ensure that GAO reaches the free example. According to Eq. (3), the full gradient estimation at a given
vicinity of the optimum solution at least once. The gradient of the Amax position requires 2 Amax values for each dimension of the problem.
descriptor is locally estimated using centered finite differences at each Thus, to estimate the gradient in a 5-dimensional space, 10 TFM images
iteration step: have to be computed. Since the aim is to perform a statistical analysis,
( ) ( ) the procedure has to be extensively tested. For this reason, a surrogate
( ) f xj + h − f xj − h
∇f xj = , (3) model is implemented to bypass the time-consuming computation of the
2h theoretical propagation times and the TFM images (Steps 1 and 2 in
where h is an infinitesimal perturbation vector, and xj the position at Fig. 2).
which the gradient is computed. We first describe the sample considered in the numerical study,
The literature provides a great variety of possible decaying learning which is representative of the welding sample in the experimental study.
rates [22–25]. In the examples covered in this paper, the learning rate is The KRR surrogate model is detailed next and used for a sensitivity
chosen so as to benefit from the speed of a linear decrease and from the analysis. Then, some insights on the Particle Swarm Optimization (PSO)
accuracy of a cosine decrease: are provided. Finally, the robustness of the procedure is statistically
analyzed and discussed.

3
C. Ménard et al. NDT and E International 116 (2020) 102327

3.1. Case study for simulations assume an isotropic propagation), and Nz = 101 pixels (sampling step of
λ/10) in the axial direction.
The numerical case study is an anisotropic and homogeneous Inconel
182 weld with an orthotropic crystalline structure. The inspection
configuration is presented in Fig. 3: a linear array is immersed in water 3.2. Surrogate model by Kernel Ridge Regression
and horizontally placed over the sample featuring 3 side-drilled holes of
2 mm diameter located at different depths. The water column height is The Ridge Regression (RR) is a linear regression algorithm that in­
set to 20 mm. The array operates around 5 MHz and is composed of 64 troduces a regularization term in the interpolation problem to decrease
elements of 0.5 mm width with a 0.6 mm pitch. Elements are excited the variance of the regression, with the trade-off of an increased bias. In
with a Gaussian-type pulse with 80% bandwidth at − 6 dB and sampled the current study, its implementation is achieved using Scikit-learn, a
at 100 MHz. Echoes are simulated in CIVA using the method of sepa­ Machine Learning Python library [28]. Let us consider a regression
ration of variables [26], without accounting for attenuation. A 2D im­ problem on a data set formed by D variables, which are sampled by P
aging problem with longitudinal waves is assumed, therefore reducing points by dimension on a regular grid. The total number of samples in
the optimization problem to 5 material parameters, to describe the the data set is then M = PD . In this database, let us consider the
phase velocity VQL of the Quasi-Longitudinal (QL) waves [27]. In our observable parameter y = (y1 ,…,yi ,…yM ), such that y ∈ RM×1 (e.g., the
case, the density ρ = 8.26 g.m− 3 of the medium is assumed to be known: maximum amplitude of the TFM image). We want to predict the value of
√̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅
√̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅ y due to variations of some control parameters X = (x1 ,…xi ,…xM ), with
Γ11 + Γ33 + (Γ11 − Γ33 )2 + 4Γ213 X ∈ RM×D and xi ∈ R1×D (e.g., the 4 elastic constants and the inclination
VQL = √̅̅̅̅̅ . (5) angle of the crystalline structure). A standard linear regression would
2ρ
use the ordinary least square method to compute the weighting vector
The Christoffel components are given by w ∈ RD×1 , so that y = Xw. However, if the number of samples M is too
⎧
⎨ Γ11 = C11 sin2 (θi − θ) + C55 cos2 (θi − θ) small, the variance of w may be too large, thus making the prediction for
Γ = C55 sin2 (θi − θ) + C33 cos2 (θi − θ) , (6) unseen data unreliable. For this reason, the RR controls the variance of w
⎩ 33 by introducing the penalization term μ ∈ R. The RR then consists in
Γ13 = (C13 + C55 )cos(θi − θ)sin(θi − θ)
minimizing the cost function J(w):
where C11 , C13 , C33 and C55 are the 4 unknown elastic constants and θis
J(w) = ||Xw − y||2 + μ||w||2 . (7)
the unknown tilt angle of the crystalline structure around the Y axis. The
θi parameter is the incidence angle of the wave in the medium. The The use of the Euclidean norm forces weighting coefficients in w to
actual material properties are given in Fig. 3, along with the related decrease when they relate to less contributing variables of the problem.
slowness curves of the QL waves of interest. It however does not set them to zero, unlike the Lasso Regression [29].
In order to quicken the procedure, the ROI is limited to an area of By solving ∇J(w) = 0, it can be shown [17] that the w coefficients are:
10 × 10 mm2 centered around reflector no2. The image resolution is ( )− 1
Nx × Nz = 5151 pixels with Nx = 51 pixels in the lateral direction w = XT X + μI XT y. (8)
(sampling step of λ/5, where λ is the wavelength at 5 MHz when we
The estimation of a new sample y’ located at x’ is then:

Fig. 3. Scheme of the numerical setup (a):

immersion inspection using an array trans­
ducer placed 20 mm above the surface;
Inconel-based V-shape weld mockup with 3
side-drilled holes at depths 20 mm (1), 30
mm (2) and 40 mm (3) from the surface; the
holes are spaced 20 mm apart along the
azimuthal axis and are of 2 mm diameter.
Slowness curves (b): quasi-longitudinal wave
inside the weld with an orthotropic crystal
structure (green) and longitudinal wave in­
side an isotropic medium with velocity 5500
m s− 1 (orange). The 0∘ angle corresponds to
the Z axis in the setup schema and the 90∘
angle to the X axis. In Table (c): actual values
of the unknown properties. (For interpreta­
tion of the references to color in this figure
legend, the reader is referred to the Web
version of this article.)

4
C. Ménard et al. NDT and E International 116 (2020) 102327

( )− 1
y’ = wT x’ = yT XXT + μI Xx’. (9) ∑
SStot = (yi − y)2 . (13)
i
The Kernel Ridge Regression associates the Ridge Regression with the
“kernel trick” to transform the linear regression into a non-linear one. The aim is then to find the kernel design that provides the best co­
Let us consider any variable change applied to the i-th parameter xi → efficients of determination for both the regression on the training data
Φ(xi ). Let us also define an inner product Kkl = Φ(xk )Φ(xl ), where k and l and the prediction on testing data. The closer the R2 score to 1 (by lower
refer to training sample indices. One can show [30] that Eq. (9) becomes values), the better the regression/prediction. The input data X are
standardized to have a null mean value and a unit variance for each
y’ = yT (K + μI)− 1 κ, (10) dimension of the problem. This way, none of the dimensions has a
variance larger than the others that could dominate the R2 score. The
where a component κk = xk x’. Consequently, the variable change is of output y, obtained by simulation, exhibit values that are of unit order.
no concern to perform the regression, only the inner product is required. The best kernel design, compared with linear and polynomial, is found
This property is called the “kernel trick”. The choice of the kernel de­ to be a Gaussian one:
pends on the addressed problem. ( )
The choice of the database is a matter of trade-off between the K(xk , xl ) = exp − σ ||xk − xl ||2 (14)
desired sampling step of the database and the computational cost of the
KRR. In the present case, a regular grid is chosen for its simplicity. A with σ = 3.162⋅10− 4 and with a regularization parameter μ = 1 for the
preliminary study of sections of the 5D space revealed that the variations ridge penalization. The corresponding training and testing scores are
of the Amax estimator are qualitatively rather smooth along each R2train = R2test ≃ 1.0. The respective exploration range of the σ and μ
dimension. Moreover, in most cases a loose grid is sufficient as the KRR hyper-parameters were set to [10− 6 , 101 ] and [10− 7 , 102 ] with a 0.25
benefits from the prior knowledge provided by the kernel choice. It is logarithmic step size.
found that a database sampling of 10 GPa for the elastic constants and 0. In order to validate the best KRR model, a comparison between
7◦ for the tilt angle is sufficient for describing the Amax variations with actual unseen samples (i.e., new values of Amax obtained from simu­
satisfactory accuracy. Complexity-wise, the Scikit-learn implementation lations with CIVA) and their KRR estimates was conducted. The com­
uses a Singular Value Decomposition to compute the pseudo-inverse in parison was made using a set of 5000 new samples generated with a
Eq. (10), with a computational complexity about O(M2 ), where M is the Latin hypercube in the weld parameters space. We recall that Amax is
number of samples [31]. A database sampled with M = PD = 16 807 the maximum amplitude of the TFM image in a 10 × 10 mm2 region that
samples is chosen, i.e., P = 7. The limits of the database space are given only contains the reflector no2. The actual and estimated Amax values
in Table 1. They are set in order to provide a sufficiently large explo­ are normalized and plotted in Fig. 4. The regression is accurately ach­
ration area while including the actual material properties. ieved with a mean absolute error equal to 2.4 ± 2.8% for predicting on
The choice of the hyper-parameters of the KRR (ridge regularization unseen samples and can now be substituted to CIVA for generating new
term, kernel type and kernel parameters) is a trade-off between the ac­ [x, y] samples.
curacy of the regression and the reliability of its prediction to unseen
data (bias vs variance). The aforementioned database of 16 807 samples
is used in a threefold cross-validation procedure, where 80% of the 3.3. Sobol sensitivity analysis
samples are randomly chosen to train the KRR, and the remaining 20%
are left for testing the prediction. This method helps avoiding overfitting Thanks to the KRR, we now have a model M such that y = M (x),
and selection bias. The best kernel is then chosen by trials and measuring with x ∈ R1×D the position vector of y in the D-dimensional space. A
the corresponding errors: for each threefold cross-validation, the design sensitivity analysis can then be conducted using Sobol indices [32]. The
and hyper-parameters of the kernel are modified and a score is given to aim of this analysis is to quantify the contribution of each input
quantify the model error. To quantify the accuracy of the regression and parameter x(p) ∈ R (1 ≤ p ≤ D), or combinations of input parameters, on
the reliability of the prediction, the coefficient of determination (noted
R2) metric is used during the training and testing phases. Let us define
the mean of the observed data used for the training/testing phase: y =
1
∑m0
m0 i=1 yi , where m0 < M is the number of samples. The coefficient of
determination for the prediction is given by
SSres
R2 = 1 − , (11)
SStot

where the residual sum of squares of the prediction is

∑
SSres = (yi − y’i )2 , (12)
i

while the total sum of squares is

Table 1
Upper and lower boundaries of the weld parameter space to compute the KRR-
based surrogate model.
Parameters Lower bounds Upper bounds

C11 (GPa) 230 290 Fig. 4. Comparison between the actual values of the image quality descriptor
C13 (GPa) 90 150 Amax (blue solid line) with the estimated values given by the Kernel Ridge
C33 (GPa) 230 290 Regression surrogate model (red dots), using 5000 new samples computed by
C55 (GPa) 70 130 simulations and randomly generated. (For interpretation of the references to
θ (∘) 6 10 color in this figure legend, the reader is referred to the Web version of
this article.)

5
C. Ménard et al. NDT and E International 116 (2020) 102327

the variance of the output y ∈ R. In our case, this means identifying the images are uniformly and randomly generated by varying the input
most relevant material properties for governing the variations of the parameters inside the database delimited by the parameter values in
Amax descriptor. The analysis is performed using the Python’s library: Table 1. Out of these images, the coordinates of the Amax pixel are
SALib [33]. retrieved and displayed on a synthetic view in Fig. 6. This figure displays
Three kinds of Sobol indices are derived. The higher the index value, a close-up around the side-drilled hole no2, the latter being represented
the higher the contribution of the input(s). The first-order indices do by the black circle. Two cases are considered: (a) the positions of the
indicate the contribution of an input parameter x(p) when it is the only Amax pixel when only the three most contributing parameters (C13 , C33
one to vary: and C55 ) are varying, and (b) the positions when only the two less
[ ( ⃒ )] contributing parameters (C11 and θ) are varying. The red marker locates
Var E y⃒x(p) the Amax pixel position in the image computed with the actual material
Sp = . (15)
Var(y) properties. Its Amax value is set to 0 dB and serves as reference ampli­
The second-order indices indicate the effect of two parameters x(p) tude. The green markers represent the Amax pixels when
− 3 dB ≤ Amax ≤ 0 dB. The blue markers are the Amax pixels for which
and x(q) varying at a time, with p ∕ = q:
− 6 dB ≤ Amax < − 3 dB. The black markers are the pixels when −
[ ( ⃒ (p) (q) ) ( ⃒ (p) ) ( ⃒ )]
Var E y⃒x , x − E y⃒x − E y⃒x(q) 12 dB ≤ Amax < − 6 dB. In Fig. 6a, the sensitivity of the position of
Sp,q = . (16)
Var(y) the echo to the three most contributing parameters is easily observable:
the markers create a scatter plot around the actual position of the echo
Finally, total-effect indices provide information on the overall and indicate position errors up to 5 mm. Moreover, it can be observed
contribution of a single input parameter. They take into account the that the set of points is globally arranged in three separated sub-
first-order indices and all higher orders with which the parameter of domains, which demonstrates the existence of a relation between
interest is involved: Amax and the correct position of the echo: the highest Amax (green
∑ ∑
STp = Sp + Sp,q + Sp,q,r + … + S1,2,…,p,…,D. (17) markers) are located in the sub-domain including the actual echo posi­
q∕
=p q<r tion. In Fig. 6b, the variations of the echo position with respect to var­
iations of the two less contributing parameters are also visible. The
q,r∕
=p

variation pattern is explained by C11 only contributing to the wave speed

In each equation, p, q and r refer to one component of the x vector. variations along the X axis of the crystalline structure, which is tilted by
The results of the sensitivity analysis are provided in Fig. 5: for θ. The variation range of the latter is small enough around the actual
readability, only indices with values greater than 5⋅10− 3 are shown, value for the markers to form this pattern. The position error of the echo
which is why C11 is not present. It is observed that the most contributing is up to 4 mm and the colored labels are again delimited in space.
parameters to the optimization problem are C55 and C13 . The main part Therefore, even if they exhibit tiny contributions to the variations of the
of their contributions comes from their first-order index. The remaining Amax value, C11 and θ do have an effect on the location of the defect
arises from their interactions with other parameters. The elastic echo. For this reason, these two parameters cannot be discarded from
component C33 contributes in a lesser way to the variance of Amax, the optimization procedure.
although its contribution is mainly related to its interactions with C55 It is interesting to discuss the distribution of the green markers,
and C13 . Finally, both θ and C11 exhibit little contributions to the which are to represent the most optimized echo position with respect to
problem, especially the latter. From these observations, one can deduce Amax. The pattern shows that an individual optimization may lead to a
that only C55 , C13 and C33 do have a major contribution to the variance wrong echo positioning, as the standard deviation of the distribution is
of Amax, while θ and C11 could be discarded. It is also noted that the sum quite large. However, if we consider the green markers as a swarm and
of the indices is not equal to unity: this is due to the 5 input parameters take its average position, we now obtain a fair estimate of the actual
not being fully uncorrelated [34]. defect position, while reducing the risk of reaching an outlier. Thus,
However, the correct location of the echo has also to be considered. during the optimization procedure, the GAO needs to be averaged over
In order to illustrate the effect of the 5 parameters on the echo position, several trials in order to give a correct estimate of the echo position.
It should also be noted that the results of this sensitivity analysis are
strongly case-dependent. For a given material, the values of Sobol
indices may vary depending on the parameters of the inspection
configuration, such as the aperture of the array, its frequency, or its tilt
angle for oblique incidence imaging. In the present case, the indices
would also be different if the analysis were performed with flaws located
away from the transducer central axis, and especially with flaws located
beyond the array aperture. Wherefore, the aforementioned results
cannot be taken as ground truth: they are dependent on this specific
study, for which we only intend to focus the optimization around the
flaw no2 in the mockup.

3.4. Particle swarm optimization

For reminder, the Particle Swarm Optimization (PSO) algorithm is

implemented to assess the good convergence of the GAO algorithm. The
PSO is a global optimization tool proposed by Kennedy and Eberhart
[18,19] that was applied to the learning phase of neural networks. A
swarm of particles is randomly initialized and evolves in the solution
space. The displacement of each particle is dependent on its previous
inertia, on the best recorded position during its past displacements, and
Fig. 5. First-order, second-order and total-effect Sobol indices of the optimi­ on the best position of the whole swarm. Let us consider the swarm of
zation problem. Only the indices whose value is above 5⋅10− 3 are displayed. particles represented by the position matrix Xj in the D-dimensional
The others are considered too low to contribute to the Amax variance. space, at optimization step j. The swarm obeys

6
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 6. Illustration of bad material properties estimation on the positioning of the maximum amplitude of the echo with randomly generated images: (a) variations of
the most contributing parameters (C13 , C33 and C55 ) and variations (b) of the less contributing parameters (C11 and θ). Position of the Amax pixel in the image
computed with the actual material properties (red marker); positions of Amax pixels for − 3 dB ≤ Amax ≤ 0 dB (green markers); positions for − 6 dB ≤ Amax < −
3 dB (blue markers); positions for − 12 dB ≤ Amax < − 6 dB (black markers). (For interpretation of the references to color in this figure legend, the reader is referred
to the Web version of this article.)

Xj+1 = Xj + Vj+1 , (18)

Table 2
Hyper-parameter values for both GAO and PSO algorithms.
where Vj+1 is the inertia of the particles at step j + 1:
( ) ( ) GAO PSO
Vj+1 = αVj + β Qp − Xj + γ Qs − Xj (19)
J 15 α 0.5
Qp contains the best particle positions: best positions found by each η0 (Cij ) 10 GPa β 0.5
particle on its path. Qs is the best swarm position: best position ever η0 (θ) 0.5 ∘
γ 0.5
found by the whole swarm. The α, β and γ coefficients are chosen in
order to penalize or not the displacement of a particle with respect to Qp
and Qs . V0 is usually set to be null. This algorithm is based on the properties. These values are averaged over the 5 dimensions of the
collaboration of individuals to reach a common goal and mimics the problem. The average error of the GAO is 7.5 ± 6.4% while it is 7.3 ±
auto-organization behavior of a flock of birds. 7.8% for the PSO. When compared to the actual material properties,
Both implementation and customization of this algorithm are easily both GAO and PSO give fair estimates. These non-negligible errors were
achieved and the hyper-parameters tuning phase does not require much expected as the optimization problem is ill-posed. We modify 5 material
effort. PSO is able to handle optimization problems with many di­ properties to account for the wave velocity in the medium during the
mensions and still provides a fast convergence. In this study, the PSO is optimization. To illustrate that the non-uniqueness of the solution is not
used in its simplest form, thanks to the optimization problem having an issue, the slowness curves of the QL waves are presented along with
only 5 dimensions. Several extensions can be found in the literature to their corresponding errors in Fig. 7. The error with respect to the actual
enhance the PSO: cluster of swarms [35], particle confinement [36], slowness is inferior to 1.5% with either GAO and PSO. Consequently, the
dynamic inertia coefficient [37], hybridization [38], etc. optimization is validated when estimating the unknown anisotropic
In the following, the convergence of the swarm is controlled by the homogeneous material properties of the inspected structures.
standard deviation along each dimension of the position matrix Xj and The 100 trials results are averaged to compute the optimized TFM
by the standard deviation of the associated yj , while it searches for the images. As the anisotropy is homogeneous, the elastic parameters
maximum Amax value. The optimization stops when one of the super­ coming from the optimization procedure around defect no2 can be used
vised standard deviation reaches a selected threshold. to compute a large TFM image involving the three defects and the
backwall, as shown in Fig. 8. For comparison purposes, the images
3.5. Results and discussion relying on an isotropic medium assumption and on the actual material
properties are also exhibited. The isotropic velocity used to compute the
The statistical analysis was conducted with 100 trials of the opti­ image is 5500 m s− 1, which corresponds to the elasticity constants C11 =
mization algorithms (GAO and PSO) and using the KRR surrogate model C33 = 260 GPa, C13 = 100 GPa and C55 = 80 GPa. With this sound
previously defined. In the case of the GAO, the maximum number of speed, the backwall echo is almost well positioned with respect to its
restarts is set to 10. In the case of the PSO, the swarm is composed of 100 actual position. Both GAO and PSO provide high quality optimized
particles and the optimization stops when the standard deviation of Xj is images. The lateral spreading of the defect echoes is significantly
inferior to 1 GPa for the Cij and 0.1◦ for θ, or when the standard devi­ reduced and the positioning of the defect and backwall echoes is correct,
ation of Yj is inferior to 0.1. The values of the hyper-parameters for the with respect to the reference image computed with the actual parame­
two algorithms are indicated in Table 2. For PSO, the values of α, β, and γ ters. To better appreciate the image improvement in terms of position
were chosen based on the study in [39]. The mean error on the material and width of the echoes, the horizontal and vertical echodynamic curves
property estimates is computed with respect to the actual material are also displayed in Fig. 8. These curves represent the evolution of the

7
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 7. Slowness curves of the Quasi-longitudinal waves and corresponding errors (a): slowness with the actual material properties (green), hidden by GAO results
(blue) and PSO results (red); slowness estimation error with GAO (dashed blue) and with PSO (dashed red). Corresponding actual and optimized material properties
values (b). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

maximum amplitude along the X axis (horizontal echodynamic curve) for the previous numerical study. The inspection configuration remains
and along the Z axis (vertical echodynamic curve). For each case, the unchanged compared to the setup in Fig. 3 (immersion linear array
amplitudes are normalized with respect to the Amax value around side- centered on defect no2 with a water column height of 20 mm) and the
drilled hole no2. The horizontal echodynamic curve is extracted without ROI for the optimization procedure is again a 10 × 10 mm2 area centered
taking into account the backwall echo. Orange dashed lines highlight the on defect no2. The second experimental case is an isotropic ferritic steel
position of defect no2 around which the optimization was performed. covered by an anisotropic 316L stainless steel layer (ρ = 7.96 g.m− 3). As
We observe that both GAO and PSO optimized images look very shown in Fig. 12 the inspection is also conducted in immersion, and the
similar to the reference image computed using the actual material optimization procedure consists in maximizing the image amplitude in a
properties. There is a small shift of the backwall and flaws echoes with 10 × 10 mm2 area centered on defect no1 (1.5 mm diameter), which is
respect to their actual positions. However the increase of the image the closest to the cladding layer.
quality is clearly visible. One also observes that the backwall between The anisotropy of both structures is assumed to be homogeneous and
− 25 mm and − 15 mm along the X dimension is not visible. This is due to only 5 parameters are of interest during the optimization. The elastic
the crystalline structure being tilted of 8◦ around the Y axis. parameters to vary are the same as for the numerical study (C11 , C13 , C33
Quantitatively, the correlations between GAO and reference images, and C55 ). For the V-shape weld, the fifth parameter is the angular shift θ
and between PSO and reference image, are very strong. The value of the of the crystal structure around the Y axis, while for the cladded
correlation coefficient is almost 1.0 in both cases. With respect to the component the fifth parameter is the anisotropic layer thickness, noted
image computed with an isotropic assumption, the amplitude gain of h. All the experimental results presented hereafter are obtained after
both GAO and PSO images is +13 dB for defect no2, i.e., the same averaging 100 GAO trials or 100 PSO trials with 100 particles, and using
amplitude gain obtained with the actual anisotropy parameters (+13 the same values of hyper-parameters as in the numerical case (Table 2).
dB). In terms of the Amax descriptor, both algorithms perform well to The data acquisitions were performed with a multichannel Multi­
produce optimized images. Moreover the positioning errors, given in X++™ (Eddyfi Technologies, Canada) and a 5 MHz linear array (Ima­
Fig. 9, illustrate that both algorithms yield a satisfactory estimate of the sonic, France) composed of 64 elements of 0.5 mm width and with a 0.6
actual position of the flaw. The errors are computed for each one of the mm pitch. The array elements were excited with a rectangular-type
100 trials and averaged. With the GAO, the error is inferior to 1.1 ± 0.9 electric pulse whose time length is 0.2 μs, and signals were digitized
mm for lateral positioning and inferior to 0.7 ± 0.5 mm for depth with a sampling frequency of 100 MHz.
positioning. With the PSO, the error is inferior to 0.4 ± 0.4 mm for A short remark should be added concerning the computation times.
lateral positioning and to 1.1 ± 0.7 mm for the depth positioning. Once With the CIVA software, a TFM image of size 51 × 101 pixels requires
again, these errors are due to the non-uniqueness of the solution and to 0.83s to be computed (50 images per minute). Therefore, with the
the fact that no constraint is set about the echo position during the aforementioned hyper-parameters of the algorithms, one GAO trial lasts
optimization. As previously mentioned, averaging several trials of around 35 min and one PSO trial lasts around 42 min. A complete run of
optimization helps at estimating the actual position of the defect in the 100 GAO then requires 2.5 days, while a complete 100 PSO run needs 3
structure, though negligible uncertainties remain. days. These values are given only for indication purpose as no efforts are
Thereby, the optimization procedure appears to be robust and reli­ made here to reduce the computational cost.
able. It is able to correct the degradation of a TFM image inside an un­
known structure and provide fair estimates of the unknown material
4.1. Optimization of images in an Inconel V-shape weld
properties in the case of homogeneous anisotropy. The accuracy of the
results is satisfactory, given that the image quality description is based
As in the previous section dedicated to simulations, Fig. 10 displays
on only one single scalar criterion. The good convergence of the GAO is
the TFM images computed with isotropic reconstruction model (a), with
assessed when compared to the PSO results: estimation and defect
the actual anisotropy (b), with GAO (c) and with PSO (d). Once more, we
localization are obtained with almost the same values and errors, and
observe that the optimization procedure is able to improve the image
are close to the actual material properties and the associated TFM image.
quality. The image degradation due to the isotropic assumption is cor­
rected and the final images are very similar to the one computed with
4. Experimental results with homogeneous anisotropic samples
the actual material properties. The Amax increase around defect no2,
compared to its value in the isotropic assumption image, is of +10 dB for
Now that the robustness and the reliability of the optimization have
the images obtained with GAO and PSO, which is the same as the image
been discussed, two experimental cases are considered. The first case is
computed with the actual anisotropy (+10 dB). The positioning error of
the Inconel 182 V-shape weld (ρ = 8.26 g.m− 3) that served as template
defect no2 is inferior to 0.6 ± 0.4 mm along the X axis and inferior to 0.8

8
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 8. TFM images with simulated signals in the Inconel V-shape weld: isotropic reconstruction model (a), actual properties (b), gradient-based optimization (c),
particle swarm optimization (d). The black line at the bottom of the images indicates the backwall interface. The red circles represent the side-drilled holes. The
horizontal and vertical curves are echodynamics corresponding to the maximum amplitude along each axis. The orange dashed lines highlight the coordinates of the
maximum amplitude (Amax estimator) for defect no2. Each image is normalized by the maximum amplitude measured in the 10 × 10 mm2 region around side-drilled
hole no2 (i.e., the amplitude of the echo of defect no2 is 0 dB for the four images). (For interpretation of the references to color in this figure legend, the reader is
referred to the Web version of this article.)

± 0.5 mm along the Z axis with GAO. Respectively, the error along X is SNR = 20log10 (Amax/Anoise), where Anoise is the RMS value of the
inferior to 0.5 ± 0.9 mm and to 0.3 ± 0.6 mm along Z with PSO. Once noise (amplitudes below − 12dB with respect to Amax) around the defect
again, the optimization provides high quality images of the inside of the echo; and the Echo Spreading Estimator, defined by ESE = Amax/S− 6dB ,
structure, with uncertainties similar to those obtained with simulated where S− 6dB is the surface of the echo at − 6 dB. From a set of 5000
signals. The values of the estimated weld parameters, along with cor­ images, the aforementioned estimators were extracted. Then, the cor­
responding slowness curves and their errors, are given in Fig. 11. The relation coefficients of the variations of these estimators across the set of
estimation error is of 5.5 ± 7.8% with GAO and of 6.7 ± 8.4% with PSO, points are computed with respect to the Amax criterion. Between Amax
and the errors of the associated slowness curves are inferior to 1.0% and and SNR, the correlation coefficient is 0.99, while between Amax and
3.1% for GAO and PSO, respectively. These results are again comparable ESE, it is 0.94. The variations of SNR and ESE are consequently almost
to those obtained with simulated signals and validate the procedure to perfectly correlated to Amax, and cannot bring more relevant informa­
yield fair estimates of unknown properties in structures with homoge­ tion to reduce the optimization errors. In other words, Amax and SNR
neous anisotropy. can be interchanged when describing the quality improvement in the
From these results, one may ask whether the above estimation errors experimental images.
can be further reduced by using another estimator of image quality that The optimization gives good final results but the echo of defect no3,
includes more information, such as the noise level or the surface of the the deepest one, is always degraded. The image of this defect after
echo at − 6 dB. This is why a complementary study was conducted by optimization around defect no2 was correct with the simulated data (see
defining two other parameters: the Signal-to-Noise Ratio defined by Fig. 8), but not from the experimental point of view because the material

9
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 9. Histogram of the positioning error of the echo of defect no2 with respect to the reference echo computed with the actual anisotropy parameters: GAO (a) and
PSO (b). The blue and orange bars indicate the mean error of positioning along the X and Z axes, respectively, after 100 optimization trials. The dark cyan and dark
red bars indicate the standard deviation of the mean error of positioning along the X and Z axes, respectively. (For interpretation of the references to color in this
figure legend, the reader is referred to the Web version of this article.)

Fig. 10. TFM images with experimental signals in the Inconel V-shape weld: isotropic reconstruction model (a), actual properties (b), gradient-based optimization
(c), particle swarm optimization (d). The optimization is conducted around defect no2. Each image is normalized by the maximum amplitude measured in the 10 ×
10 mm2 region around side-drilled hole no2 (i.e., the amplitude of the echo of defect no2 is 0 dB for the four images).

10
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 11. Slowness curves of the quasi-longitudinal waves and corresponding errors in the experimental V-shape weld case (a): actual material properties (green),
hidden by GAO results (blue) and PSO results (red); slowness estimation error with GAO (dashed blue) and with PSO (dashed red). Corresponding actual and
optimized material properties values (b). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

is not strictly homogeneous. A way of improvement would be to consider along the Z axis with GAO, and errors inferior to 0.1 ± 0.1 mm along the
the medium as a dual layer material. The image of the first two flaws X axis and to 0.53 ± 0.27 mm along the Z axis with PSO. In this case,
would be optimized by varying the properties of the first layer. The third GAO gives a better estimation of the defect position than PSO.
flaw and the backwall would then be addressed in a second optimization The interest of this case is to show that the optimization procedure is
considering the other layer. able to manage unknown parameters of two different natures: four relate
to the elastic behavior of the material, while the fifth relates to the ge­
4.2. Optimization of images in a cladded ferritic steel ometry of the structure. The estimation of the cladding parameters and
the associated slowness curves, together with their corresponding er­
In this experimental case, a ferritic steel is cladded by an anisotropic rors, are provided in Fig. 14. Using GAO, the estimation error of the
layer of 11.5 mm thick. The experimental setup is schematized in material properties is 3.2 ± 6.2% and the slowness estimation error is
Fig. 12. The optimization procedure is applied with the ROI centered less than 1.5%. With PSO, the errors are 5.9 ± 11.9% for the estimation
around flaw no1 at 15.5 mm from the surface. The final images are the cladding parameters and less than 2.3% for the associated slowness
displayed in Fig. 13: computed with an isotropic reconstruction model curve. This efficiency of the optimization is again demonstrated,
(a), with the actual anisotropy (b), with GAO (c) and with PSO (d). These although this last case exhibits higher error rates when estimating the
results once again show the high quality of the final images obtained cladding parameters. This may be explained by the fact that the actual
through the optimization procedure. The values of Amax around defect material properties are subject to measurement uncertainties and may
no1 in the images obtained with GAO and PSO are 10 dB higher than the not be perfectly accurate. Therefore, we can conclude that this test case
value in the image computed with the isotropic assumption, while the also validates the optimization procedure.
Amax value is 9 dB higher in the image computed with the actual ma­ It is worth noticing that the homogeneity assumption to describe the
terial properties. The position of defect no1 is accurately estimated, with whole anisotropic layer works well for the cladded sample, contrary to
errors inferior to 0.1 ± 0.1 mm along the X axis and to 0.17 ± 0.15 mm the V-shape weld case where a dual layer medium could be a better
description of the structure. Also, a shadowing effect is observed for
defects 2–5, and an optimization around one of these defects requires a
different experimental setup (array with a larger aperture or not
centered on the line of defects).

5. Conclusion

The feasibility of an adaptive procedure for TFM imaging in un­

known anisotropic homogeneous structures was investigated in the case
of point-like defects. A numerical study showed that the Gradient Ascent
Optimization proposed herein yields good results, with very satisfactory
accuracy both for imaging and characterization, given the actual ma­
terial properties and the associated image. Its robustness and reliability
were evaluated, although the optimization problem was ill-posed due to
non-uniqueness of the solution. The positions of the echoes were in good
agreement with the actual positions of the reflectors, despite the opti­
mization being only constrained by the Amax image quality descriptor.
The good convergence of the algorithm was also verified by comparing
with the results of a Particle Swarm Optimization.
Fig. 12. Scheme of the experimental setup: immersion inspection with a Experiments conducted with two nuclear steel components validated
transducer array placed 20 mm above the surface; ferritic-based steel mockup the procedure and exhibited good final results. It was however observed
with 5 side drilled holes at depths 15.5 mm (1), 23.5 mm (2), 31.5 mm (3), that the homogeneity assumption can be too strong in some cases.
39.5 mm (4) and 47.5 mm (5) from the surface and of 2.0 mm diameter. Further developments are needed to extend the method to more complex
Anisotropic cladding of 11.5 mm thickness and assumed to be homogeneous. materials, exhibiting inhomogeneous anisotropic properties, and to

11
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 13. TFM images with experimental signals in the cladded ferritic steel component: isotropic reconstruction model (a), actual properties (b), gradient-based
optimization (c), particle swarm optimization (d). The optimization is conducted around defect no1. Each image is normalized by the maximum amplitude
measured in the 10 × 10 mm2 region around side-drilled hole no1 (i.e., the amplitude of the echo of defect no1 is 0 dB for the four images).

12
C. Ménard et al. NDT and E International 116 (2020) 102327

Fig. 14. Slowness curves of the quasi-longitudinal waves and corresponding errors in the experimental cladded component case: actual material properties (green),
hidden by GAO results (blue) and PSO results (red); slowness estimation error with GAO (dashed blue) and with PSO (dashed red). Corresponding actual and
optimized material properties values (b). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

welds with various chamfer geometries. It is also planned to extend this Annual Review of Progress in Quantitative Nondestructive Evaluation, Boise, USA.
AIP Conf. Proc., vol. 1650; July 2014. p. 1657–66. 2015.
study to the TFM imaging of crack-like defects or lacks of fusion that are
[7] Zhang J, Drinkwater BW, Wilcox PD, Hunter AJ. Defect detection using ultrasonic
typically found on the weld chamfers. arrays: the multi-mode total focusing method. NDT Int 2010;43(2):123–33.
[8] Gengembre N, Lhémery A. Pencil method in elastodynamics: application to
ultrasonic field computation. Ultrasonics 2000;38:495–9.
CRediT authorship contribution statement
[9] Pudovikov S, Bulavinov A, Pinchuk R. Innovative ultrasonic testing (UT) of nuclear
components by sampling phased array with 3D visualization of inspection results,
Corentin Ménard: Conceptualization, Software, Writing - original JRC-NDE 2010 8th international conference on NDE in relation to structural
draft, Validation, Methodology. Sébastien Robert: Conceptualization, integrity for nuclear and pressurized components. Berlin, Germany. In: JRC
Scientif. Techn. Reports JRC 64886; April 2010. p. 910–7. 2011.
Validation, Methodology, Supervision. Roberto Miorelli: Validation, [10] Juengert A, Dugan S, Homann T, Mitzscherling S, Prager J, Pudovikov S,
Supervision, Formal analysis. Dominique Lesselier: Supervision, Schwender T. Advanced ultrasonic techniques for nondestructive testing of
Validation. austenitic and dissimilar welds in nuclear facilities, 44th annual review of progress
in quantitative nondestructive evaluation. Provo, USA. AIP Conf. Proc., vol.
110002; July 2017. 1949, 2018.
Declaration of competing interest [11] Chassignole B, Recolin P, Leymarie N, Gueudré C, Guy P, Elbaz D. Study of
ultrasonic characterization and propagation in austenitic welds: the MOSAICS
project, 41st annual review of progress in quantitative nondestructive evaluation,
The authors declare that they have no known competing financial boise, USA. AIP Conf. Proc., vol. 1650; July 2014. p. 1486–95. 2015.
interests or personal relationships that could have appeared to influence [12] Tant KMM, Galetti E, Curtis A, Mulholland AJ, Gachagan A. A transdimensional
the work reported in this paper. bayesian approach to ultrasonic travel-time tomography for nondestructive testing.
Inverse Probl 2018;34(9). 095002.
[13] Fan Z, Mark AF, Lowe MJS, Withers PJ. Nonintrusive estimation of anisotropic
Acknowledgment stiffness maps of heterogeneous steel welds for the improvement of ultrasonic array
inspection. IEEE Trans Ultrason Ferroelectrics Freq Contr 2015;62(8):1530–43.
[14] Zhang J, Hunter AJ, Drinkwater BW, Wilcox PD. Monte Carlo inversion of
This study is part of European project ADVISE, that has received ultrasonic array data to map anisotropic weld properties. IEEE Trans Ultrason
funding from the Euratom research and training programme 2014–2018 Ferroelectrics Freq Contr 2012;59(11):2487–97.
under grant agreement No 755500. The authors would like to thank the [15] Stojakovic D. Electron backscatter diffraction in materials characterization. Process
Appl Ceram 2012;6(1):1–13.
ADVISE project partners, with special credits to Electricité De France
[16] Mark AF, Li W, Sharples S, Withers PJ. Comparison of grain to grain orientation
(EDF) for providing mock-ups to all consortium partners. and stiffness mapping by spatially resolved acoustic spectroscopy and EBSD.
J Microsc 2017;267:89–97.
[17] An S, Liu W, Venkatesh S. Face recognition using Kernel Ridge regression,
Appendix A. Supplementary data
conference on computer vision and pattern recognition, Minneapolis, USA. In: Ieee
Conf. Comput. Vis. Pattern Recog; June 2007. p. 1–7. 2007.
Supplementary data to this article can be found online at https://doi. [18] Kennedy J, Eberhart R. Particle swarm optimization, int. Conf. On neural networks,
org/10.1016/j.ndteint.2020.102327. perth, Australia, November 1995. Ieee Int. Conf. On neural networks Proc., vol. 4;
1995. p. 1942–8.
[19] Clerc M, Kennedy J. The particle swarm - explosion, stability, and convergence in a
References multidimensional complex space. IEEE Trans Evol Comput 2002;6(1):58–73.
[20] Lopez Villaverde E, Robert S, Prada C. Ultrasonic imaging of defects in coarse-
[1] Holmes C, Drinkwater BW, Wilcox PW. The post-processing of ultrasonic array data grained steels with the decomposition of the time reversal operator. J Acoust Soc
using the total focusing method. Insight Non-Destruct Test Condit Monit 2004;46 Am 2016;140(1):541–50.
(11):677–80. [21] Gardahaut A, Jezzine K, Cassereau D, Leymarie N. Simulation of ultrasonic wave
[2] Prada C, Fink M. Eigenmodes of the time reversal operator: a solution to selective propagation in welds using ray-based methods, 12th Anglo-French Physical
focusing in multiple-target media. Wave Motion 1994;20:151–63. Acoustics Conference, Fréjus, France, January 2013. J. Phys.: Conf. Ser., vol. 498;
[3] Fink M, Cassereau D, Derode A, Prada C, Roux P, Tanter M, Thomas JL, Wu F. 2014. p. 1–7.
Time-reversed acoustics. Rep Prog Phys 2000;63(12):1933–95. [22] Loshlilov I, Hutter F. SGDR: stochastic gradient descent with warm restarts, 5th
[4] Le Jeune L, Robert S, Dumas P, Membre A, Prada C. Adaptive ultrasonic imaging international conference on learning representations. France, Toulon, April 2017
with the total focusing method for inspection of complex components immersed in arXiv 2017:1608.03983.
water, 41st Annual Review of Progress in Quantitative Nondestructive Evaluation, [23] Ruder S. An overview of gradient descent optimization algorithms. arXiv 2016:
Boise, USA. AIP Conf. Proc., vol. 1650; July 2014. p. 1037–46. 2015. 1609.04747.
[5] Hoyle E, Sutcliffe M, Charlton P, Rees J. Virtual source aperture imaging with auto- [24] Bengio Y. Practical recommendations for gradient-based training of deep
focusing of unknown complex geometry through dual layered media. NDT E Int architectures. In: Neural networks: Tricks of the trade. second ed. Springer Berlin
2018;98:55–62. Heidelberg; 2012. p. 437–78.
[6] Robert S, Calmon P, Calvo M, Le Jeune L, Iakovleva E. Surface estimation methods [25] Smith LN. No more pesky learning rate guessing games. arXiv 2015:1506.01186.
with phased-arrays for adaptive ultrasonic imaging in complex components, 41st

13
C. Ménard et al. NDT and E International 116 (2020) 102327

[26] Darmon M, Leymarie N, Chatillon S, Mahaut S. Modelling of scattering of [32] Saltelli A, Ratto M, Andres T, Campolongo F, Cariboni J, Gatelli D, Saisana M,
ultrasounds by flaws for NDT, ultrasonic wave propagation in non homogeneous Tarantola S. Global sensitivity analysis, the Primer. John Wiley and Sons; 2008.
media. Springer Berlin Heidelberg; 2009. p. 61–71. chap. 4.
[27] Royer D, Dieulesaint E. Elastic waves in solids I, Free and guided propagation. [33] Herman J, Usher W. SALib: an open-source Python library for sensitivity analysis.
Springer-Verlag Berlin Heidelberg; 2000. chap. 4. J Open Sourc Software 2017;2(9).
[28] Pedregosa F, Varoquaux G, Gramfort A, Michel V, Thirion B, Grisel O, Blondel M, [34] Li G, Rabitz H, Yelvington PE, Oluwole OO, Bacon F, Kolb CE, Schoendorf J. Global
Prettenhofer P, Weiss R, Dubourg V, Vanderplas J, Passos A, Cournapeau D, sensitivity analysis for systems with independent and/or correlated inputs. J Phys
Brucher M, Perrot M, Duchesnay E. Scikit-learn: machine learning in Python. Chem A 2010;114(19):6022–32.
J Mach Learn Res 2011;12:2825–30. [35] El Dor A, Clerc M, Siarry P. A multi-swarm PSO using charged particles in a
[29] Tibshirani R. Regression shrinkage and selection via the Lasso. J Roy Stat Soc B partitioned search space for continuous optimization. Comput Optim Appl 2012;
1996;58(1):267–88. 53:271–95.
[30] Saunders C, Gammerman A, Vovk V. Ridge regression learning algorithm in dual [36] Eberhart R, Simpson P, Dobbins R. Computational intelligence PC tools. Academic
variables, International Conference on Machine Learning. Madison, USA. In: Press Professional; 1996.
Proceedings of the 15th International Conference on Machine learning; July 1998. [37] Chatterjee A, Siarry P. Nonlinear inertia weight variation for dynamic adaptation
p. 515–21. 1998. in particle swarm optimization. Comput Oper Res 2006;33:859–71.
[31] Géron A. Hands-on machine mearning with Scikit-learn, Keras, and TensorFlow. [38] Talbi EG. A taxonomy of hybrid metaheuristics. J Heuristics 2002;8:541–64.
second ed. O’Reilly Media; 2019. chap. 4. [39] Trelea IC. The particle swarm optimization algorithm: convergence analysis and
parameter selection. Inf Process Lett 2003;85:317–25.

14