A multiband approach for accurate numerical simulation of frequency dependent ultrasonic wave propagation in the time domain and nearfield acoustical holography in wedge propagation

Finite element (FE) simulations are popular for studying propagation and scattering of ultrasonic waves in nondestructive evaluation. For a large number of degrees of freedom, time domain FE simulations are much more efﬁcient than the equivalent frequency domain solution. However, unlike frequency domain simulations, time domain simulations are often poor at representing the speed and the attenuation of waves if the material is strongly damping or highly dispersive. Here, the authors demonstrate efﬁcient and accurate representation of propagated and scattered waves, achieved by combining a set of time domain solutions that are obtained for a set of frequency ranges known as bands, such that, in combination, the authors’ multiband solution accurately repre-sents the whole wave spectrum. Consequently, high accuracy is achieved, at minor computational cost, using a modest number of bands. The multiband technique is implemented for ultrasonic wave propagation in highly attenuating polyethylene material, using three frequency bands, and can yield a reduction in empirical acoustic properties fractional error compared with respective time domain simulations, in propagation duration, of a factor of 1.4, and in full-width-half-maximum, of a factor of 10. Last, the accuracy of this approach is further exempliﬁed in a wave scattering simulation.


I. INTRODUCTION
Nondestructive evaluation (NDE) is widely used in industry to detect potential defects in engineering components such that the likelihood of critical damage or failure is minimised. The development of accurate and reliable NDE methods increasingly uses numerical model simulations. The work presented here arose from a need to perform simulations of ultrasound NDE of high-density polyethylene (HDPE) pipes. The simulations were performed using finite element (FE) analysis, but this was not straightforward because HDPE is a highly attenuating material. The methodology presented here was the result of developments to enable accurate FE simulations for this problem, but it is presented here as a general methodology as it has a wide value for applications with damping materials. FE simulations can be performed in either time or frequency domains. It is usual that time domain calculations are performed explicitly, such that current solutions are obtained directly from previous solutions, and under certain preferred conditions-central difference time marching, with diagonal mass representation 1 -this can be done by spatially local calculations without ever needing to assemble or invert a full system matrix. With no loss of accuracy or generality, the efficiency is further improved through the use of parallel processing. A particularly effective approach is to use graphics processors, such as is done by the Pogo FEA (Ref. 2) software. High efficiency is of particular significance when modelling wave propagation in large regions where a fine FE mesh is required, necessitating field solutions for a large quantity of elements; multiple millions can be used in many current applications.
The principle FE alternative, frequency domain finite element (FDFE) simulations, are solved implicitly and for a large number of harmonic solutions. This is highly computationally expensive, particularly for large modelled geometries, requiring the assembly and inversion of system matrices. 1 Further to this, some large models cannot be solved practically at all using FDFE because the system matrix is too large for the computer memory, and inversion is not possible within a practical time scale, even though the matrix is sparse.
The problem for the modeller arises because explicit time domain finite element (TDFE) solutions are often poor at accurately representing the variation with frequency of acoustic properties-specifically the dispersion relation, which comprises attenuation and phase velocity-and therefore for many such applications, neither TDFE or FDFE are suitable.
Consequently, we developed the multiband finite element (MBFE) approach to simulation, with which a full, a) Electronic mail: jack.egerton09@imperial.ac.uk accurate dispersion relation can be modelled by blending the solutions from multiple separate TDFE solutions that are obtained at different frequency bands. MBFE couples both the accuracy and efficiency of FE analysis of waves with dispersion relations that can exhibit frequency dependence of general form. The multiband approach is developed to recover the lost accuracy of explicit TDFE without incurring the heavy computational penalty of FDFE, or the limitations in model size of the latter.
Our example application material, HDPE, is viscoelastic. The properties of materials such as this can be described by extrapolation and inference-often from multiple decades below ultrasonic frequencies-of frequency-dependent material property parameters, such as storage and loss moduli, 3,4 provided that the modelling constraints described above do not apply. We also achieve an accurate description of ultrasonic wave propagation and scattering in the MBFE methodology without need for extrapolation across many frequency decades and the consequent inference of material properties, and without such modelling constraints as model size, found in similar FDFE models, using acoustic properties-attenuation and phase velocity.
In our preceding study of the ultrasonic viscoelastic acoustic properties of HDPE pipe material, 5 we provide a unique, accurate, and reliable description of the viscoelastic properties of HDPE, covering the range of frequencies and temperatures needed for practical ultrasound NDE of this material. The properties are presented there in the desired manner as spectra of attenuation and phase velocity. The variation in properties with frequency is the modelling challenge which we aim to address here; while this will take the specific example of HDPE properties, the approach is general and can be used for any dispersive medium.
The multiband technique may be used with dispersion relations of general frequency dependence, because it places no constraint on the frequency dependence of the acoustic properties-attenuation and phase velocity-used to define the dispersion relation. The multiband technique is here applied to the example of viscoelastic HDPE, in which bulk longitudinal waves propagate. More generally, such application of FE analysis-to the accurate description of viscoelastic or highly attenuating wave propagation, and to other waves for which attenuation and phase velocity have frequency dependence-is important to many areas of research, including NDE, 6 medical imaging, 7 and seismology. 8 This article is set out as follows. First, in Sec. II, is the background theory required to formulate MBFE; next, in Sec. III, is the outline of the MBFE method; in Sec. IV are the results and analysis of MBFE applied to wave propagation and scattering examples in HDPE pipe material; and in Sec. V conclusions are drawn.

II. MATERIAL MODEL
The theory used to define the ultrasonic acoustic properties of viscoelastic HDPE pipes is provided here to form the basis for exemplification of the MBFE model.
The spectrum of a one-dimensional harmonic plane wave can be expressed as where x is location, f is frequency, and k x ðf Þ, or tersely k(f), is the dispersion relation where v p ðf Þ is phase velocity and aðf Þ is attenuation. For this entire study, the following units are used, ½v p ¼ ms À1 and ½a ¼ m À1 . In our preceding HDPE pipe acoustic properties study 5 we presented a description of ultrasonic viscoelastic bulk wave dispersion that conforms to theoretical viscoelastic constraints and that is necessarily causal, using the Kramers-Kronig relation. The variation in properties with frequency is the modelling challenge which we aim to address here; while this will take the specific example of HDPE properties, the approach is general and can be used for any dispersive medium. Consequently, attenuation and phase velocity vary with frequency using and where 0 < a 1 ; 1 < y < 2, and positive phase velocity at zero frequency, 0 < v p1 , are all known empirical coefficients. It is shown in our prior analysis 5 that from Eq. (1), general forms of attenuation and phase velocity are and v p;gen f ð Þ ¼ 2pfx where / S 0 and / S are phases of the initial and propagated waves.

III. METHOD
The method for obtaining the MBFE model is described here, followed by its validation using wave propagation simulations that well match those of the empirical acoustic properties obtained above. The simulation technique is also applied to a scattering example.

A. MBFE
The underlying concept of the MBFE methodology is to solve the simulation multiple times in the time domain, in each case covering only a subset of the frequency bandwidth, and then assemble the multiple results to obtain the outcome for the broader bandwidth problem. Each of the narrowband solutions needs to be narrow enough to achieve acceptable accuracy, given the frequency-varying properties, and at the same time wide enough to avoid the need for large numbers of contributing bands.
The wave characteristics used in the described implementation of the MBFE model are chosen to be suited to real inspection requirements. Piezoelectric ultrasonic transducers are designed with a Q-factor-the degree of damping of an object oscillating at resonance-that results in a single pulse duration that is generally on the order of 3-4 cycles. These short duration pulses have broad frequency bandwidths, and consequently exhibit greater variation in frequencydependent dispersion within their non-negligible frequency range. FE simulations in this study are therefore conducted with cycles, n cyc ¼ 4, such that meaningful comparison can be made between accurately represented frequencydependent dispersion, using the MBFE approach, and otherwise using TDFE.
A commonly used pulse window is a half-period sinusoid given by a Hann function. The amplitude spectrum of this function has large peaks, known as side lobes, either side of the central peak at pulse centre frequency, f 0 . This is unsuitable when general frequency dependent dispersion is assumed, as it is necessary to describe viscoelastic media such as HDPE, because during wave propagation the frequency content of the main lobe can be attenuated more than the unphysical side lobes, resulting in dominance of these erroneous frequencies in propagated pulses. This is avoided through selection of an envelope that generates an amplitude spectrum with smaller side lobes, such as the Blackmann-Harris window, 9 P n ð Þ ¼ a 0 À a 1 cos 2pn N À 1 þ a 2 cos 4pn N À 1 À a 3 cos 6pn N À 1 ; for samples n ¼ 1; ::; N, where a 0 ¼ 0:35875; a 1 ¼ 0:48829; a 2 ¼ 0:14128; and a 3 ¼ 0:01168. This is the chosen window for all FE simulations in this article. Using an input pulse with accurate frequency content relaxes the filtering requirements in post-processing of the simulated wave; if large side lobes are present the high-and low-pass filters used necessarily have corner frequencies close to f 0 , resulting in greater undesired filtering of the frequency content of the main spectrum.
The low accuracy of the frequency dependent acoustic properties obtained using TDFE result from the necessary reduction of Eq. (2) to a fixed frequency solution. This fixed frequency representation is also used in the MBFE procedure, but accuracy is improved because each band is relatively narrow. The fixed frequency for this approximation is justifiably chosen to be that which the wave contains most of-the modal value of its amplitude spectrum-which is f 0 for approximately single-peak pulses, such as the tone bursts used here. The resulting simplified dispersion is with corresponding terse variables, k 0 , v p0 , and a 0 . Rayleigh damping is a damping model available in many commercial explicit TDFE software packages that is appropriate to the dispersion in Eq. (8), and that does not compromise the efficiency of the explicit time marching scheme, which is described using the equation of dynamic equilibrium 10 where ½M; ½C; ½K, and ½F are mass, damping, stiffness, and loading matrices, respectively. The following, u, _ u; and € u are particle displacement and its first and second time derivatives. The damping matrix can be expressed in terms of both the mass and stiffness matrices 10 where ½C; ½M, and ½K are the damping, mass, and stiffness coefficients and C M and C K are the Rayleigh coefficients, which can be defined in terms of attenuation for a fixed frequency using 11 and where the centre angular frequency is x 0 ¼ 2pf 0 . Assuming a harmonic displacement solution to Eq. (9) uðtÞ ¼ u 0 expðix 0 tÞ; where the angular frequency is x ¼ 2pf , and using Eq. (10), the equilibrium equation becomes implying the following transformations from an undamped to damped system for density and Young's modulus The low accuracy of TDFE results from the following. Describing the general frequency dependence of wave acoustic properties in HDPE using power law viscoelastic dispersion is highly accurate, as shown in Ref. 5. Conversely, it is inaccurate to reduce these acoustic properties to a fixed frequency description. Much of the accuracy of the power law description, or other dispersion relations, is recovered using MBFE. The favoured approach to conducting MBFE, known as "k-blending," is outlined here with relevant equations: (1) Divide the amplitude of the spectrum of the initial signal into a set number of bands, n band , containing equal frequency content (equal area under jS init j vs f curve). (2) Perform separate TDFE simulations with centre frequencies equal to the frequency at the centre of each band. (3) Window and filter all waveform outputs such that distinct initial and latent signals are isolated. (4) Calculate attenuation and phase velocity from these separated signals. (5) Combine these acoustic properties using the dispersion relation of Eq. (2). (6) Blend together the dispersion relations of all bands using weighting functions, yielding k for the whole bandwidth. (7) Calculate an accurate description of the transfer function, H(f), that maps an initial pulse to a propagated pulse.
The window used in step (3) is a Tukey window, which is a rectangular window with Hann (sinusoidal) leading and trailing edges. Unlike for the above example of toneburst generation, where significant attenuation and corresponding dominance of amplitude spectrum side lobes results, windowing of waveforms is not subject to this complication, therefore Hann edges are no issue. A rectangular function is accurate because it does not alter pulse amplitude within the window range. However, it does sharply truncate the signal, which can alter the pulse frequency content. By contrast, use of Hann edges in Tukey windowing smooths truncation, minimising this effect. The corner frequencies of the filters in step (3) are distant from significant frequency content and have first-order roll-off.
The attenuation and phase velocity in step (4) are those of Eqs. (5) and (6). The blending described in step (6) involves linear weighting functions for the nth band, b n , that equal 1 at band centres where TDFE simulations are most accurate, equal 0.5 at adjacent band edges, and equal 0 at adjacent band centres. In the outer halves of the first and last bands all frequency content is used-equivalently b n ¼ 1 there. Using this blending, the whole bandwidth can be described with high use of accurate frequency content and low use of inaccurate content. Frequency content beyond the first and last band is alternatively linearly extrapolated for comparison of accuracy.
Below we describe and formulate general methods we have considered for the blending of attenuations and phase velocities, as well as blending of dispersion relations. Blending of propagated spectra, attenuations, and phase velocities, dispersion relations, and transfer functions are all inequivalent. For example, the blends of attenuations and phase velocities are not the same as the blends of dispersion relations, as follows: and v p;blend ¼ X while the dispersion blending used in MBFE is such that The mapping of initial pulse to propagated pulse described in step (7) includes both the amplitude and the phase of the spectra and is defined as follows: To ensure the k-blending MBFE technique is most suitable, other approaches to MBFE are also developed for comparison of accuracy, generality, and to a lesser degree, simplicity. These unfavoured approaches are listed here.
The first omitted approach takes the transfer function of each band evaluated at solely the band centre frequencies then assigns these fixed values, obtained from each TDFE simulation, to all frequencies within their corresponding band. The second approach is to linearly interpolate the transfer function between the band centre values. The next approach is logarithmic interpolation of the transfer function. Of blending approaches, first the propagated spectra of each band are blended, then the phase velocities and attenuations are blended, followed by dispersions, and finally the transfer functions. Of these, blending the dispersion relations (k-blending) proves accurate, general, and simple to implement and use.

B. FE meshing and geometries
The FE models described here represent the propagation and scattering of plane longitudinal waves in HDPE pipe material in two dimensions. In Fig. 1, the propagation FE simulations are in two dimensions; at the top and bottom model boundaries, displacement normal to these boundaries is set to zero; and HDPE is approximated as homogeneous and isotropic. A plane wave is approximated by a 100 mm linear source oriented 90 to the axial direction. The model is set within the FE software to have no displacements normal to the direction of wave propagation. The wave is monitored at two locations, one 10 mm from the source in the propagation direction, the next 100 mm beyond that.
The sum of all point sources along the source line approximates a plane wave, which is most suitable for HDPE where shear waves are attenuated on the order of 10 times that of the already highly attenuated longitudinal waves. 12 The simulated region is bounded by a large absorptive region that has decreasing stiffness towards the boundary, known as the stiffness reduction method (SRM). 13 This highly attenuates waves propagating in this region such that boundary reflections are minimal and the resulting wave is not complicated by undesired signals. The SRM was primarily designed for use with elastic materials. In such cases, the damping of the SRM decreases to zero towards the modelled region. We have adjusted the implementation of the SRM to decrease damping towards the viscoelastic damping of HDPE at the inner edges of the SRM. This ensure damping continuity at the boundary.
The optimal frequency for ultrasonic inspection of HDPE pipes is selected based on the following compromise. Higher frequencies can increase spatial resolution. Also, higher frequencies produce more directive beams. Last, attenuation in HDPE, and many other media, increases with frequency. Because of this, MBFE and TDFE simulations of both wave propagation and scattering are conducted with median propagation paths of x s ¼ 100 and 53.5 mm, respectively, all at frequency, f ¼ 2.25 MHz. Frequencies similar to this are often used in current industrial inspection of HDPE pipe joints.
The chosen FE mesh comprises structured triangular plane strain elements. The mesh is describable as fine 14 with element size dx k 0 20 ; where the wavelength at centre frequency, f 0 , is k 0 ¼ v p0 =f 0 and v p0 is the phase velocity at centre frequency. This mesh refinement, and our other chosen refinement, which is half the spatial step of this, are sufficiently fine for our purpose of demonstrating the improvements in accuracy of MBFE over TDFE simulations, for given mesh refinements. Given that for all simulations the time step is 2.5 times as refined as the spatial step, the time step contributes negligible computational error, while ensuring numerical stability. For increased FE simulation accuracy in both TDFE and MBFE, shown in Fig. 2, at the four square-oriented corners of two adjacent triangular source or monitor elements, four nodes are displaced parallel to lines that project to the centre of the squares that comprise two triangular elements. After summation of the displacements at the four nodes, where two nodes are common to both elements, this excitation approximates that of a purely longitudinal (pure-L) wave point source. 14 Following the study of the propagation and attenuation of plane waves, we have also demonstrated the MBFE methodology on a simulation of scattering from a discontinuity, in a study relevant for NDE. Shown in Fig. 3, the scattering FE simulations are two-dimensional and the scatterer is a cylindrical void with a radius of 3 mm. This is an approximation to a pertinent defect example, where a void may form around a liquid or solid manufacturing inclusion in the HDPE pipe joint. The FE mesh comprises triangular elements with the size stated in Eq. (22). Before and after scattering off the cylindrical void, the wave is detected at a pure-L monitor.

IV. RESULTS AND ANALYSIS
Provided here are the results and analysis of MBFE applied to our wave propagation and scattering examples. Cartoon of two adjacent square-oriented triangular source or monitor elements with four nodal displacements, depicted by arrows, which are directed parallel to lines that project to the centre of the squares that comprise the two triangular elements. Fig. 4 are 4-cycle tonebursts propagated 100 mm in the FE geometry depicted in Fig. 1. The solid curve is our benchmark propagated pulse, u, predicted using Eq. (2) with select acoustic properties obtained from the procedure described in our preceding HDPE pipe acoustic properties study. 5 Dot-dashed curve is the TDFE pulse, u TD , dx ¼ k 0 =40, and a time increment of dt ¼ 1= ð100f 0 Þ. Dot-dashed curve is the MBFE pulse with n band ¼ 3, u MB;3 .

Shown in
Qualitatively, it is clear that there is significantly improved matching to u for u MB;3 compared with u TD , specifically in both the amplitude and the phase of the pulse. Also the period is visibly much shorter for u TD than it should be, because the low frequencies are attenuated too much and the high frequencies too little with frequency independent attenuation.
In summary, using as few as three bands in MBFE results in a significantly improved representation of ultrasonic wave propagation relative to that of TDFE simulation.
The TDFE and MBFE attenuations are exemplified in Figs. 5 and 6 for n band ¼ 3 and n cyc ¼ 4, d x ¼ k 0 =40, and d t ¼ 1=ð100f 0 Þ. The attenuation of TDFE is approximately equal to the empirical attenuation at f 0 , resulting in attenuation being too high at frequencies below f 0 and too low at frequencies above f 0 .
The blending of attenuation in Fig. 6 shows close agreement with the desired empirical solution, obtained in our acoustic properties study, 5 over the frequency range covered by the bands. The bands are defined such that equal frequency content of the initial pulse exists in each band. However, the dispersion blending solution between bands cannot be an ideal linear approximation to the empirical dispersion because other sources of computational error exist in FE simulation. An example of this is the error associated with approximating the modelled geometry with elements of finite extent, where a coarser mesh decreases the FE accuracy; the accuracy increases as the mesh is refined. Similar results are expected for refinement of the time increment, dt.
In summary, using as few as three bands in MBFE results in a significantly improved representation of ultrasonic attenuation relative to that of TDFE simulation.
In Fig. 7, the phase velocity of TDFE is approximately equal to the empirical phase velocity at f 0 , resulting in velocity being too high at frequencies below f 0 and too low at frequencies above f 0 . For increasing frequency, the TDFE velocity error relative to the frequency independent ideal solution, v p;1 , increases because the number of elements per wavelength, k=40, decreases. In Fig. 8, as with MBFE attenuation, MBFE phase velocity shows close agreement to the empirical curve and improved accuracy at frequencies above and below f 0 compared with any of the three individual TDFE solutions that this blend is obtained from.
In summary, using as few as three bands in MBFE results in a significantly improved representation of ultrasonic phase velocity relative to that of TDFE simulation.
To further quantify the improved accuracy of MBFE over TDFE, for a given number of frequency bands, n band , mesh spacing, dx, and time increment, dt, the full-widthhalf-maxima (FWHM) and propagation durations (PDs) of the FE pulses are obtained as fractional differences from those of the empirical propagated pulse, u, for a 4-cycle toneburst, shown in Figs. 9 and 10, respectively. The FWHM corresponds to a 6 dB drop in amplitude on a logarithmic scale and the PD can infer wave propagation distance using knowledge of the wave velocity in the medium. The TDFE solutions are threshold lines for comparison with 3, 7, and 15 bands. The two alternative MBFE techniques, with and without extrapolation of the solution beyond the edge bands, are also compared.
In Fig. 9, the FWHM discrepancy is where mahðuÞ is the maximum amplitude of the Hilbert transform of the longitudinal wave displacement, and in Fig.  10, the PD discrepancy is  9. Full-width-half-maximum discrepancies between the empirical propagation of a 4-cycle toneburst and the TDFE and MBFE approximations to this. The dotted line is the TDFE fractional error with dx ¼ k 0 =20. For the same mesh refinement, the diamonds are the MBFE FWHM fractional errors and the crosses are the fractional errors for MBFE with extrapolation. Dotdashed line is the TDFE fractional error with dx ¼ k 0 =40. For the same mesh refinement, the squares are the MBFE fractional errors and the circles are the fractional errors for MBFE with extrapolation.
where t mahðu P0 Þ is the time that the maximum amplitude of the Hilbert envelope of the initial pulse, u P0 , occurs, and the other times correspond to the maxima for the other pulses. In Fig. 9, the TDFE solution with less mesh refinement has unexpectedly low fractional error. This is because, as is shown in Figs. 4 and 5 for dx ¼ k 0 =40, the pulse does not have the required frequency down shift, and therefore no pulse broadening in the time domain. However, an equal and opposite error occurs with the velocity, as also shown in Figs. 4 and 7, where high frequencies are significantly more dispersive than for the empirical velocity, resulting in erroneous pulse broadening.
These two errors counteract one another to produce a coincidentally accurate FWHM, where FWHM ¼ 0:018. By increasing the mesh refinement of the TDFE solution, and therefore decreasing numerical error, the fractional uncertainty in FWHM shown in Fig. 9 actually increases to FWHM ¼ 0:096. This is because the two numerical errors described above that have equal and opposite contributions to pulse dispersion in the time domain reduce at different rates with increased mesh refinement.
Because of these various contributing sources of error in FWHM, MBFE has high accuracy for three bands but there is no consistent trend towards increasing accuracy with increasing bands for FWHM. With three bands the FWHM errors for the coarse and fine mesh refinements of MBFE without extrapolation are FWHM ¼ 0:074 and 0.046. With extrapolations these errors reduce to FWHM ¼ 0:010 and 0.009, which are approximately a factor of 2 more accurate than the erroneously high accuracy of the coarse refined TDFE solution and approximately a factor of 10 more accurate than the TDFE solution with less computational error.
In summary, using as few as three bands in MBFE results in significantly improved representations of ultrasonic pulse FWHM relative to those of TDFE simulation, with the exception of one anomalously accurate TDFE result, explained above.
In Fig. 10, the coarse and fine mesh TDFE PD errors are t ¼ 0:047 and 0.038. The MBFE errors asymptote towards the highest accuracy achievable for a given mesh refinement. With three bands the PD errors for the coarse and fine mesh refinements of MBFE without extrapolation are t ¼ 0:042 and 0.035. With extrapolations these errors reduce to t ¼ 0:040 and 0.033, which are approximately a factor of 1.1 more accurate than the coarse refined TDFE solution and approximately a factor of 1.4 more accurate than the TDFE solution with less computational error.
In summary, using as few as three bands in MBFE results in significantly improved representations of ultrasonic pulse PD relative to those of TDFE simulation.
The absolute reduction in error from TDFE to MBFE is greater for FWHM than PD because FWHM is influenced by attenuation and phase velocity accuracy while PD error is only influenced by phase velocity; for viscoelastic media with dispersion relations that adhere to the Kramers-Kronig causality constraint, as with Eqs. (3) and (4), attenuation varies more than phase velocity with frequency and is therefore less well represented by the TDFE frequency independent solution. Extrapolation provides the largest reduction in MBFE FWHM error because of the high frequency dependence of attenuation, while mesh refinement provides the largest reduction in MBFE PD error because of the relatively large reduction in velocity accuracy when the number of elements per wavelength is lower.
Consequently, MBFE has been quantitatively validated using the example of ultrasonic wave propagation in HDPE that has a known general functional form for frequency dependent dispersion. MBFE has also been shown to provide significant improvements over TDFE for as few as three bands.

B. Wave scattering
Wave scattering is a second case where dispersion characteristics, including the influence of the scatterer, follow general frequency dependence describable by MBFE. The scattered pulse shapes of TDFE simulations will be compared with those of MBFE with three bands. The convergence of FIG. 10. PD discrepancies between the empirical propagation of a 4-cycle toneburst and the TDFE and MBFE approximations to this. The dotted line is the TDFE fractional error with dx ¼ k 0 =20. For the same mesh refinement, the diamonds are the MBFE PD fractional errors and the crosses are the fractional errors for MBFE with extrapolation. Dot-dashed line is the TDFE fractional error with dx ¼ k 0 =40. For the same mesh refinement, the squares are the MBFE fractional errors and the circles are the fractional errors for MBFE with extrapolation. the 3 band solution to a 15 band solution will be depicted for a high mesh refinement MBFE example with extrapolation at the outer band edges.
The following analysis is of simulations using the cylindrical void scattering geometry shown in Fig. 3, where the wave is monitored at the same location before incidence on the scatterer and after scattering. In Fig. 11, the solid curve is the 15 band MBFE scattered pulse, u MB;15 , the dashed curve is the TDFE scattered pulse, u TD , and the dotted-dashed curve is the band MBFE scattered pulse, u MB;3 . The TDFE pulse has an amplitude that is approximately double that of the 15 band MBFE pulse. Also, its period is too short, which will result in a FWHM that is too low. The phase is approximately onequarter of a cycle wrong. The PD is approximately equal to that of the 15 band MBFE pulse. The 3 band MBFE solution has fully converged to the 15 band MBFE pulse, within uncertainty bounds much smaller than the computational error caused by the finite mesh size.
In summary, using as few as three bands in MBFE results in a significantly improved representation of ultrasonic wave scattering relative to that of TDFE simulation.

V. DISCUSSION AND CONCLUSIONS
Conventional TDFE simulations provide insufficient accuracy when describing ultrasonic waves in materials where wave propagation results in a large reduction in wave amplitude, because they are limited to restrictive damping models, such as Rayleigh damping. FDFE is discussed as an alternative approach to such simulation, but this method is inherently highly computationally expensive for large model sizes and, using the usual implementation in many commercial software packages, FDFE can also be limited in its description of different frequency dependent dispersion relations. Accurate application of FDFE to, for example, viscoelastic waves often requires prior accurate knowledge of the frequency dependence of storage and loss moduli, rather than just the frequency dependence of attenuation and phase velocity. The MBFE approach is beneficial because such material properties as storage and loss moduli require significant materials testing to obtain. The inference of these acoustic properties using specific relations to material properties can result in increased uncertainty or error, potentially caused by such a mathematical description using unsuitable approximations to variations in, or unsuitable operating ranges in, the parameters with which the acoustic properties vary. For example, these material properties must usually be inferred using extrapolation from data obtained at frequencies significantly below ultrasonic frequencies. In wave simulation, use of accurate and reliable fully parameterised acoustic properties inherently bypasses these potential complications and sources of uncertainty or error. MBFE is consequently devised and proposed as an appropriate approach to modelling ultrasonic viscoelastic waves, and other such dispersion relations with general frequency dependence, only using prior knowledge of acoustic properties.
The MBFE approach is validated using a wave propagation example and shown to significantly increase accuracy over the TDFE example. For HDPE pipe material, while the large variation of its attenuation with frequency results in improved accuracy, even for a high number of MBFE bands, the smaller variation with frequency of its phase velocity leads to higher achievable absolute accuracy; using MBFE with three bands is highly beneficial for sizing accuracy, for example using a 6 dB drop in amplitude, and also beneficial for location accuracy because of its improved PD accuracy.
Scattering from a cylindrical void is demonstrated as an alternative case where MBFE is implemented to simulate waves with transfer functions that have general frequency dependence. The MBFE 3 band solution is fully converged to the 15 band pulse shape, within uncertainty bounds, while the TDFE solution features significant discrepancy.
MBFE is shown to improve accuracy over TDFE without incurring high computational cost. Furthermore, MBFE can be significantly more computationally efficient than FDFE for similar accuracy, or potentially more efficient with higher accuracy than an FDFE simulation implemented in many current commercial software packages. FIG. 11. FE simulated longitudinal wave pulses propagated in HDPE and scattered from a 3 mm cylindrical void. The input pulses are 4-cycle tonebursts with Blackman-Harris envelopes. The centre frequency of the input pulse is f 0 ¼ 2:25 MHz and the mesh spacing is dx ¼ k 0 =40 and the time increment is dt ¼ 3=ð400f 0 Þ. The pulse amplitudes are normalised to the maximum of the amplitude of the Hilbert transform of the input pulse, mahðu 0 Þ. The solid curve is the 15 band MBFE scattered pulse, u MB;15 , the dashed curve is the TDFE scattered pulse, u TD , and dot-dashed curve is the 3 band MBFE scattered pulse, u MB;3 .