Detecting Delamination via Nonlinear Wave Scattering in a Bonded Elastic Bar

In this paper we examine the effect of delamination on wave scattering, with the aim of creating a control measure for layered waveguides of various bonding types. Previous works have considered specific widths of solitary waves for the simulations, without analysing the effect of changing the soliton parameters. We consider two multi-layered structures: one containing delamination"sandwiched"by perfect bonding and one containing delamination but"sandwiched"by soft bonding. These structures are modelled by coupled Boussinesq-type equations. Matched asymptotic multiple-scale expansions lead to coupled Ostrovsky equations in soft bonded regions and Korteweg-De Vries equations in the perfectly bonded and delaminated region. We use the Inverse Scattering Transform to predict the behaviour in the delaminated regions. In both cases, numerical analysis shows that we can predict the delamination length by changes in the wave structure, and that these changes depend upon the Full Width at Half Magnitude (FWHM) of the incident soliton. In the case of perfect bonding, we derive a theoretical prediction for the change and confirm this numerically. For the soft bonding case, we numerically identify a similar relationship using the change in amplitude. Therefore we only need to compute one curve to determine the behaviour for any incident solitary wave, creating a framework for designing measurement campaigns for rigorously testing the integrity of layered structures.


Introduction
Solitary waves are of significant interest from both a mathematical perspective as well as in physical and engineering applications.They often arise as solutions to nonlinear equations such as the KdV equation (and its extensions) in shallow water [1][2][3][4], the Benjamin-Ono equation for internal waves of stratified fluids [5,6], the nonlinear Schrödinger equation in optics [7,8], and flexural waves in the beam equation [9], to name a few.From a purely mathematical perspective, there are many studies into the existence and behaviour of solitons, for example as solutions to the Boussinesq equation [10,11].Boussinesq-type equations are of interest in this study, in the context of solid mechanics.It is well-known that they can describe long longitudinal bulk strain waves in elastic waveguides, such as rods and metal plates (see e.g.[12][13][14][15][16][17]). Practical experiments have confirmed that longitudinal bulk strain solitons exist in these waveguides, which validated theoretical findings [18][19][20][21].
Indeed, layered waveguides with bonding between the layers can be modelled by the socalled "doubly dispersive equation" (DDE), which can be derived using nonlinear elasticity theory for long longitudinal waves in a bar of rectangular cross-section [12,13].The DDE for a bar of rectangular cross-section σ = 2a × 2b has the form where ρ is the density, E is the Young's modulus, v is the Poisson's ratio, while l, m, n are the Murnaghan's moduli, and a and b are geometric and physical parameters.
The case when the interlayer bonding is missing over part of the structure, known as delamination, is important for a wide range of applications in non-destructive testing for structural damage in the multi-layer beam-like structures found throughout civil and mechanical engineering.The governing mathematical model then takes the form of a scattering problem and, for a perfectly bonded waveguide (represented in experiments by cyanoacrylate), we find a series of Boussinesq equations, with continuity conditions on the interface between sections [22].Incident solitons fission into multiple solitons with dispersive radiation, agreeing with theoretical predictions [22], numerical simulations [23,24] and experimental observations [20,25,26].
In the case of an imperfect "soft" bonding (represented in experiments by polychloroprene), a model based upon a series of anharmonic coupled dipoles can be used to derive coupled regularised Boussinesq (cRB) equations to model long nonlinear longitudinal bulk strain waves in a bi-layer, assuming sufficiently soft bonding [27].In this case, when the materials in the layers are sufficiently close, an incident solitary wave evolves in the bonded region into a solitary wave with a one-sided, co-propagating oscillatory tail, known as a radiating solitary wave.In the delaminated regions, the solitary wave detaches from its tail and this can be used as a measure of delamination [28].More recently, in the limiting case of a semi-infinite delamination, we find the emergence of Ostrovsky wave packets in bonded regions [29].The Ostrovsky equation originally arose in the context of shallow-water waves, where the rotation of the Earth is considered [30], and the evolution of wave packets generated from an initial pulse has been extensively studied [31].
In this paper we aim to use theoretical predictions and numerical simulations to establish a prediction for the delamination length based upon changes in the wave during its propagation.We will consider a range of initial conditions by varying the Full Width at Half Magnitude (FWHM) of the incident wave, whereas previous studies have only considered a single fixed incident soliton [23,28].Our aim is to find a relationship between the generated delamination curves for different values of FWHM, so that only one curve needs to be computed, significantly reducing the computation time.This wider range of predictions allows for the design of measurement campaigns for detecting and measuring delamination in layered waveguides.We will consider a multi-layered symmetric structure with perfect bonding, as well as a two-layered structure with soft bonding.In both cases, we will consider delamination 'sandwiched' by bonding.These structures are illustrated in figures 1 and 2, and are inspired by an existing experimental set-up [20].
The paper is structured as follows.In Section 2 we introduce the equations describing longitudinal wave propagation in both the perfectly bonded and soft bonded cases.We also introduce the weakly-nonlinear solution for the perfectly bonded case so that we can create a measure of the delamination length using theoretical predictions.In Section 3, we begin by illustrating the evolution of incident solitary waves in the cases of both a perfectly bonded waveguide and a soft bonded waveguide.Next, for the perfectly bonded case, we use theoretical predictions to determine the length of the delamination region for a variety of incident solitary waves, tested via numerical simulations and measuring the amplitude of the transmitted wave with reference to the incident wave.This gives rise to a relationship between FWHM of the incident soliton and the delamination length, allowing for efficient computation for other incident waves.A similar result is also found for the soft bonded case, analysing the decrease in amplitude after the wave propagates through a delaminated region.The theoretical prediction is difficult in this case, so we rely on numerical observations and instruction from the previous case.In Section 4 we conclude our discussions.
2 Problem set-up 2.1 Perfectly bonded case Fig. 1: Bi-layer structure with an initial perfect bonded region for x 0 < x < x 1 , a delaminated region for x 1 < x < x 2 and a perfect bonded region for x 2 < x < x 3 .We assume that the materials in both layers are identical.
We consider the scattering of a long longitudinal strain solitary wave in a perfectly bonded layered bar with delamination in the centre, as illustrated in Figure 1.Note that we only illustrate two layers in this figure, but the setup can accommodate any number of symmetric layers, as we assume that the materials of the layers are identical and that the bonding is the same between all layers.This problem is described by the regularized non-dimensional Boussinesq equations [22] u where i = 1, 3 represent the perfect bonded regions and i = 2 represents the delaminated region.We have the coefficients α i , β i and c i , which can theoretically vary between sections (representing a waveguide with different materials in each section), but for our purposes we will assume that the sections are of one and the same material.The parameter ε is the small wave parameter.The Boussinesq equations are complemented with continuity conditions, namely continuity of longitudinal displacement and continuity of normal stress where σ (i) is defined by our original equation (3) when written in the form We consider α i = 1 for all i, β 1,3 = 1 and where n represents the number of layers in the structure and k is defined by the geometry of the waveguide.Referring to Figure 1, the cross section has width 2a and the height of each layer is 2b/n.In terms of these values, as there are two layers, n = 2 and k = b/a.In our numerical simulations we will consider various n and k values.

Weakly-nonlinear solution
In order to find theoretical predictions for the evolution of the solitary waves, we construct a weaklynonlinear solution and use theoretical results for the derived equations.For brevity, we only provide a summary of the results below, more details can be found in Refs.[23,24,28,29].We seek a weaklynonlinear solution for the strains x of the form where ξ = x − c i t, η = x + c i t and X = εx.Substituting the respective weakly-nonlinear solutions into the differentiated form of (3), then applying space-averaging (see [23,28,29]) yields the leading order solutions To determine "initial conditions" for the equations derived in each section, we substitute (7) into the continuity conditions ( 4) -( 5) to find values for T and R at the interface, in terms of the previous transmitted wave.This gives rise to transmission and reflection coefficients in terms of c i .As we have assumed that the waveguide is one and the same material (so c i = 1 for all i) we have full transmission and no reflection.

Theoretical predictions
We can rewrite the transmitted wave equation ( 8) in the form For a sufficiently rapidly decaying initial condition U 0 (χ) on the infinite line, the solution to ( 10) is related to the spectral problem for the Schrödinger equation where λ is the spectral parameter.Finding the evolution of the scattering data for the discrete and continuous spectra and using these to reconstruct the solution to the KdV equation is known as the Inverse Scattering Transform (IST) [32].We can use the results from the IST to create theoretical predictions for the solitons in the delaminated region, as well as in the second bonded region.We assume that there is an incident soliton in the first region, which is a travelling wave solution and thus will move in time but retains its shape.To illustrate the theoretical predictions we consider the second region, where we have β 2 defined as in (6) and the initial condition for (10) in this region then takes the form (12) In this case the solution will consist of either one soliton, or a series of solitons, characterised by eigenvalues in the discrete spectrum, and accompanying dispersive radiation determined by the continuous spectrum.In some cases we may see the fission of the initial soliton, which is when more than one soliton is generated, in particular when β 2 ̸ = 1.
The discrete eigenvalues of (11) take the form where for n = 1, 2, . . ., N .Recalling ( 6), the number of solitons, N , generated in the delaminated region is given by the largest integer satisfying the inequality We can see from ( 14) that, for β = 1 we will have one soliton, while for β < 1, we will have more than one soliton and as β becomes smaller, more solitons will be generated.This corresponds to either an increase in layers in the waveguide, or a change in geometry.As τ → ∞, the solution will evolve into a train of solitary waves, ordered by their heights, propagating to the right and some dispersive radiation (a dispersive wave train) propagating to the left (in the moving reference frame), i.e.
where χ n is the phase shift.In the context of our problem, if there is an infinite delamination then the solitons will separate and rank order, while for finite delamination the solitons will only separate for a large delamination.This allows us to create a measure of the delamination length, by comparing the measured signal at the end of the bar to the theoretical prediction.We introduce the incident solitary wave for T (1) , the exact travelling wave solution of (8), as where v is the phase speed.Solitary waves are often measured in experiments in terms of their Full Width at Half Magnitude (FWHM), so we rewrite this as and hence we obtain This allows us to generalise the FWHM based measure to any size of incident solitary wave.

Imperfect bonding case
The second case we consider is when we have a two layered waveguide with "soft" bonding between the layers.This is illustrated in Figure 2. The Fig. 2: Bi-layer structure with an initial soft bonded region for x 0 < x < x 1 , a delaminated region for x 1 < x < x 2 and a soft bonded region for x 2 < x < x 3 .We assume that the materials in both layers are similar, that is, their material properties differ by O(ε).
longitudinal displacement in the bonded regions is described by the regularized non-dimensional equations for x i−1 < x < x i , while in the delaminated regions we have Boussinesq equations As with the perfectly bonded case, these equations are complemented with continuity conditions at the interfaces between the sections.We have continuity of longitudinal displacement and continuity of normal stress for i = 1, 2, where σ u and σ w are defined by ( 19) and ( 20) as respectively.We will consider the case here where the materials in the layers are similar, namely c−1 = O(ε).We can construct a weakly-nonlinear solution to this system of equations, as was done in [28], however we cannot obtain any direct theoretical predictions from this approach as the derived coupled Ostrovsky equations are not solvable via the Inverse Scattering Transform.Therefore, we will explore this case numerically to determine a measure of delamination.

Numerical results
We now aim to use the derived weakly-nonlinear solution and the theoretical predictions of Section 2 to introduce a measure of the delamination length in terms of the change in wave structure.
In this section, we first demonstrate the effect of delamination on the transmitted soliton for both the perfectly bonded and soft bonded waveguides in Section 3.1.We then introduce a measure of the delamination length for the perfectly bonded case in Section 3.2, and for the soft bonded case in Section 3.3.In both cases we consider how these measures scale with respect to the incident soliton in order to rapidly recompute results for a wide range of initial conditions.We will use the finite difference scheme from Ref. [24] to solve the original Boussinesq system and a semi-analytical method using a pseudospectral scheme for the derived KdV equations for the perfectly bonded case, similar to the one used for coupled Ostrovsky equations in [28].In all cases our simulations will use a grid spacing of ∆x = 0.01 and a time step of ∆t = 0.01 for the finite-difference scheme.For the semi-analytical method we take ∆X = 5 × 10 −4 and ∆ξ = 0.1, corresponding to N = 131, 072.

Examples of scattering
Firstly, we demonstrate the effect of delamination on the propagation of an incident solitary wave, in both scenarios described in Section 2. For the perfectly bonded case, let us assume a spatial domain x ∈ [−100, 1000] and we have a delamination starting at x = 0 of length D. We present three cases: no delamination, a delamination of length D = 50 and a delamination of length D = 300.These results, as well as a comparison between the weakly-nonlinear solution and the directly computed solution, are shown in Figure 3.We can see from Figure 3a that, in the case of no delamination, the soliton propagates without change in shape or structure.When delamination is introduced in figures 3b and 3c, the soliton fissions into two solitons with accompanying dispersive radiation.Comparing the case of D = 50 to D = 300, we can see that the second soliton has become more distinct from the radiation and the first soliton has tended towards its theoretically predicted amplitude.Figure 3d shows that there is a reasonable agreement between the direct numerical simulation and the semi-analytical method when D = 300, with a slight phase shift and amplitude change.This will be reduced for smaller values of ε (as we have constructed a series in increasing powers of ε), however we have qualitatively the same structure.
Similarly, we consider the case of soft bonding and an incident soliton.We assume a spatial domain x ∈ [−500, 1000], with a homogeneous (delaminated) section for x ∈ [−500, −400] to generate the same wave in both layers, bonded sections x ∈ [−400, 0] and x ∈ [D, 1000], with a delaminated region for x ∈ [0, D].We again present three cases including when there is no delamination, a shorter delamination length of   The results are shown in Figure 4 for the upper layer, where the lower layer is qualitatively similar.In the case of no delamination shown in Figure 4a we have a solitary wave with a one-sided oscillatory tail, known as a radiating solitary wave.As the delamination length increases, the solitary wave begins to lose amplitude and expel energy into its tail through an exchange of energy between the layers.These are clear signs of the presence of delamination in figures 4b and 4c, with structural changes that can be detected and the decay in amplitude can be quantified to give a measure of the length of delamination.Note that in this case there is no comparison between the simulation schemes since the semi-analytical scheme is not applicable for the soft bonded case.

Measure of delamination length for perfect bonding
We now expand upon the observations from the previous subsection by introducing a measure of delamination based upon the theoretical predictions from Section 2. We then generalise this measure to apply for different incident soliton widths.
If we assume that the solitons are well-separated, representing the case of infinite delamination, the amplitude of the soliton can be found from the IST.We have where A 1 is the amplitude of the incident soliton, A 3 is the amplitude of the lead soliton in the second bonded region, and k 2 , k 3 are the eigenvalues corresponding to the lead soliton amplitude in the second and third regions, as determined by the IST.As the delamination length is reduced, the amplitude in the third region will tend towards the initial amplitude, A 1 .
Denoting the calculated numerical solution as A num from the simulation, we introduce a measure of the amplitude of the lead soliton in the third section of the bar in comparison to the incident soliton as This corresponds with the measure used in [24], where it was assumed that FWHM = 5.We now compute the solution using the semi-analytical pseudospectral scheme for a wide range of values of FWHM with the aim of determining a general rule for the delamination length.The delamination length is chosen to be D ∈ [0, 300] and we measure the delamination in multiples of FWHM.
The results are plotted in Figure 5a for n = 3, k = 3, and in Figure 5b for n = 4 and k = 3.We can see that, as the delamination length increases, the measure σ increases until it tends to a value of 1, and this behaviour is replicated for different values of FWHM.For larger FWHM it may not reach this limit for the chosen delamination length.We can also see a similar behaviour for different values of n and k.To generalise this approach, we consider the equation for the phase speed v in terms of FWHM (18).We can see that v in inversely proportional to the square of FWHM.Thus, fixing our reference as FWHM = 5, we let σ be a function of delamination length, parametrised by FWHM, and we introduce the scaling The resulting plots are shown in Figure 6a for n = 3, k = 3, and in Figure 6b for n = 4, k = 3.We see that the scaled versions are very closely aligned, with the only disagreement stemming from the restriction on delamination length for larger values of FWHM.Therefore, after computing one curve for the smallest value of FWHM, we can reproduce all subsequent curves for any value of FWHM.This allows for the highly efficient computation of the delamination curves and for a wide range of experiments with incident waves of different amplitude.

Finite delamination with soft bonding
We now consider the soft bonded case outlined in Section 2.2.This case was also studied in Ref. [28], but for only one value of FWHM.In this work we compute the solution for a wide range of FWHM, using the finite-difference scheme in [24].The constructed weakly-nonlinear solution consists of coupled Ostrovsky equations, in contrast to the KdV equations in the previous case [28].Therefore we cannot use the IST to predict the amplitude of the waves in the bonded regions as coupled Ostrovsky equations are not integrable via the IST.The incident soliton in this case evolves into a radiating solitary wave, that is a solitary wave with a one-sided oscillatory tail.
To determine the change in amplitude, as a measure of the delamination length, we denote the amplitude of the soliton or wave packet in each region as A L , where L is the region index, and we introduce the measure Figure 7(a) presents the results for ζ, computed over a wide range of FWHM values, for the upper layer.A similar agreement is seen for the lower layer, but the results are omitted for brevity.Following the idea from the perfectly bonded case, we introduce a quadratic scaling.We choose a reference FWHM value (in this case we choose our lowest value of FWHM, namely FWHM = 5) and then calculate a scaling of the form where we have introduced We now summarise the scaling used to convert our nondimenesional variables to dimensional material parameters that can be compared to experimental data, in order to confirm whether a delamination length of 200 is reasonable.Referring to the dimensional form of the DDE in (1), with parameters (2), we introduce the scaling to nondimensional form via where We can therefore find the corresponding material length given the nondimensional length.Let us assume a PMMA bar of 10mm × 10mm crosssection, then using the parameters for PMMA from [20] we find that, for ε = 0.1, a delamination length of x = 200 in nondimensional units corresponds to a length of approximately x = 520mm, which is significant given the experimental materials are usually around 600mm long in total.Therefore, restricting our considerations to values of delamination less than 200 nondimensional units is reasonable in the context of practical applications.

Conclusion
In this paper we have considered the scattering of a bulk strain solitary wave in a delaminated bilayer with either perfect or soft bonding between the layers.The longitudinal displacements are modelled by either Boussinesq equations (perfectly bonded or delaminated sections) or coupled Boussinesq equations (soft bonded sections), with continuity conditions on the interface.Incident solitary waves undergo fission in delaminated regions in the perfectly bonded structure, providing a clear indicator of delamination.We construct theoretical estimates for the amplitude of the solitons after a delaminated region, using the Inverse Scattering Transform.A measure is introduced using the theoretical and observed values to predict the delamination length based upon amplitude changes.This is then extended for incident waves of different Full Width at Half Magnitude, and a quadratic scaling is introduced and verified by numerical results.The significance of this result is that we now only need to compute a single curve in order to perform a wide range of experiments, which significantly reduces computation times and allows for further experiments (with different incident solitons) to be performed rapidly.This was confirmed for various configurations of the waveguide.
In the case of a soft bonded waveguide with delamination, theoretical estimates cannot be derived using the Inverse Scattering Transform.In this case we resort to direct computation of the solution and a comparison between the amplitude after delamination and the corresponding amplitude for the non-delaminated case.A similar quadratic scaling can be found, which has a good agreement up to a delamination length of 200 in nondimensional units or 520mm in physical units.This is consistent for both layers of the waveguide.
This work facilitates experimentation with a wide range of initial condition parameters, and provides a framework for detecting delamination in perfectly bonded and soft bonded waveguides with similar materials in the layers.The case with distinctly different materials in the layers is more complex, and some preliminary studies have been conducted into quantifying delamination [29].

Declarations
• Funding: Jagdeep S. Tamber would like to thank Nottingham Trent University for funding through their PhD studentship scheme.• Competing interests: The authors have no relevant financial or non-financial interests to disclose.• Ethics approval: Not applicable.
• Consent to participate: Not applicable.
• Consent for publication: Not applicable.
• Availability of data, code and materials: The datasets generated during this study can be reproduced using equations throughout the paper and the cited numerical methods.The codes used to generate the datasets are available from the corresponding author on reasonable request.
Comparison of schemes for D = 300.

Fig. 3 :
Fig.3: The solution at t = 900 in the final section of the perfectly bonded waveguide, for various delamination lengths, and comparison of the direct numerical (blue, solid line) and semianalytical (red, dotted line) simulations.Parameters are ε = 0.1, FWHM = 5.0, n = 2 and k = 2.The finite-difference method uses a computational domain of [−100, 1000] and for the semi-analytical pseudospectral method we have N = 131, 072.

Fig. 5 :
Fig. 5: Plots of the change in amplitude of the lead transmitted soliton, in comparison to the incident soliton and theoretical predictions, as measured by σ, for various values of FWHM.Here we take ε = 0.1.

Fig. 6 :
Fig. 6: Plots of the scaled delamination curves for the change in amplitude of the lead transmitted soliton, in comparison to the incident soliton and theoretical predictions, as measured by σ.Here we take ε = 0.1.