Open codes

$$\tag{Content}\label{Content}$$

Content

• \( \ref{KWSINC} \): a package for arbitrary source and receiver positioning in finite-difference grids and seismic data resampling with Kaiser-windowed sinc functions .
• \( \ref{FDFDMATRIX3D} \) : a package for 3D finite-difference frequency-domain seismic wave modeling in the VTI visco-acoustic approximation.

$$\tag{KWSINC}\label{KWSINC}$$

KWSINC : a package for arbitrary source and receiver positioning in finite-difference grids and seismic data resampling with Kaiser-windowed sinc functions.

This report reviews how to perform seismic data resampling in time and implement source and receiver positioning in finite difference schemes with Kaiser-windowed sinc (KWS) function as proposed by\citet{Hicks_2002_ASR}. We provide a package of Fortran codes that implement 1D, 2D and 3D KWS functions in a fortran module and interface them with finite-difference time-domain (FDTD) and finite-difference frequency-domain (FDFD) codes for numerical validation of the source and receiver positioning. I found over years that this approach is tedious to implement and can be a source of errors and bugs, which prompts me to write this report with as many details as possible. The KWS method implements a monopole or dipole source/receiver at arbitrary position in a uniform finite difference grid by approximating a delta function and its derivative by a KWS function, namely a truncated limited-bandwidth representation of a Delta function.
In our codes, we first tabulate the values of the KWS function in the domain alpha in ]-0.5;+0.5], where alpha is the algebraic distance of the closest grid point from the source/receiver in a local coordinate system whose origin is the position of the source/receiver. This tabulation allows us to avoid recomputing the KWS coefficients for each source/receiver position each time we need to perform a wave simulation (for example, at each iteration of Full Waveform Inversion (FWI)). We also review how to process the free surface by mirroring the coefficients that are above the free surface with a sign -1. The algorithm is validated for different scenarios in a homogeneous half space against analytical solution. Our implementation rigorously follows the indexing convention defined by Hicks (2002). We also provide a code that subsample/oversample seismograms with KWS function and validate it with a numerical example. This can be useful when we want to store the values of the wavefields at the boundaries of the computational domain at the Nyquist limit before their recomputation backward in time during gradient computation during FWI application.

Contacts: Stephane Operto (email: operto@geoazur.unice.fr)

References:

• G. J. Hicks, Arbitrary source and receiver positioning in finite-difference schemes using Kaiser windowed sinc functions, Geophysics, 67, 156--166, 2002.
• S. Operto, Seismic data resampling and source/receiver positioning in finite-difference simulatio with Kaiser-windowed sinc function (Hicks' method), WIND project, Technical report n 11, 2022 (provided as the documentation of the KWSinc package).

Download KWSInc package (.tgz, 1.2 Mb)

Download KWSInc documentation (.pdf, 4Mb)

Back to \ref{Content}

$$\tag{FDFDMATRIX3D}\label{FDFDMATRIX3D}$$

FDFDMATRIX3D: a package for 3D finite-difference frequency-domain seismic wave modeling in the VTI visco-acoustic approximation.

The numerical scheme relies on the so-called 27-point mixed-grid stencil . This stencil has a compact support covering the 27 grid points of the 8 cubic cells surrounding the central point. This compact support is crucial to minimize matrix fill-in when sparse direct solver are used to solve the system after a LU decomposition of the impedance matrix. This compact support is also useful to minimize communication when the system is solved with iterative solver based upon domain decomposition preconditioner. This compact support is obtained by building the stiffness matrix with second-order accurate stencils. Accuracy is obtained by linearly combining with optimal weights different stiffness matrices build on different coordinate systems and by building a consistent mass matrix with optimal weights. The weights are estimated through a classical dispersion analysis in homogeneous media. They can be non-adaptive or adaptive . In the first case, dispersion in homogeneous media is minimized jointly for several values of G (number of grid points per wavelength) and each row of the impedance matrix is built with the same weights. In the adaptive case, the weights are functions of 1/G. These functions are computed in discrete form by minimizing dispersion separately for each value of G sampling the domain of the function (1/G belongs to [0;0.25] in our implementation). Then, each row of the impedance matrix is built with the weights corresponding to the local wavelength. We have shown that the wavelength-adaptive stencil provides a more uniform accuracy in heterogeneous media.

This package provides a MATLAB code that computes the non-adaptive and adaptive weights of the stencil by dispersion analysis, sequentiel FORTRAN codes that build the impedance matrix and sparse right-hand sides (point sources) and a MPI-PARALLEL FORTRAN code calling the sparse multifrontal solver MUMPS that computes monochromatric wavefields from the output files of the two previous programs (sparse impedance matrix and right-hand sides).

The simulations can be performed for five benchmarks (homogeneous model, linear model, 3D SEG/EAGE overthrust and salt models and GO_3D_OBS model ). The accuracy of the FDFD wavefields can be assessed agaisnt analytical solutions (homogeneous and linear models) and the wavefields computed with the convergent Born series (CBS) method. The velocity grids and the CBS wavefields are made available to perform the simulations. Codes to compute analytical solutions in homogeneous and linear models are provided in the package.

Contacts: Hossein Agahmiry (email: aghamiry@geoazur.unice.fr), Laure Combe (email: combe@geoazur.unice.fr), Stephane Operto (email: operto@geoazur.unice.fr)

References:

1. H. Aghamiry, A. Gholami, L. Combe and S. Operto, Accurate 3D frequency-domain seismic wave modeling with the wavelength-adaptive 27-point finite-difference stencil: a tool for full waveform inversion, Geophysics, https://doi.org/10.1190/geo2021-0606.1.
2. S. Operto, R. Brossier, L. Combe, L. Métivier, A. Ribodetti, and J. Virieux. Computationally-efficient visco-acoustic finite-difference frequency-domain seismic modeling in vertical transversely isotropic media with sparse direct solvers. Geophysics, 79(5), 2014, https://doi.org/10.1190/geo2013-0478.1.
3. R. Brossier, V. Etienne, S. Operto, and J. Virieux. Frequency-domain numerical modelling of visco-acoustic waves based on finite-difference and finite-element discontinuous galerkin methods. In D. W. Dissanayake, editor, Acoustic Waves, pages 125--158, SCIYO, 2010.
4. S. Operto, J. Virieux, P. Amestoy, J-Y. L'Excellent, L. Giraud, and H. Ben Hadj Ali, 3D finite-difference frequency-domain modeling of visco-acoustic wave propagation using a massively parallel direct solver: A feasibility study, Geophysics, 72(5):SM195--SM211, 2007, https://doi.org/10.1190/1.2759835.

Simulations in the 3D SEG/EAGE Salt model. Left column: CBS wavefield. Middle column: FDFD wavefield computed with wavelength-adaptive 27-point stencil. Right column: Differences. First two rows are depth slices and last two rows are vertical sections. See Aghamiry et al. (2022) for more details.

Download here FDFDMATRIX3D package (.tgz, 14Mb)

Download version 1.1 of FDFDMATRIX package (May 2022) FDFDMATRIX_v1.1 package (.tgz, 13Mb)

Download here documentation of FDFDMATRIX3D (.pdf, 18Mb)

Download here documentation of FDFDMATRIX3D_v1.1 (.pdf, 17Mb)

Updated March 7th 2022: Scripts were modified such that the environment variable $fdfdpath=../../../BIN by default if it was not defined by the user. In that case, it is assumed that the user runs the demo from the tree structure of the package. Scripts psimage.sh were corrected such that the path to the CBS solutions is set to ../CBS.

Update May 2022 (v1.1): scripts psimage.sh for plot were corrected by removing undefined environement variables. Benchmarks for validating against analytical solutions arbitrary source & receiver positionings with Kaiser-windowed sinc functions in homogeneous infinite media and half spaces are added.

Back to \ref{Content}

Number of visitors