Friday, February 17, 2017

Don’t ignore seismic attenuation

Based on agilelibre

So good remind that is nice to refer and link this nice article and illustration about the low-pass diffusion filtering effect of earth.

It is trick that the inherent process of wave propagation from a punctual source on 2D or 3D implies in energy loss by dispersion in dimensionality know as inverse-square law.  In 3D $ 1/r $.

Monday, November 7, 2016

Regex patterns for ANP well files

For managing, renaming and moving multiple files at once. You can use:

  1. rename (Perl utility not in all Linux distros), find and grep. Some simple reminds:
    • lists all files (f) in the folder "folder" and its subfolders: find folder/ -type f
    • shows the matches of file names by the grep regex Perl expression: find folder/ -type f | grep -P '(\d-\w+-(?:\d+|\d+\w)-ES)' --color
    • remove duplicate well names in file names appended with '-' (like 1-RSM-1A-ES-1-RSM-1A-ES.lis)  find arranging/ -type f -exec rename -n 's/[-](\d-\w+-(?:\d+|\d+\w)-ES)//g' {} ";"
  2. use python 3+ glob and re modules

For well data files from BDEP ANP (pre-98 or after) bellow is a useful REGEX. It finds any file containing a BDEP ANP standard well name mixed almost anything. It splits the name in 4 groups.
First a number from 1 to 9. Second a group of characters minimum 2. Third the number of the well. Forth the state where the well was drilled, two characters. A separator between these fields may exist many times or may not even exist.

Python style raw string:



Tuesday, October 22, 2013

Equivalent Staggered Grid Acoustic Wave Equation

Equivalent Staggered grid  (Di Bartolo et al., 2012) comes from the first order system of equations that come from the continuity equation and Newton's Second Law. First the vectorial linear equation of the motion for the particle velocity field $\mathbf{v}$ is given by:

\rho( \mathbf{r}) \frac{\partial \mathbf{v}(\mathbf{r},t)}{\partial t} + \mathbf{\nabla} p(\mathbf{r},t) = \mathbf{f}(\mathbf{r},t)

where $\mathbf{r}$ is the position vector, $t$ is the time, $\rho$ os the density of the medium, $\mathbf{f}$ is the density of external body force applied and $\mathbf{\nabla} $ is the gradient operator.
The equivalent scalar equation for pressure $p$ is given by:

\frac{1}{\kappa( \mathbf{r})} \frac{\partial p(\mathbf{r},t)}{\partial t} + \mathbf{\nabla} \cdot \mathbf{v}(\mathbf{r},t) = \frac{\partial i_v(\mathbf{r},t)}{\partial t}

where $\kappa$ is the adiabatic compression modulus of the medium and $i_v$ , which represents a source, is a volume density injection (an air gun, for example). The $\mathbf{\nabla}$ applies the divergence to the vectorial velocity field $\mathbf{v} $

The Staggered Grid schemes have the great advantage of being able to model regions of water-rock interface, allowing application in the presence of materials with any Poisson's ration (including water), even when density contrasts are considered. In summary are able to model coupling between elastic and acoustic media.

Di Bartolo et. al. 2011 gives an overview of the Staggered Grid methods and its advantages against non-Staggered Methods. Finally its own formulation defined as Equivalent Staggered Grid is based on the following versions of equations (1) and (2), manipulated to remove the vectorial field $\mathbf{v} $.

\rho( \mathbf{r}) \mathbf{\nabla} \cdot \left( \frac{1}{\rho (\mathbf{r})} \mathbf{\nabla} p(\mathbf{r},t) \right) - \frac{1}{c^2(\mathbf{r})} \frac{\partial^2 p(\mathbf{r},t)}{\partial^2 t} = -s(\mathbf{r},t)


s(\mathbf{r},t) = \rho( \mathbf{r}) \frac{\partial^2 i_v(\mathbf{r},t)}{\partial^2 t}-\rho( \mathbf{r}) \mathbf{\nabla} \cdot \left(\frac{1}{\rho(\mathbf{r})} \mathbf{f}(\mathbf{r},t)\right)

and with the propagation velocity

c(\mathbf{r}) = \sqrt{\frac{\kappa(\mathbf{r})}{\rho(\mathbf{r})}}

For practical implementation the source term $s(\mathbf{r},t)$ is simplified to $ \rho( \mathbf{r}) \frac{\partial^2 i_v(\mathbf{r},t)}{\partial^2 t} $ using a null density force $ \mathbf{f}(\mathbf{r},t) $

So the equation (3) can be written as :

\rho( \mathbf{r}) \mathbf{\nabla} \cdot \left( \frac{1}{\rho (\mathbf{r})} \mathbf{\nabla} p(\mathbf{r},t) \right) + s(\mathbf{r},t) &=& \frac{1}{c^2(\mathbf{r})} \frac{\partial^2 p(\mathbf{r},t)}{\partial^2 t} \nonumber \\
c^2(\mathbf{r}) \rho( \mathbf{r}) \mathbf{\nabla} \cdot \left( \frac{1}{\rho (\mathbf{r})} \mathbf{\nabla} p(\mathbf{r},t) \right) + c^2(\mathbf{r}) \rho( \mathbf{r}) \frac{\partial^2 i_v(\mathbf{r},t)}{\partial^2 t} &=& \frac{\partial^2 p(\mathbf{r},t)}{\partial^2 t} \\
\kappa(\mathbf{r}) \mathbf{\nabla} \cdot \left( \frac{1}{\rho (\mathbf{r})} \mathbf{\nabla} p(\mathbf{r},t) \right) + \kappa(\mathbf{r}) \frac{\partial^2 i_v(\mathbf{r},t)}{\partial^2 t} &=& \frac{\partial^2 p(\mathbf{r},t)}{\partial^2 t}

From (6) the finite differences marching scheme is defined using the grid properties are $ c(\mathbf{r}) $ and $ \rho( \mathbf{r}) $ and $i_v(\mathbf{r},t)$ can be any combination of punctual sources. The schema is of 2nd order in time using central differences for the second derivative of pressure related to time

For each dimension the divergence term in (6) is written as:

$$ \mathbf{\nabla} \cdot \left( b (\mathbf{r}) \mathbf{\nabla} p(\mathbf{r},t) \right) = \frac{\partial}{\partial x}\left( b (\mathbf{r}) \frac{\partial p}{\partial x}\right)$$

Translated to finite differences and equivalent staggered grid takes the form:

$$ \frac{\partial}{\partial x}\left( b (\mathbf{r}) \frac{\partial p}{\partial x}\right)^n_{i} $$
with $b$ called buoyancy the inverse of density.

=\frac{9}{8} \left[ \frac{b_{i+\frac{1}{2}}}{\Delta x} \left( \frac{9}{8} \frac{ p^n_{i+1}-p^n_i }{\Delta x} - \frac{1}{8} \frac{ p^n_{i+2}-p^n_{i-1} }{3\Delta x}\right) -\frac{b_{i-\frac{1}{2}}}{\Delta x} \left( \frac{9}{8} \frac{ p^n_{i}-p^n_{i-1}}{\Delta x} - \frac{1}{8} \frac{ p^n_{i+1}-p^n_{i-2} }{3\Delta x}\right) \right]
-\frac{1}{8} \left[ \frac{b_{i+\frac{3}{2}}}{3\Delta x} \left( \frac{9}{8} \frac{ p^n_{i+2}-p^n_{i+1} }{\Delta x} - \frac{1}{8} \frac{ p^n_{i+3}-p^n_{i} }{3\Delta x}\right) -\frac{b_{i-\frac{3}{2}}}{3\Delta x} \left( \frac{9}{8} \frac{ p^n_{i-1}-p^n_{i-2}}{\Delta x} - \frac{1}{8} \frac{ p^n_{i}-p^n_{i-3} }{3\Delta x}\right) \right]

better form by Andre

=\frac{9}{64 \Delta x ^2} \left[ b_{i+\frac{1}{2}} \left( 27 \frac{ p^n_{i+1}-p^n_i }{3} - \frac{ p^n_{i+2}-p^n_{i-1} }{3}\right) -b_{i-\frac{1}{2}} \left( 27 \frac{ p^n_{i}-p^n_{i-1}}{3} - \frac{ p^n_{i+1}-p^n_{i-2} }{3}\right) \right]
-\frac{1}{64 \Delta x ^2} \left[ \frac{b_{i+\frac{3}{2}}}{3} \left( 27 \frac{ p^n_{i+2}-p^n_{i+1} }{3} - \frac{ p^n_{i+3}-p^n_{i} }{3}\right) -\frac{b_{i-\frac{3}{2}}}{3} \left( 27 \frac{ p^n_{i-1}-p^n_{i-2}}{3} - \frac{ p^n_{i}-p^n_{i-3} }{3}\right) \right]

=\frac{3}{64 \Delta x ^2} \left[ b_{i+\frac{1}{2}} \left( 27 ( p^n_{i+1}-p^n_i ) - p^n_{i+2}+p^n_{i-1} \right) -b_{i-\frac{1}{2}} \left( 27 (p^n_{i}-p^n_{i-1}) - p^n_{i+1}+p^n_{i-2} \right) \right]
-\frac{1}{576\Delta x ^2} \left[ b_{i+\frac{3}{2}} \left( 27 (p^n_{i+2}-p^n_{i+1} ) - p^n_{i+3}+p^n_{i} \right) -b_{i-\frac{3}{2}} \left( 27 (p^n_{i-1}-p^n_{i-2}) -p^n_{i}+p^n_{i-3} \right) \right]

Then opening $b$ using linear interpolation

=\frac{3}{128 \Delta x ^2} \left[ (b_{i+1}+b_i) \left( 27 ( p^n_{i+1}-p^n_i ) - p^n_{i+2}+p^n_{i-1} \right) -(b_{i-1}+b_i) \left( 27 (p^n_{i}-p^n_{i-1}) - p^n_{i+1}+p^n_{i-2} \right) \right]
-\frac{1}{1152\Delta x ^2} \left[ (b_{i+2}+b_{i+1}) \left( 27 (p^n_{i+2}-p^n_{i+1} ) - p^n_{i+3}+p^n_{i} \right) -(b_{i-2}+b_{i-1}) \left( 27 (p^n_{i-1}-p^n_{i-2}) -p^n_{i}+p^n_{i-3} \right) \right]

 Di Bartolo, L. A new family of finite-difference schemes to solve the heterogeneous acoustic wave equation - Geophysics 2012

Friday, August 9, 2013

Git and Git Hub reminds

Some short reminds:

Going back to a previous push :

git reset --hard ce82c8108042e555b57013 
where ce82c8108042e555b57013 is the SHA key 

git push -f origin

Getting all (remote) branches to (local)

git fetch origin

Cloning a (remote) branch to (local)

git checkout -b  ''  origin/''


Thursday, August 1, 2013

Simple explicit 2D scalar wave equation with convergence

Easy as a cake, after we learn.

Explicit 2D scalar wave equation

The scalar wave equation, a fluid medium with variable velocity, is defined by (\ref{ExpA}) Where $U$ and $V$ represents the pressure (or scalar displacement) and velocity fields and S(t) the source, and can be simplified in (\ref{ExpB}) for 2D case:

\frac{\partial^2 U}{\partial t^2} \frac{1}{V^2} &=& \nabla^2 U + S(t) \label{ExpA} \\
\frac{\partial^2 U}{\partial t^2} &=& V^2 \left( \frac{\partial^2 U}{\partial x^2}+ \frac{\partial^2 U}{\partial z^2} +S(t) \right)

From equation (\ref{ExpB}) we can approximate derivatives by Taylor series, manipulating we can reach backward differences in time (first order time) and centered in space N'th order; using $Dx=Dz=Ds$, remember that $j=jDs$, $k=kDs$ and $n=nDt$.

U_{jk}^{n+1} &=& \left( \frac{\Delta t V_{jk} }{\Delta s} \right) ^2 \left( \sum_{a=-N}^N w_a U_{j+a k}^n + \sum_{a=-N}^N w_a U_{j k+a}^n +S_{jk}^n {\Delta s}^2 \right) + 2 U_{jk}^{n} - U_{jk}^{n-1} \\
U_{jk}^{n+1} &=& \left( \frac{\Delta t V_{jk}}{\Delta s} \right) ^2 \left[ \sum_{a=-N}^N w_a \left( U_{j+a k}^n + U_{j k+a}^n \right) \right] + S_{jk}^n {\Delta t}^2 V_{jk}^2 + 2 U_{jk}^{n} - U_{jk}^{n-1} \\
U_{jk}^{n+1} &=& \left( \frac{\Delta t V_{jk}}{\Delta s} \right) ^2 \left[ \sum_{a=-N}^N w_a \left( U_{j+a k}^n + U_{j k+a}^n \right) + S_{jk}^n {\Delta s}^2\right] + 2 U_{jk}^{n} - U_{jk}^{n-1} \label{ExpC}

For forth order space, we have $N=2$ and $w$ is:
$$ w = \frac{1}{12} [-1, 16, -30, 16, -1] $$

To make the approximate finite differences solution (or series) above converge to the correct solution  $ R $ must always be bounded (the stability criteria):

\[ R = \frac{\Delta t V_{max}}{\Delta s} \]

Where $V_{max}$ is the maximum velocity.
Based on work and simplified from Chen [2] the forth order schema requires :

$$ \Delta t \leq \frac{2 \Delta s}{ V_{max} \sqrt{2/12(1+16+30+16+1)}} \implies \Delta t \leq \frac{ \Delta s \sqrt{3}}{ V_{max} \sqrt{8}} = \frac{ \Delta s \sqrt{3}}{ V_{max} \sqrt{8}}$$

Must also be highlighted that the closer $\Delta t$ is to the limit above better is the solution.

  • Last line in equation  (\ref{ExpC}) describes the time step forward implementation (second order in time and forth in space) for a 2D grid.
  • To guarantee stability the Laplacian should not be contaminated, just multiply by R.
  • Use an epsilon $\epsilon$ (like 10^-4) to make $\Delta t$ closer to its limit of convergence. Like $\Delta t = \frac{ \Delta s \sqrt{3}}{ V_{max} \sqrt{8}}*(1-\epsilon) $ 


Simple cookbook example from Fatiando branch (under development) here. The animation bellow shows a receiver recorded at a diffraction point from a Ricker wavelet source .


[1] Alford R.M., Kelly K.R., Boore D.M. (1974) Accuracy of finite-difference modeling of the acoustic wave equation Geophysics, 39 (6), P. 834-842

[2] Chen, Jing-Bo (2011) A stability formula for Lax-Wendroff methods with fourth-order in time and general-order in space for the scalar wave equation Geophysics, v. 76, p. T37-T42

Wednesday, May 22, 2013

Digging out end up in Javascript

Ok, I accept this has nothing to do with the previous post ... but anyway since I found this code, I did in college years ago about linear momentum, and translated it to JavaScript using Html5 canvas and found it cool, here I post it just for fun...

This text is displayed if your browser does not support HTML5 Canvas.

Tuesday, May 14, 2013

Miss behavior MathJax and Blogger

If you face any problem seeing equations (numbering), try to visualize the post you want individually.
Unfortunately the numbering of equations gets messed up when viewing multiple posts at once.

Sorry :S For now no workaround!