Three-dimensional Ultrasonic Reverse-Time Migration Imaging of Submarine Pipeline Nondestructive Testing in Cylindrical Coordinates

: Submarine pipelines are a safe and energy-e ﬃ cient mode of gas transport. However, due to the complex manufacturing process and harsh operating environment, submarine pipelines are subject to fatigue cracks under long-term cyclic loading. A comprehensive and high-precision characterization strategy for submarine pipelines can e ﬀ ectively prevent potential safety hazards and have signi ﬁ cant economic and social repercussions. As a ma tt er of fact, pipeline defects cannot be reliably detected with current traditional 2D methods. On the other hand, in ultrasonic testing, cylindrical geometry increases the complexity of the 3D wave ﬁ eld in the submarine pipeline space and signi ﬁ cantly in ﬂ uences the accuracy of the detection results. In this paper, we put forward a novel method for 3D ultrasonic image testing that is suitable for cylindrical coordinates. In order to accurately simulate the ultrasonic signal received from pipelines, we generalize the 3D staggered-grid ﬁ nite-di ﬀ erence method from Cartesian coordinates to cylindrical ones and simulate the full wave ﬁ eld in the 3D pipeline space. Then, signal processing is performed on the ultrasound simulation records, and 3D reverse-time migration imaging of submarine pipeline defects can be e ﬀ ectively achieved using the reverse-time migration method and cross-correlation imaging conditions. The results obtained from simulations and real ﬁ eld data show that the proposed method provides high-quality 3D imaging of defects in pipelines, taking into account multiple sca tt ering and mode conversion information at the bo tt om of the defects.


Introduction
Submarine pipelines can connect subsea oil and gas resources with the entire onshore oil and gas production management system by the fastest, safest, and most economical route, which is called the "lifeline" of offshore oil and gas engineering [1]. However, due to the complex manufacturing processes and severe service environments [2], submarine pipelines are prone to crack voids, inclusions, and other defects, which greatly affects their mechanical properties and results in their premature failure. In real engineering applications, the initial imperfection is always introduced onto the pipes during the manufacture and installation procedures [3]. Submarine pipelines are subjected to fatigue loading in the harsh environment of the seabed operation conditions [4]. When the pipes are subjected to external pressure, failure first occurs in the cross-section with the most severe initial defects [5]. Additionally, once a crack develops, corrosion may promote its propagation, resulting in the cracking of the submarine pipeline [6]. When the pipeline is damaged and leakage occurs, it will not only cause economic losses, but also may bring about safety problems, endangering the safety of both people and the marine environment [7]. Therefore, before the submarine pipelines are put into service, a comprehensive high-precision inspection of them, including the rapid and accurate detection of defects, may effectively prevent safety hazards, which would guarantee a beneficial economic and social effect [8].
Several methods are available to inspect pipelines, including magnetic particle testing [9], magnetic flux leakage testing [10], acoustic wave testing [11], and penetration testing [12]. Ultrasonic testing may be used to detect the structure and shape of internal defects in submarine pipelines due to its high sensitivity, light equipment, and because it poses no harm to humans or the environment [13]. Ultrasonic testing can effectively detect the flaw size, the crack location [14], the elastic properties of materials [15], and the layup stacking sequence of composite materials [16,17] through the propagation of surface waves, guided waves, and body waves [18]. Then, by using ultrasonic imaging methods, defects in submarine pipelines can be visualized in an intuitive way [19][20][21]. The traditional ultrasonic methods of pipeline inspection include the synthetic aperture focusing technique (SAFT) [22][23][24], time-of-flight diffraction (TOFD) [25][26][27], and the total focusing method (TFM) [28,29]. Those methods have their own advantages and disadvantages. The SAFT is able to quickly provide media images; however, it has limitations in detecting vertical interfaces, lower boundaries of defects, and structures involving high-impedance contrasts. Additionally, there may be artifacts in SAFT-generated images due to surface waves, multiple reflections, and the mode conversions of the wave field originating at interfaces, resulting in wrong conclusions about the defect s location [30]. The TFM, on the other hand, only considers the direct ray path of ultrasound, without considering the mode conversions and multiple scattering arising from the interaction of ultrasonic waves with defects, which reduces the accuracy of the image [31,32]. Moreover, multiple wave reflections from the pipeline s lower side make the signals more complicated to analyze and lower the signal-to-noise ratio, often leading to artifacts in the reconstruction. As a consequence, it is difficult to image structures with vertical boundaries or complex geometry [33,34]. In order to circumvent such difficulties, modern ultrasonic testing techniques often involve full acoustic or elastic wave-equation modeling. In this paper, we address the 3D ultrasonic imaging of a pipeline using the reverse-time migration (RTM) technique, which has received much attention in the field of geophysics [35][36][37][38][39]. In recent years, RTM methods, originating from seismic imaging, have gained popularity in ultrasonic nondestructive testing applications [40][41][42][43]. The RTM method is a pre-stack imaging technique based on full-wave extrapolation. It does not suffer from the presence of oblique structures and offers a high resolution for the imaging of complex structures. In contrast to ray-based methods, RTM includes the effects of multiple scattering and mode conversions, as well as multiple wave reflections from the defect s lower side [33]. This allows one to gain more information, enabling the imaging of vertical interfaces and boundaries and providing higher quality images of interior defects [41,42]. However, most RTM studies on submarine pipelines are conducted using Cartesian coordinates [41,44], whereas, in any realistic case, the pipeline has an irregular cavity shape. In particular, for cylindrical structures, the pipeline s walls should be approximated as a staircase boundary using 3D Cartesian coordinates. Thus, it cannot accurately delineate the pipeline cavity, and the resulting grid scattering and dispersion analysis is affected by artifacts that reduce the imaging accuracy. To avoid those unwanted effects, we use cylindrical coordinates to model submarine pipelines in ultrasonic testing, which are more suitable for discretizing grids, and thus ensure higher accuracy [45]. In particular, using cylindrical coordinates makes the subdivision grids suitable to accurately represent the submarine pipeline cavity structure [46][47][48].
To address the shortcomings of traditional 2D methods and improve the accuracy of the ultrasonic characterization of submarine pipelines, we propose a 3D ultrasonic simulation and RTM imaging of submarine pipelines based on cylindrical coordinates. In our approach, the simulation of the 3D ultrasonic wave field and the wave field characteristic analysis of the submarine pipeline is performed by setting up a double free-surface and absorption boundary. Then, the ultrasonic RTM method is used to achieve high-quality 3D imaging of defects in pipelines. Numerical examples and real field data are used to prove the reliability and effectiveness of our method, paving the way for potential applications in practical ultrasonic testing.

Methodology
To perform accurate numerical simulations of ultrasonic waves in a 3D pipeline, we have developed a variable staggered-grid time-domain finite-difference numerical simulation method (FDM) in cylindrical coordinates. In the following paragraphs, we provide a brief description of the elastic wave equation, grid discretization, and boundary conditions in cylindrical coordinates.

Equations of Motion and Grid Discretization in Cylindrical Coordinates
For isotropic media, the first-order velocity-stress equation in the cylindrical coordinates (r, θ, z) can be expressed as follows [49]: , ,  (1) and (2) using a central differencing scheme with a staggered grid both in spatial and temporal domains [50]. As shown in Figure 1, within cell (i, j, k), the normal stresses on the half nodes (i + 1/2, j, k), (i, j + 1/2, k), and (i, j, k + 1/2). Then, according to the staggered grid with discretized points of velocity and stress components in Figure 1, Equations (1) and (2) are dissected in a finite-difference staggered-grid scheme. More details of the 3D time-domain FDM can be found in Liu et al. [49] for cylindrical coordinates. τ τ τ represent the shear stresses. The normal stresses are placed at the corner points around the staggeredgrid cell, the particle velocity components are located at the cell edges, and the shear stresses are sampled at the center of the cell faces.

Boundary Conditions
In finite-difference forward modeling of submarine pipelines, the boundary conditions are particularly important to properly simulate the propagation of ultrasonic waves. Figure 2 illustrates a schematic diagram of the implementation of boundary conditions. The absorbing boundary is placed along the R-and Z-direction of the model area, while the double free-surface boundary is along the inner and outer θ-direction of the pipeline. In our study, for the sake of computational efficiency, we introduce an improved vacuum formulation [51] into the cylindrical coordinates to set the double free-surface boundary condition in the θ-direction. At variance with the acoustic-elastic boundary method, which requires setting the free-surface boundary condition individually for each case, the improved vacuum formulation (IVF) is adaptable to an irregular free-surface. In order to eliminate interference from the region outside of the model, an absorbing boundary condition is assumed on the exterior of the model. We employ the split-field perfectly matched layer (S-PML) [52] to optimize absorption and minimize computational costs. As shown in Figure 2, the S-PML is placed in the R-direction and Z-direction of the pipe to absorb the body waves propagating toward the model boundary. The attenuation factor d(xi) of the absorbing boundaries is given by the following expression: where R is the theoretical reflection coefficient, L is the thickness of the absorbing boundary, and xi is the distance between grid points and model boundaries in the i direction.

Reverse-Time Migration Method
Reverse-time migration (RTM), originally proposed by Whitemore, is based on the two-way wave equation, on which reverse-time extrapolation is performed on the time axis [35]. RTM was first applied to the field of seismic imaging, providing higher imaging accuracy than previous methods. Besides accuracy, RTM has no inclination limitation and may be applied to arbitrary complex velocity models [36,37]. The RTM algorithm consists of the following three steps: 1. forward propagation of the source wave fields; 2. backward propagation of the receiving wave fields; and 3. imaging using imaging conditions. In this study, we used cross-correlation imaging conditions, followed by source normalization, ultimately leading to the following expression [53]: is the field at the receiver. After that, the Laplace filtering method is employed to eliminate the low-frequency artifacts caused by the imaging conditions [54].

Signal Processing
The wavefields obtained from pipeline ultrasound inspection records are characterized by complex features, such that the preprocessing of records is required before imaging to enhance RTM accuracy. In this study, wavelet extraction and dynamic balance in the channel are employed. Wavelet extraction methods may be classified into two categories: deterministic [55] and statistical [56]. In this study, the statistical wavelet extraction method is employed to extract wavelets from reflected waves. To this aim, we have to first select the reference trace, and then search for wavelets of other tracks within the travel time range of the reflected wave using this reference track as a guide. Finally, we normalize the wavelets of all traces, and then stack all traces in the record. Figure 3 illustrates the steps of the proposed method. The first step is to construct a geophysical 3D submarine pipeline model and implement the staggered-grid finite-difference method in cylindrical coordinates. Then, the SPML-absorbing and the IVF double free-surface boundary conditions are implemented in cylindrical coordinates and a 3D ultrasonic wave field simulation of the submarine pipeline is performed. RTM calculations represent the third and final step. RTM itself involves three steps, i.e., the forward propagation of source wave fields, the backward propagation of receiver wave fields, and the imaging step using the imaging condition. The procedure is carried out as follows: At first, a source wavelet is placed on the pipeline surface and used to excite the propagation of ultrasonic waves. The source wave fields during the forward propagation from T = 0 to T = max are calculated. Then, the recorded signal at the boundaries is time-reversed and simultaneously propagated back into the simulation domain to obtain the source wave fields from T = max to T = 0. Then, signal processing is performed and used as a signal at the receiver position. The receiver wave fields are then calculated using the FDM. Finally, the source wave fields are cross correlated with the receiver wave fields at each time step to construct the image. The final ultrasonic 3D RTM imaging is obtained using Laplace filtering [54].

Survey Layout
Our 3D cylindrical model of a submarine pipeline is shown in Figure 4. It includes a total of six survey lines arranged in the model region. Geophones are placed at the bottom, middle, and top of the pipeline s outer wall. The geophones are piezoelectric ceramic ultrasonic probes. Three survey lines are circumferentially placed (red lines in Figure 4), with a geophone spacing of 1.0 mm (Line-θ1, Line-θ2, and Line-θ3) and three are placed axially (blue lines in Figure 4), with a geophone spacing of 1.0 mm (Line-Z1, Line-Z2, and Line-Z3). The source is located at the center of the pipeline s outer wall.

Modeling
To analyze the propagation characteristics of the 3D wave field during the ultrasonic testing of pipelines, and to study the wave field characteristics resulting from different defects, we have designed three models (see Figure 5). Model-1 is a combined model of horizontal slag inclusion and hole defects, where the thickness of the slag is 2.0 mm, the width is 51.45 mm, the diameter of the hole is 2.4 mm, and the distance between the hole and the slag is 40 mm. Model-2 is a vertical slag inclusion model with a vertical thickness of 1.0 mm, a width of 51.45 mm, and a length of 60.0 mm. Model-3 describes another vertical slag inclusion with a slag width of 1.0 mm and a length of 60.0 mm. The model parameters are listed in Table 1.   Table 1.

Forward Modeling Results
In this study, we employ a spatial fourth-order and temporal second-order variable staggered-grid FDM. In our 3D cylindrical submarine pipeline model, the wall thickness is 45.0 mm, the inner diameter is 250.0 mm, the arc length is 102.97 mm (with a 20° rounding angle), and the axial length is 120.0 mm. The model size in the R-, θ-, and Z-directions is 45.0 mm × 102.97 mm (20°) × 120.0 mm. The radial step in the R-direction and the axial step in the Z-direction are Δr = Δz = 0.2 mm. The azimuthal step in the θ-direction is Δθ = 0.046°. The corresponding arc length increases with the increase in wall thickness (0.2 mm at the inner wall and 0.236 mm at the outer wall). The time step of Δt = 0.01 µs. To ensure the consistency with the actual ultrasonic source, we assume a point-like force source. As shown in Figure 6, the source wavelet is a ricker wavelet with a frequency of 0.7 Mhz, a signal width of 0.33-1. 16 Mhz, and a wavelet delay of 1.71 µs. In order to visually analyze the propagation of the wavefields in Model-1, we consider the 3D wavefield snapshots at T = 12 µs, and the 2D wavefield snapshot profiles from three slice directions at T = 6, 12, and 18 µs. Figure 7a shows the 3D wavefield snapshot, where the yellow color indicates the defect s location; Figure 7b shows a slice diagram, where green indicates a ROZ slice located at the midpoint position in the θ-direction, purple denotes a ROθ slice located at the midpoint in the Z-direction, and orange indicates a θOZ slice located at the pipeline s outer wall; and Figure 7c shows snapshots of the 2D wavefield in three slices.
In Figure 7, at 6 µs, we only see the direct wave P, because the wavefront has not yet encountered the defect. At 12 µs, since the Fresnel condition is not satisfied, the diffraction wave Php is generated when the P wave passes through the hole. The diffraction wave Pap1 is produced when the P wave passes through the upper interface of the slag inclusion. At 18 µs, we see the reflection wave Pp1, which is generated when the P wave reaches the free-surface boundary of the pipeline s outer wall. When the reflected wave Pp1 reaches the slag inclusion defect, it generates a reflected diffracted wave Ppap. Overall, the kinematics and dynamics of each wave item are consistent with the wave field law, which confirms the accuracy of the wave field simulations.

Model-1
Based on the observation system in Figure 4, we were able to generate the wave field records of the two sets of survey lines through forward modeling. The simulation record of Model-1 is illustrated in Figure 8. Trace refers to the number of geophones in the survey line. As can be seen from the plots, the direct wave P, the reflected wave Pp1, the multiple Pp2, the diffraction wave Php, and the diffraction wave Pap1 carry most of the energy. Nevertheless, the diffraction waves Ppap generated by the Pp1 are also evident. According to the simulations, the arrival time of the direct wave P is 9.

Model-2
The results from Model-2 are consistent with those from Model-1. In particular, we see in Figure 9 that increasing the propagation length of the P wave leads to a reduced energy of the wave field and a corresponding reduction in the energy received by the geophone on both sides of the pipeline. Nevertheless, both the primary reflection waves Pmp and secondary reflection waves Pmp2 resulting from the vertical slag inclusion defect and the reflected waves Pmm reflected from the inner wall of the pipeline have high energy, which can be precisely identified. In addition, the diffraction waves PmD1 and PmD2 from the vertical slag defect endpoints are also clearly visible.

Model-3
In addition, the results of Model-3 are consistent with those of Model-1. We see in Figure 10 that the direct wave P and the reflected wave Pap from the vertical slag inclusion are weak in Line-θ1 and Line-θ3, while the reflected wave Pp1 and the reflected wave Ppap from the inner wall of the pipeline through the vertical slag are stronger. Nevertheless, the direct and reflected waves are clearly visible along Line-θ2, with energy decreasing with distance.

Signal Processing Results
The signal processing of wavefield records is one of the key steps in RTM imaging, and it ultimately determines the overall quality of the imaging. Here, we present the results from the signal processing of the original recordings of Model-3. Looking at Figure  10a, we can see that there is a large amount of information in the unprocessed original record, which interferes with the defect reflection wave. We thus start with removing the direct wave P (see Figure 11) and then proceed with the wavelet extraction method to suppress the free-surface boundary reflection wave. Finally, the dynamic balance in the track is used to obtain the final record, which is characterized by an improved signal-tonoise ratio.  Figure 12 illustrates the results of RTM imaging for Model-2 after signal processing. As shown in Figure 12a, the unprocessed RTM image has a very low signal-to-noise ratio and contains many artifacts, which makes it impossible to identify the slag inclusion defect accurately. From Figure 12b, we can see that removing the direct wave allows us to eliminate the yellow-dashed-frame artifact on the outer wall of the pipeline. After the suppression of the free-surface boundary reflection wave, we obtain the image shown in Figure 12c, where the red-dotted-frame artifact has also been removed. Finally, using dynamic balance within the track, we improve the signal-to-noise ratio and obtain the final image of Figure 12d. When compared with the initial model in Figure 12a, it is evident that the position and shape of the vertical slag inclusion defect are essentially the same, and no other artifacts are present.  Figure 13 shows the three-dimensional ultrasonic RTM images obtained in cylindrical coordinates denoised with a Laplace filter and with the noise of the receiver point removed. The results for a hole defect are shown in Figure 13a, where the location and shape are clearly visible. The same is true for the boundaries and the position of the slag inclusion interface, as shown in Figure 13b. Figure 13c shows that, for a vertical slag inclusion defect, the upper and lower boundaries can be well imaged, as well as its bottom boundary. The position of the boundaries corresponds to those of the real model, and, although the lack of reflection point removes the information in the central region, the four boundaries are sufficient for determining the location of the slag inclusions.

Experimental Setup and Observing System
Ultrasonic seismic physics simulation laboratory equipment is used to collect the actual data underwater for pipelines with defects. Figure 14 illustrates the major components of the data acquisition equipment. The ultrasonic data acquisition system consists of a computer, an ultrasonic pulse transmitter, a water pool, high-speed data-acquisition, and a probe-motion double-3D-coordinate automatic positioning control system. The experimental procedure is as follows: The computer sets the parameters for the sampling points, the starting and ending position of the receiver probe, and so on. Then, the 3D positioning device moves between the sampling start and end points. Once the ultrasonic receiving transducer reaches a sampling point, the ultrasonic pulse generator transmits a fixed-length synchronous signal. Finally, the ultrasonic signal received by the ultrasonic receiving transducer is sent to the computer for processing.
The experimental specimen is shown in Figure 14c. The submarine pipeline length is 550 mm, the outside diameter is 219 mm, and the wall thickness is 45 mm. In order to simulate slag inclusions in actual an ultrasonic testing, cement is injected into the crack to simulate a low velocity body. The slag inclusion on the outer wall of the pipeline measures 60 mm in arc length, 3 mm in width, and 22.5 mm in depth. As shown in Figure 14d, for this slag entrapment model, we designed an observation system on the outer wall of the pipeline, with a total of seven lines. The spacing between each line is 1 cm, the length of each line is 150 mm, the channel spacing is 1 mm, and the sampling time is 60 µs. The source is located at the end point, and the offset distance is 2 cm. The ultrasonic pulse generation receiver frequency is 0.5 MHz.  Figure 15 illustrates the actual filtered data. The direct wave P can be clearly seen, as well as the surface wave R, and slag inclusions diffraction wave Pap. Since water is used as a coupling agent, the ultrasonic excitation propagates in the pipeline, and also in the water, thus producing a direct wave Pw, whose speed is much slower than the propagation in the steel pipe (medium speed), but with greater energy. We have also constructed a pipeline model with the same slag inclusion defect, matching the size and physical parameters of the test. The simulation data (Figure 16a) are compared to the experimental data ( Figure 16b). The figure illustrates that the diffracted wave Pap position of the defect is basically the same. In the experimental data, we can see the effects of t noise and absorption by the water layer, resulting in a low excitation frequency, incomplete wave field separation, and R energy covering part of the effective wave field. However, the diffraction wave Pap of the slag inclusions can still be identified, thereby confirming the reliability of our cylindrical FDM simulations in providing theoretical guidance for the ultrasonic nondestructive testing of pipeline defects. Based on the actual ultrasound Line-1~Line-7 data, cross-correlation imaging conditions are used for RTM imaging after data processing. The imaging results after noise processing are shown in Figure 17b. The image clearly shows the lower interface and vertical boundaries of the defect, and the wave field energy is concentrated at the real cementfilled defects. The results confirm the reliability of the method proposed in this paper, and that ultrasonic testing of submarine pipelines is feasible.

Conclusions
In this study, we have proposed a staggered-grid FDM in cylindrical coordinates to match the natural cylindrical symmetry of submarine pipeline cavities. A realistic pipeline model has been designed by including a double free-surface boundary and an absorbing boundary. Using this scheme, we have simulated the ultrasonic wave fields resulting from three types of defects in pipelines, i.e., a hole, a vertical slag inclusion defect, and a horizontal slag inclusion. Numerical and experimental examples are provided to verify the reliability and accuracy of the method. Signal processing has been performed using the ultrasonic information at hand. In particular, by using cross-correlation imaging conditions, 3D RTM imaging of a pipeline s space in cylindrical coordinates has been realized. Compared to traditional 2D ultrasonic testing methods, our scheme is capable of providing accurate, high-quality 3D imaging of pipeline defects with high resolution and accuracy. The improvement in accuracy comes from taking into account the converted wave on the inner wall of the pipelines and the information coming from the multiple waves reflection at the bottom of the defects. Numerical and experimental results indicate that the method is effective, and that it may be potentially applied to the practical ultrasonic nondestructive testing of submarine pipelines.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.