# Simulation of Wave Propagation Using Finite Differences in Oil Exploration
## 2 Materials and Methods  
### Description of the Model  
The mesh was divided into N\({}_{x}\) by N points. It defined \(\Delta x\) and \(\Delta z\) as the distance between the points or nodes of the mesh. If \(\Delta t\) was the increment in time, then \(t=k\Delta t\), where k was the step-in time with \(k=1,2,\ldots,n\). With this discretization, one can describe a scheme of infinite difference to approximate Equation (5) by rewriting it as follows:  
\[\frac{\partial^{2}p}{\partial t^{2}}-\rho c^{2}\left[\frac{\partial}{\partial x }\left(\frac{1}{\rho}\frac{\partial p}{\partial x}\right)+\frac{\partial}{ \partial z}\left(\frac{1}{\rho}\frac{\partial p}{\partial z}\right)\right]= \delta(r)(t) \tag{6}\]  
The second derivative concerning time can be approximated using the equation in differences centered on the second-order [3]:  
\[\frac{\partial^{2}p}{\partial t^{2}}=\frac{p(i,j,k+1)-2p(i,.k)+p(i,j,k-1)}{ \partial t^{2}} \tag{7}\]  
where \(p(i,j,k)=p\) (\(i\Delta x,j\Delta z,k\Delta t\)). Similarly, the spatial derivative is also approximated by an equation in finite differences:  
\[\left\{\frac{\partial}{\partial x}\left[\frac{1}{\rho}\left(\frac{\partial p }{\partial x}\right)\right]\right\}_{i,j}^{k} =\frac{\rho^{-1}(i+\frac{1}{2}j)}{\Delta x}\left[\frac{\partial p }{\partial x}\right]_{i+\frac{1}{2}j}^{k}-\frac{\rho^{-1}(i-\frac{1}{2}j)}{ \Delta x}\left[\frac{\partial p}{\partial x}\right]_{i-\frac{1}{2}j}^{k} \tag{8}\] \[=\frac{\rho^{-1}(i+\frac{1}{2}j)}{\Delta x}\left[\frac{p(i+1,j,k) -p(i,j,k)}{\Delta x}\right]-\frac{\rho^{-1}(i-\frac{1}{2}j)}{\Delta x}\left[ \frac{p(i,j,k)-p(i-1,j,k)}{\Delta x}\right]\]  
Moreover, this is similar for the component in z, which is given by:  
\[\left\{\frac{\partial}{\partial z}\left[\frac{1}{\rho}\left(\frac{\partial p }{\partial z}\right)\right]\right\}_{i,j}^{k} =\frac{\rho^{-1}(i_{j}+\frac{1}{2})}{\Delta z}\left[\frac{\partial p }{\partial z}\right]_{i+\frac{1}{2}j}^{k}-\frac{\rho^{-1}(i_{j}-\frac{1}{2})}{ \Delta z}\left[\frac{\partial p}{\partial z}\right]_{i,j-\frac{1}{2}}^{k} \tag{9}\] \[=\frac{\rho^{-1}(i_{j}+\frac{1}{2})}{\Delta z}\left[\frac{p(i,j+1,k)-p(i,j,k)}{\Delta z}\right]-\frac{\rho^{-1}(i_{j}-\frac{1}{2})}{\Delta z} \left[\frac{p(i,j,k)-p(i,j-1,k)}{\Delta z}\right]\]By substituting Equations (7)-(9) into (6), simplifying and taking the value \(\rho^{-1}\) (\(i\pm 1/2\), \(j\)) as the average value between two nodes, one has to solve it for future time steps:  
\[ \begin{split} p(i,j,k+1)=& 2p(i,j,k)-p(i,j,k-1)\\ &+\frac{\rho(i,j)c^{2}(i,j)\Delta t^{2}}{2\Delta x^{2}}x\Big{\{} \Big{[}\rho^{-1}(i+1,j)+\rho^{-1}(i,j)\Big{]}[p(i+1,j,k)-p(i,j,k)]\\ &-\Big{[}\rho^{-1}(i,j)+\rho^{-1}(i-1,j)\Big{]}[p(i,j,k)-p(i-1,j, k)]\\ &+\bigg{(}\frac{\Delta x}{\Delta z}\bigg{)}^{2}\Big{[}\rho^{-1}(i,j+1)+\rho^{-1}(i,j)\Big{]}[p(i,j+1,k)-p(i,j,k)]\\ &-\Big{[}\rho^{-1}(i,j)+\rho^{-1}(i,j-1)\Big{]}[p(i,j,k)-p(i,j-1,k)]\Big{\}}+\delta^{i,j}_{is,j}f(k)\end{split}  \tag{10}\]  
When performing simulations that involve partial differential equations, it became evident that, aside from the rounding and truncation errors of the series (bearing in mind that the expressions of finite differences were derived from a second-order Taylor series development), there were errors inherent in the process due to the challenges of accurately adapting the physical model to describe the natural physical system and the limitations of the boundary conditions. Convergence and stability tests ensured that the finite difference approximations used in discretization tests were correct. There are ways to calculate the convergence of a code when analytical solutions are known. Once the value of \(\Delta x\) is known, the minimum number of nodes can be determined to present a numerical solution for this system. The Courant-Friedrichs-Lewy (CFL) condition is necessary to achieve a numerical solution for differential equations in partial derivatives' convergence [20].