Wave Motion 94 (2020) 102526
Contents lists available at ScienceDirect
Wave Motion journal homepage: www.elsevier.com/locate/wamot
Numerical simulation of ultrasonic time reversal on defects in carbon fibre reinforced polymer ∗
M. Lints b,a , , A. Salupere a , S. Dos Santos b a b
Department of Cybernetics, School of Science, Tallinn University of Technology, Akadeemia tee 21, 12618 Tallinn, Estonia INSA Centre Val de Loire, Blois Campus, 3 Rue de la Chocolaterie CS 23410, F-41034 Blois Cedex, France
article
info
Article history: Received 15 January 2019 Received in revised form 18 December 2019 Accepted 28 January 2020 Available online 31 January 2020 Keywords: Contact acoustic nonlinearity Resonance excitation Ultrasound Delayed time reversal Carbon fibre reinforced polymer
a b s t r a c t An ultrasonic non-destructive testing method is under development for detecting small scale damage in carbon fibre reinforced polymer. This method relies on signal processing to detect the presence of nonclassical nonlinear signature of microcracking or delamination in the material. Using reciprocal time reversal with nonlinear spectroscopy, the ultrasound can be focused on defect which is in unknown location in test object and may be far from the sending or the receiving transducer. With spectral analysis of defect and delayed time reversal input, the defect can be excited to larger displacement amplitudes by using its resonance. The study is conducted using finite element simulations in two dimensions for carbon fibre composite plate. © 2020 Published by Elsevier B.V.
1. Introduction Carbon Fibre Reinforced Polymer (CFRP) has gained widespread use in modern manufacturing from non safety-critical applications in advanced technologies to recently use in safety critical parts of consumer products. Still, its higher cost of manufacturing compared to metal parts and the possible complex damage mechanisms necessitate the development of robust, simple and affordable Non-Destructive Testing (NDT) methods. Microcracking and delamination are among many defects that are hard to detect using cheap ultrasonic apparatus. The recent advances of Structural Health Monitoring (SHM) applications place further constraints on recently developed ultrasonic NDT methods such as scalability, simple and low-weight apparatus and global method of testing the object, as opposed to simple scanning of material for defects. Nonclassically nonlinear defects in linear materials have been studied for detection methods for decades: McCall and Guyer have developed a theory for hysteretic nonlinear materials [1] and with Meegan and Johnson have conducted respective experiments [2] on the influence of small scale damage on classical and nonclassical nonlinear effects of solid materials. One type of nonlinear damage is a crack which can close intermittently during excitation, and can be used to recognise and image defects by [3–6], including resonance effects [7,8]. Previous authors [9–12] have examined the nonlinear spectroscopy of complex materials for detection of small-scale defects which can be difficult to find using traditional linear ultrasound methods. In case of CFRP, the nonlinear methods could prove useful for detecting smaller than wavelength defects like delamination and microcracking, which are typical failure modes of the laminate. Ciampa and Meo [13] have studied detecting the nonlinear signature of damaged complex solids using reciprocal Time Reversal (TR). They used a simple setup consisting of one broadcasting and one receiving transducer. By assuming hysteretic behaviour of damaged material, and extracting the third order nonlinear signature (by ∗ Corresponding author at: Department of Cybernetics, School of Science, Tallinn University of Technology, Akadeemia tee 21, 12618 Tallinn, Estonia. E-mail addresses:
[email protected] (M. Lints),
[email protected] (A. Salupere),
[email protected] (S. Dos Santos). https://doi.org/10.1016/j.wavemoti.2020.102526 0165-2125/© 2020 Published by Elsevier B.V.
2
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
employing Nonlinear Elastic Wave Spectroscopy — NEWS), they were able to detect the damaged location of a composite test object. This paper studies an ultrasonic NDT method based on reciprocal Time Reversal (TR) signal processing used to focus the ultrasonic waves, using a simple test setup with advanced signal processing. Since it relies on the internal reflections of the material to act as additional ‘‘virtual transducers’’ used in the focusing process, it is advantageous in case of materials with complex geometry or internal structure, such as composites [13–15]. Furthermore, for potential SHM applications, the method should be (semi-)global, meaning that the broadcasting or receiving transducer need not be located exactly on the defect. The main idea of this paper is exciting a nonlinear contacting defect resonance by delayed TR [16] (which is a novel method based on modification of the reciprocal TR process). To do this, the main steps are to firstly excite the medium by a broadband chirp and conduct PI analysis to extract the nonlinear component [17] (although other NEWS methods could be used). Then the received signals from simple chirp and PI focused waves are compared spectrally, to find the potential resonance frequency of the defect (assuming that the PI focused wave excites the defect resonance frequencies more than simple chirp). It is then straightforward to transmit the PI focusing signal with the period of suspected defect resonance period. This is done by delayed TR principles. The paper is structured as follows. Second section describes the simulation model and contains the theoretical background of the used signal processing steps used, including reciprocal TR, TR-NEWS and delayed TR. These can be used for detection, localisation and excitation of the defect globally in a single-channel ultrasonic testing configuration. The third section describes the simulation results, analysis of these results, spectroscopy and the effect on resonance enhancement. 2. Theoretical background 2.1. Reciprocal time reversal The traditional direct TR methods usually require large transducer arrays, where acoustic waves from the source transducer are captured by an array, time reversed and then transmitted from previously receiving array to the initial source, to reconstruct the initial focusing, as developed by Fink [18]. In case of a reflecting or scattering medium with low absorption, it is possible to replace the short recording of the transducer array with a long recording of a single transducer, making the method suitable for single-channel setup [15,19]. Furthermore, the medium heterogeneity and the time reversal window length T are linked to the focusing amplitude and quality improvement. A conceptually different TR method, called reciprocal TR has both ultrasonic transmissions happen in the same direction without exchanging the roles of the transducers [10]. It is an efficient method for spatio-temporal focusing of ultrasonic waves in complex medium [14]. The testing equipment consists of one dedicated transmitting and another receiving transducer (or alternatively a laser vibrometer). In the simplest form, the ultrasound focusing happens under the receiving transducer/vibrometer. The focusing process relies on the reciprocity of linear waves in undamaged medium. Since the wavefield spreads throughout the medium during the process, any local deviations from linearity can manifest in imperfectly focused signal, deviations from wave reciprocity (loss of symmetry) [20], difference in PI signals [17] or by other effects detectable with various NEWS methods [3,7]. Additional possibilities to attain a single spatio-temporal focusing include deconvolution method, related to TR [21], but is not covered here. For simplest reciprocal TR, a delta-pulse signal is transmitted. It is received during a time window T (TR window) to record the initial arrival of the signal and its various echoes. If this received signal is time reversed and re-sent without changing the roles of the transducers, then the delays of the echoes cancel out to produce a spatio-temporally focused signal at the receiver location. This can be further improved by using a coded signal as an initial excitation, which contains more wave energy, such as a linear chirp [22,23]. Then an additional cross-correlation step is needed. The steps for reciprocal TR are as follows [16]:
• A linear chirp-coded signal is transmitted from receiver c(t) = A · sin (ψ (t)), where ψ (t) is a linearly changing phase. The receiving transducer records a signal y(t) with time window T after propagation. In linear case, the received signal can be analysed as a convolution between the sent signal and the (unknown) impulse response (or Green’s function) h(t) of the medium y(t) = h(t) ∗ c(t). • Since the autocorrelation of a linear chirp is a delta signal c(t) ∗ c(t − T ) ∼ δ (t − T ), it is possible to reveal the effect of the impulse response by cross-correlating the sent and received signals Γ (t) = y(t) ∗ c(T − t) = h(t) ∗ c(t) ∗ c(T − t) ∼ h(t) ∗ δ (t − T ). • Propagating the time reversed Γ (T − t) through the material in the same direction and configuration will result in a spatio-temporally focused wave yTR = Γ (T − t) ∗ h(t) ∼ δ (t) under the receiving transducer. 2.2. Delayed time reversal The delayed TR is a further development of reciprocal TR to optimise the received spatio-temporally focused wave by using the focused yTR (t) ∼ δ (t) as a new basis [16]. In a linear wave propagation any modification of the sent Γ (t − T )
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
3
will be reflected accordingly in the focused wave yTR (t), so it is possible to construct a new excitation
Γs (T − t) =
n ∑
ai Γ (T − t + τi ) =
i=0
n ∑
ai Γ (T − t + i∆τ ),
(1)
i=0
where ai is the ith amplitude coefficient and τi the ith time delay. In case of uniform time delay the ∆τ is the time delay between samples. In linear material this would give following signal at the focusing location: ydTR (t) =
[ ∑
] ai Γc (T − t + τi ) ∗ h(t)
linearity
HHHHHH
i
=
∑
ai Γc (T − t + τi ) ∗ h(t) =
∑
ai yTR (t − τi ).
(2)
i
i
In this paper this principle is used to excite the resonance frequency of a defect. Eq. (1) can alternatively be expressed by using a modulation function instead of periodic delayed sum.
Γs (T − t) = Γ (T − t) ∗ m(t),
(3)
where m(t) is some modulation function and will produce a focused wave with shape ydTR (t) = Γ (T − t) ∗ m(t) ∗ h(t) ∼ δ (t) ∗ m(t) = m(t).
(4)
In short, delayed TR is a method for creating a prescribed signal at the one-channel system receiver using deltapeak like focused wave of reciprocal TR as a basis. This is effective in case the medium has complex wave propagation characteristics (internal reflections), but can be hindered by attenuation. Instead of having a ‘‘nice’’ signal at transmitter and chaotic at the receiver, it calculates a Γ (T −t) or Γs (T −t) which will produce a signal of required shape at the receiving side [16] (or focusing location in case of this paper). Among other possibilities this allows either probing, scanning or detection of nonlinear effects. 2.3. Nonlinear elastic wave spectroscopy Time reversal methods are useful for creating a spatio-temporally focused ultrasonic waves, increasing signal to noise ratio and increasing signal amplitude by pulse compression among other uses. However, in order to detect the presence of a defect or create a spatio-temporal focusing in the area of the defect, the TR needs to be coupled with NEWS. In context of this work, it would be called delayed TR–NEWS. 2.3.1. Pulse inversion One of the simplest nonlinearity analysis methods is Pulse Inversion [17]. It is similar to widely-known Scaling Subtraction Method [24]: the sign of the amplitude is modified instead its magnitude, to take advantage of excitation symmetry. Further development of such methods have been researched [25], but this paper considers simple PI for wave field focusing onto the defect area. For PI, in case of linear media, positive and negative chirps can be sent through the medium, yielding positive and negative received signals ypos (t) = h(t) ∗ c(t),
yneg (t) = h(t) ∗ (−c(t)) ,
(5)
which in case of linear, undamaged medium have identical results in opposite signs yd = ypos + yneg = h(t) ∗ c(t) + h(t) ∗ (−c(t)) = h(t) (c(t) − c(t)) = 0.
(6)
The above relation cannot be applied for nonlinear processes: the convolution by impulse response assumes linearity to be valid. There are many numerical and physical experiments conducted with PI showing its ability to detect defects by the difference yd ̸ = 0 [4,9,23]. Ultrasonic waves can be focused at the defect area instead of the receiver by first extracting the nonlinear signature of the defect using PI, and then focusing on this signature using reciprocal TR. The nonzero difference yd is computed (Eq. (6)) and correlated with c(t) to yield the nonlinear signature response of the medium (deviation from the linear case): yd (t) ∗ c(T − t) = ΓNL (t). This is time-reversed and re-propagated to create signal focusing at the defect region [9]. Again this global focusing of the wave on the defect area does not require any experimental tuning, but emerges from the signal processing steps. 2.3.2. Resonance excitation and defect spectroscopy The resonance of a CAN defect could be exited more during illumination by ΓNL (T − t) (which excites mainly the defect area) than by the simple chirp c(t) (which excites mostly area near the transmitter): this can be seen from the maximum recorded u2 displacement amplitude plots in Fig. 2. The received output signal spectra can be compared (chirp against PI) to reveal which frequencies are amplified by the increased focusing of the ultrasonic wavefield near defect. We can assume these frequencies to be linked to the resonant frequency of the defect. Then it is possible to use delayed TR (Eqs. (1) and (2)) to modulate the focused waves at the defect area according to the found resonance frequency, to further excite the damage.
4
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
Fig. 1. Simulation configuration and TR-NEWS focusing of waves on the defect which is far from transducers. Pseudocolor plot shows displacement (in m) in y direction. Spatio-temporal focusing results from transmission of the ΓNL signal. Drawing shows the configuration of the test, including shear wave transmitting transducer, single receiving signal point and location of 10 mm defect in the material.
2.3.3. Defect resonance excitation Finally the steps of the method can be summarised. 1. Illuminate the sample with chirps (positive c(t) & negative −c(t)). Extract the nonlinearity signature using PI yd = ypos + yneg . 2. Cross-correlate the difference with chirp, yielding nonlinear impulse response yd (t) ∗ c(T − t) = ΓNL (t). 3. Create PI focused ultrasound on the defect by re-broadcasting the time-reversed nonlinear impulse response ΓNL (T − t). 4. Conduct spectroscopy of the output signals: find the enhanced frequency by comparing the output signals from chirp excitation and PI focused excitation. This enhanced frequency can be the resonance frequency of the defect. 5. Create resonance excitation by delayed TR–NEWS procedure. Apply the found resonance time period as the delay value in Eq. (1). 6. By transmitting the Γs (T − t) signal, the defect should be excited to larger displacement amplitudes, using the same level of output power. The validity and usability of this procedure will be investigated in this paper using FEM simulations. The advantages of this method are that it is only single-channel method, capable of using laser vibrometer and it does not require explicit focusing or tuning for focusing process. It is a global method meaning little to no scanning is required for the defect detection, but the disadvantage is that the spectrum needs to be analysed in step 4, which requires human interaction. 2.4. Simulation model It is known that the analytical methods have limited use in fully examining the nonlinear wave propagation in complex media, therefore this study is conducted using numerical experiments. The formulae used in describing the delayed reciprocal TR–NEWS utilise convolution by impulse response as a model for wave propagation in linear material. In case of nonlinearities, these concepts are not applicable (with the loss of linearity and reciprocity of wave propagation being needed for detecting the nonlinear effects). The computational model (Fig. 1) is a CFRP laminate, plate of 150 mm length by 3.24 mm thickness. It consists of carbon fibre layers in alternating 90◦ and 45◦ directions with pure epoxy layers between (total of 22 pairs). Boundary conditions on the left and right side are absorbing to simulate a larger plate. The analysed defect lengths are 10 mm, 5 mm or 1 mm. The defect is contacting with initially zero gap and force, with its centre always located at 100 mm from left side of the simulation region and 1 mm from the top surface. The rest of the simulation and contact model is an extension of previous work [26,27] except that in this work the receiving transducer is not near the defect. The transmitting and receiving transducers are fixed and are never interchanged, since this work uses reciprocal TR. Transmitter has allowed bandwidth from 0 to 2 MHz for all signal transmission steps: chirp c, time-reversed correlation Γ or its modified versions ΓNL or Γs . The initial linear chirp excitation is between 0.2 to 2 MHz for 60 µs. The time window of reciprocal TR is 120 µs for 1 and 5 mm defect cases and 180 µs for 10 mm defect case. Maximum input signal stress is limited to 50 kPa for all signals and is angled at 70◦ from normal, as indicated in Fig. 1. The signal is sampled at 2 Gs/s. The time window length is quite close to minimal for this model due to the high computational cost, but sufficient. Due to the infinite plate NDT problem, the sides have absorbing boundary conditions, so they reflect limited amount of wave energy back for refocusing, making the model more realistic but pessimistic. The crack motion is governed by Lagrange plus penalty method in normal direction and penalty method in tangential direction with coefficient of friction µ = 0.6. Other parameters and procedures for this model can be found in [27]. 3. Analysis and discussion of simulation results Assuming that the PI focusing excites the defect area comparatively more than chirp input signal does, we can compare and analyse the power spectra of the respective output signals (in this case the u2 displacement component at the
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
5
Fig. 2. Comparison of chirp (c(t)) and Pulse Inversion excitation (ΓNL (T − t)) transmissions via maximal generated u2 displacement through all time. In case of a good focusing on defect area (PI), the maximal displacement amplitudes take place near the defect area. In case of simple chirp it is near the source. Here PI is effective for focusing in case of 10 and 5 mm defects and not effective for 1 mm defect.
receiving transducer) to try to detect the frequency of the crack ‘‘clapping’’. Fig. 2 shows that for 10 mm and 5 mm defects there is a large enough difference between PI and chirp in maximum excitation amplitude near the defect area to possibly detect the change in spectrum, while for the 1 mm defect there might not be enough difference to detect the nonlinear signature. The symmetry of the maximum displacement values are due to the symmetry of the medium and the single-channel test setup.
6
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
Fig. 3. Spectral densities of the output signals: comparison between chirp and PI excitations, where PI should magnify the spectrum of nonlinear contact due to its better focusing of ultrasonic power on the defect location.
3.1. Spectroscopy of the output signal In NDT applications only the input and output signals are known. Numerical simulation has the added benefit of revealing the wave motion inside the material and in the damaged region. This helps in the development and verification of the procedures. Spectra of both the output signal and crack surface motion are analysed and compared. We want to know if the output signal contains sufficient spectral information to reveal the defect resonance and we can use the additional information from spectrum of the crack (which is available in simulation only) to verify it. The spectral densities are used for comparison of the spectra [28]:
{ S(k, t) =
2|U(k,t)|2 N |U(k,t)|2 N
,
,
k = 1, 2, 3, . . . , N /2 − 1, k = N /2,
(7)
where |U(k, t)| is the magnitude of kth wavenumber of FFT of a variable (here displacement u2 at the receiver). This is suitable for analysis of the output signal. Alternative methods could be used for comparing the spectra of the signals, like spectral amplitudes, periodograms, and other methods [29]. The spectral analysis of output signals are shown in Fig. 3. For 10 mm defect, the 42 kHz frequency is amplified, for 5 mm it is 132 kHz and for 1 mm defect it is difficult to judge but we choose 425 kHz. 3.2. Spectroscopy of the defect motion s The normal gap of the slave and master surface node coordinates are calculated gN = xm 2 − x2 for spectral analysis. The spectrum of the contact gap can be visualised in 2D after normalisation:
Snorm (k, t) =
S(k, t) Ssum (t)
,
Ssum (t) =
N /2 ∑
S(k, t).
(8)
k=1
This is shown in Fig. 4 for the three defect sizes. The ‘‘clapping’’ of the defect under PI focused wave loading is evident as vertical lines. From this, it is straightforward to calculate the average temporal period of the crack motion for verifying the found resonance frequencies (since the data for the crack is only available in simulations, not NDT measurements). For 1 mm sized defect the period of its motion is difficult to conclude, which can be due to the lack of PI wave focusing on the defect area or just the small size of the defect compared to the maximum 2 MHz input signal frequency in the simulation.
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
7
Fig. 4. Normalised spectrum of crack motion for 10 mm, 5 mm and 1 mm defects showing the periodical ‘‘clapping’’ of the defect under PI focused excitation.
An alternative method to analyse the motion of the defect gap is the cumulative spectral density, which considers the summed contributions of wavenumbers from mth number to N /2, cutting out the low frequency contribution for better analysis of nonlinear effects [28]. Sc (m, t) =
N /2 ∑ k=m
Snorm (k, t),
m = 1, 2, 3, . . . , N /2.
(9)
8
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
Fig. 5. Cumulative spectrum of 10 mm crack motion for various values of m.
The value m must be chosen so that it cuts out the unneeded low-frequency content, in order to better indicate the ‘‘clapping’’ motion of the defect. Fig. 5 shows an example of the cumulative spectrum of the contact gap motion in case of the 10 mm defect. It is straightforward to select a suitable m value and then find the frequency of clapping by a second discrete Fourier transform of the cumulative spectrum (which is not shown in this work).
3.3. Crack resonance excitation results The output signal spectrum (Fig. 3) is used to find the possible resonant frequency of the defect to be excited using delayed TR–NEWS procedure. For 10 mm defect, the 42 kHz frequency is amplified, for 5 mm it is 132 kHz and for 1 mm defect 425 kHz. The frequencies are verified using the aforementioned spectral analysis of the defect gap motion (discrete Fourier transform of cumulative spectrum yielding respectively 50 kHz, 138 kHz and 423 kHz for resonance frequencies). The supposed period of the defect resonance was input to delayed TR–NEWS formula (Eq. (1)) which completes the signal processing procedure for crack resonance excitation listed in Section 2.3.3. Further optimisation could be possible in regard of the resonance frequency detection. The procedure is moderately successful in resonance excitation. Its effect on the amplitudes of the captured signals (at the receiving transducer) can be seen in Fig. 6. As the input amplitude is always limited to 50 kPa, any increase in measured amplitudes at focusing location (when comparing single PI focusing with delayed TR resonance excitation) would indicate increase of power in the system thanks to the resonance of the nonlinear defect. In case of the 10 mm defect, the resonance excitation initially gives good results over PI, but then the resonance fades. Further optimisation and tuning might be needed for spectral analysis and longer TR window might be needed. In case of 5 mm defect in CFRP plate the results are good: The resonance is excited successfully and the output signal amplitude grows. In case of 1 mm defect there is no appreciable effect, which might be due to the difficulty of high-frequency signal reaching the defect through the layered and complex media. Because in the simulation it is possible to analyse also the defect motion, the results of resonance excitation can be verified by comparing the crack gap opening during excitation when ‘‘illuminating’’ the crack with PI focused wave versus with delayed TR–NEWS (in essence periodic PI illumination with the suspected crack resonance period). Fig. 7 shows this comparison: the master and slave node distances of a crack are summed, and delayed TR–NEWS result is subtracted from PI result (yielding the measure |∆y|). It can be seen that in case of 10 mm defect the delayed TR–NEWS resonance excitation initially yields larger crack opening motion than PI focusing excitation, but this fades soon. More optimisation might be needed in spectral analysis, since the simulation is very sensitive to the correct choice of supposed resonance frequency. In case of the 5 mm defect, there are enough delayed impulses to generate a good resonance in the simulation time period. Its larger resonance frequency is easier to determine from spectral analysis and the TR window is sufficiently long to produce a good resonance. The 1 mm defect case is problematic: due to insufficient PI focusing on the defect area, the correct resonance value could not be determined nor excited. This is most probably due to the difficulty of high frequency content of the chirp reaching the defect area, caused by the layered model.
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
9
Fig. 6. Output signals: u2 amplitude from single PI focused excitation versus periodic delayed TR–NEWS excitation.
Fig. 7. Crack gap opening comparison from single PI focus excitation or resonance from periodic delayed TR–NEWS excitation. The green colour indicates times where delayed TR–NEWS resonance is creating a larger crack opening than single PI focusing, red colour means vice versa. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
4. Conclusion Simulations in complex CFRP media, using a single-channel ultrasonic testing setup limited to 2 MHz, show that it is possible to detect the nonlinearity signature globally of a CAN type defect using the TR–NEWS procedure. Reciprocal TR procedure with NEWS allows to use the nonlinearity signature to focus the ultrasonic wave energy onto a small defect area which is far from either the transmitting or the receiving transducer. Comparing the received ultrasonic signals from focused versus unfocused signals, it is possible to measure the resonance frequency. Using the found frequency of the defect, the delayed TR–NEWS can be leveraged to periodically excite the damaged region with the focused ultrasound. Resonance excitation results show potential for enhancing the wave motion more than a simple PI focused wave. However, care must be taken that firstly the defect can be resolved by PI focusing and secondly that a large enough TR time window is used in case of larger defects that require longer excitations for resonance. These results need to be verified in physical experiments for materials of various attenuation levels, since it can impact the focusing quality of the delayed TR–NEWS method. 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.
10
M. Lints, A. Salupere and S. Dos Santos / Wave Motion 94 (2020) 102526
Acknowledgments This work was supported by institutional research funding IUT IUT33-24 of the Estonian Ministry of Education and Research. References [1] K.R. McCall, R.A. Guyer, Equation of state and wave propagation in hysteretic nonlinear elastic materials, J. Geophys. Res. 99 (1994) 23887–23897. [2] D. Meegan, P.A. Johnson, R.A. Guyer, K.R. McCall, Observations of nonlinear elastic wave behavior in sandstone, J. Acoust. Soc. Am. 94 (6) (1993) 3387–3391. [3] K.E.-A. Van Den Abeele, P.A. Johnson, A. Sutin, Nonlinear elastic wave spectroscopy (NEWS) techniques to discern material damage, Part I: Nonlinear wave modulation spectroscopy (NWMS), Res. Nondestruct. Eval. 12 (2000) 17–30. [4] I.Y. Solodov, N. Krohn, G. Busse, Nonlinear ultrasonic NDT for early defect recognition and imaging, in: Proceedings of European conf. on NDT (ECNDT), 2010, 734–758. [5] B. Sarens, B. Verstraeten, C. Glorieux, G. Kalogiannakis, D. Van Hemelrijck, Investigation of contact acoustic nonlinearity in delaminations by shearographic imaging, laser doppler vibrometric scanning and finite difference modeling, IEEE Trans. Ultrason. Ferroelectr. Freq. Control 57 (6) (2010) 1383–1395, http://dx.doi.org/10.1109/TUFFC.2010.1557. [6] S. Delrue, K. Van Den Abeele, Three-dimensional finite element simulation of closed delaminations in composite materials, Ultrasonics 52 (2) (2012) 315–324, http://dx.doi.org/10.1016/j.ultras.2011.09.001. [7] K.E.-A. Van Den Abeele, J. Carmeliet, J.A. TenCate, P.A. Johnson, Nonlinear elastic wave spectroscopy (NEWS) techniques to discern material damage, Part II: Single-mode nonlinear resonance acoustic spectroscopy, Res. Nondestruct. Eval. 12 (2000) 31–42. [8] A. Klepka, L. Pieczonka, W.J. Staszewski, F. Aymerich, Impact damage detection in laminated composites by non-linear vibro-acoustic wave modulations, Composites B 65 (2014) 99–108, http://dx.doi.org/10.1016/j.compositesb.2013.11.003. [9] T. Goursolle, S. Dos Santos, O. Bou Matar, S. Callé, Non-linear based time reversal acoustic applied to crack detection: Simulations and experiments, Int. J. Non-Linear. Mech. 43 (3) (2008) 170–177, 11th International Workshop on Nonlinear Elasticity in Materials. http: //dx.doi.org/10.1016/j.ijnonlinmec.2007.12.008. [10] T.J. Ulrich, A.M. Sutin, R.A. Guyer, P.A. Johnson, Time reversal and non-linear elastic wave spectroscopy (TR NEWS) techniques, Int. J. Non-Linear Mech. 43 (3) (2008) 209–216. [11] J. Kober, Z. Prevorovsky, Theoretical investigation of nonlinear ultrasonic wave modulation spectroscopy at crack interface, NDT & E Int. 61 (2014) 10–15, http://dx.doi.org/10.1016/j.ndteint.2013.09.001. [12] K. Hajek, V. Nenakhova, Principles of a defect localisation in nonlinear ultrasonic mixing impulse spectroscopy, Adv. Mil. Technol. 9 (2) (2014) 41–48. [13] F. Ciampa, M. Meo, Nonlinear elastic imaging using reciprocal time reversal and third order symmetry analysis, J. Acoust. Soc. Am. 131 (6) (2012) 4316–4323. [14] M. Griffa, B.E. Anderson, R.A. Guyer, T.J. Ulrich, P.A. Johnson, Investigation of the robustness of time reversal acoustics in solid media through the reconstruction of temporally symmetric sources, J. Phys. D: Appl. Phys. 41 (8) (2008) 085415. [15] D. Givoli, Time reversal as a computational tool in acoustics and elastodynamics, J. Comput. Acoust. 22 (03) (2014) 1430001. [16] M. Lints, S. Dos Santos, A. Salupere, Solitary waves for non-destructive testing applications: Delayed nonlinear time reversal signal processing optimization, Wave Motion 71 (2017) 101–112, http://dx.doi.org/10.1016/j.wavemoti.2016.07.001. [17] C.-C. Shen, Y.-H. Chou, P.-C. Li, Pulse inversion techniques in ultrasonic nonlinear imaging, J. Med. Ultrasound 13 (1) (2005) 3–17, http: //dx.doi.org/10.1016/S0929-6441(09)60073-4. [18] M. Fink, Time reversal of ultrasonic fields - Part I: Basic principles, IEEE Trans. Ultrason. Ferroelectr. Freq. Control 39 (5) (1992) 555–566. [19] C. Draeger, M. Fink, One-channel time reversal of elastic waves in a Chaotic 2D-silicon Cavity, Phys. Rev. Lett. 79 (1997) 407–410. [20] T. Ulrich, M. Griffa, B. Anderson, Symmetry-based imaging condition in time reversed acoustics, J. Appl. Phys. 104 (6) (2008) 064912. [21] T. Ulrich, B. Anderson, P.-Y. Le Bas, C. Payan, J. Douma, R. Snieder, Improving time reversal focusing through deconvolution: 20 questions, Proc. Mtgs. Acoust. 16 (2012) http://dx.doi.org/10.1121/1.4764487. [22] E.C. Farnett, G.H. Stevens, Pulse compression radar, Radar Handb. 2 (1990) 10–11. [23] S. Dos Santos, Z. Převorovský, Imaging of human tooth using ultrasound based chirp-coded nonlinear time reversal acoustics, Ultrasonics 51 (6) (2011) 667–674. [24] C. Bruno, A. Gliozzi, M. Scalerandi, P. Antonaci, Analysis of elastic nonlinearity using the scaling subtraction method, Phys. Rev. B 79 (6) (2009) 064108. [25] S. Dos Santos, C. Plag, Excitation symmetry analysis method (ESAM) for calculation of higher order non-linearities, Int. J. Non-Linear. Mech. 43 (2008) 164–169. [26] M. Lints, A. Salupere, S. Dos Santos, Simulation of detecting contact nonlinearity in carbon fibre polymer using ultrasonic nonlinear delayed time reversal, Acta Acust. United Ac. 103 (6) (2017) 978–986. [27] M. Lints, A. Salupere, S. Dos Santos, Simulation of Defects in CFRP and Delayed TR-NEWS Analysis, Research Report Mech 320/17, Tallinn University of Technology, Department of Cybernetics, 2017. [28] A. Salupere, The pseudospectral method and discrete spectral analysis, in: E. Quak, T. Soomere (Eds.), in: Applied Wave Mathematics, Springer Berlin Heidelberg, Berlin, Heidelberg, 2009, pp. 301–333. [29] P. Stoica, R. Moses, Spectral Analysis of Signals, Prentice Hall, 2005.