Fast And High Accurate Algorithm And Its Implementation For Acoustic Waves And Elastic Waves In Reverse Time Migration

Reverse time migration(RTM)is a high-efficiency tools for seismic imaging.RTM has been considered to be one of the most accurate seismic pre-stack depth migration methods,especially for imaging geologically complex structures.RTM calculates the two-way wave equations numerically.It has no angle limit,and has considered complex propagating paths,such as turning,prismatic and multiple reflections.RTM can get better imaging quality than ray based Kirchhoff migration and one-way equation based miagraiton methods.However,RTM is computationally expensive for mass computation and storage to solve the two-way wave equation numerically.The bottleneck limits the widely industrial application of new methods.To ease the bottleneck,two ways can be considered.The one is using fast and accurate numerical method to solve the two-way wave equations,and thus reduce the amounts of the computation resource.The other is using high-performance hardware and optimized code to improve the efficiency of RTM.The core procedure of RTM is sovling the two-way equations.So an efficient numerical method to sovle the equations accurately and fast is important to the efficiency of RTM.We propose to use the high-order finite difference time domain method(FDTD)and the pseudospectral time domain method(PSTD)to improve the the efficiency of RTM.Traditional RTM usually uses 2nd-order central difference schema.The 2nd-order central difference schema is easy to implement but has shortcoming of high sample rate and computationally intensive.Using a high-order difference schema in space can significant reduce the spatial sample rate and grid size,thus saves mass computation time and storage.In a 3D model,using 12th-order difference in space instead of 2nd-order difference can save up to 96.3%storage and mass computation time.However,the needed computation resource of RTM is still enormous even if we use the high-order difference.The PSTD method is an effective method to further reduce the requirements.Traditional pseudospectral method uses the fast Fourier transform(FFT)algorithm to calculate the spatial derivatives,but is limited by the wraparound effect due to the periodicity assumed in the FFT.The PSTD algorithm combines the pseudospectral method with a perfectly matched layer(PML)for acoustic waves.PML is a highly effective absorbing boundary condition that can eliminate the wraparound effect.It enables a wide application of the pseudospectral method to complex models.RTM based on the PSTD algorithm has advantages in the computational efficiency compared to traditional methods such as the second-order and high order finite difference time-domain(FDTD)methods.The spatial sample rate of PSTD method is ony 2 points per minimum wavelength at the highest frequency,which is only one sixth of the 2nd-order FDTD and half of the 12th-order FDTD.For a 3D model,in the same accuracy,the PSTD method needs a grid size which is only 0.5%of the 2nd-order FDTD and 12.5%of the 12th-order.The computational requirements of RTM are enormous;to overcome this limitation,large scale CPU clusters(homogeneous hardware)have been commonly used for RTM.Besides,"heterogeneous" hardware can also be applied to RTM,such as the IBM Cell/BE processor,Field Programmable Gate Array(FPGA),and the Graphics Processing Unit(GPU).Among all platforms mentioned above,CPU and GPU have advances of stable performance,good universality and easyly programming.We use OpenMP in CPU and CUDA in GPU to speedup the RTM based on the PSTD method and the high-order FDTD method.Our RTM based on on the PSTD method and the high-order FDTD method is suitable to parallel computing.RTM based on on the PSTD method and the high-order FDTD can significant reduce the grid size and it is important to a memory-limited GPU card.RTM based on the PSTD algorithm and the high-order FDTD mdthod are particularly suitable for running on GPU because it can solve a much larger problem than the traditional FDTD method on a memory-limited GPU.To secure a better performance and generality of FFT in GPU,we present a scheme which combines 1D FFT with matrix transpositions instead of using 3D FFT directly,and get 30%performance improvement in pesudospectral derivative operater.The matrix transpositions use shared memory to improve memory access efficiency.We also apply an efficient FFT scheme proposed by Xu et al.which replaces even-sized R2C FFT with a half-sized C2C FFT,and get further 20%performance improvement in pesudospectral derivative operater.For a small amount and balanced memory swapping from computer to GPU,we save the boundaries in lieu of checkpointing scheme when we propagate the source wavefield forward and backward.The proposed RTM has an acceleration ratio of about 80 times by a Tesla K20X GPU card on a desktop computer.The simulation results of 2D and 3D models demonstrate that the proposed RTM is fast and inexpensive.