Skip to content

Commit 8bae7f2

Browse files
JDBetteridgeJDBetteridge
committed
misc: Address George's comments and tidy a bit
1 parent 76e711e commit 8bae7f2

1 file changed

Lines changed: 52 additions & 22 deletions

File tree

‎examples/seismic/tutorials/07.1_dispersion_relation.ipynb‎

Lines changed: 52 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -9,21 +9,13 @@
99
"This notebook compares the numerical dispersion of explicit finite-difference stencils for the second derivative obtained via several methods.\n",
1010
"\n",
1111
"Contents:\n",
12-
"1. **Problem setup**\n",
13-
"2. **Finite Difference Weights from Fornberg’s Algorithm** computed via the Fornberg algorithm. [1]\n",
14-
"3. **Investigation of the dispersion properties** of the weights.\n",
15-
"4. **DRP (Dispersion-Relation-Preserving) stencils (part 1)** following Tam and Webb [2], Caunt [3]\n",
16-
"5. **DRP stencils (part 2)** following Chen, Peng and Li [4]\n",
12+
"1. **[Problem setup](#1.-Problem-Setup)**\n",
13+
"2. **[Finite Difference Weights from Fornberg’s Algorithm](#2.-Finite-Difference-Weights-from-Fornberg’s-Algorithm)** computed via the Fornberg algorithm. [[1]](#References)\n",
14+
"3. **[Investigation of the dispersion properties](#3.-Investigation-of-the-dispersion-properties)** of the weights.\n",
15+
"4. **[DRP (Dispersion-Relation-Preserving) stencils (part 1)](#4.-Dispersion-Relation-Preserving-Stencils-(part-1))** following Tam and Webb [[2]](#References), Caunt [[3]](#References)\n",
16+
"5. **[DRP stencils (part 2)](#5.-Dispersion-Relation-Preserving-Stencils-(part-2))** following Chen, Peng and Li [[4]](#References)\n",
1717
"\n",
18-
"Both 3. and 4. above are achieved by computed by solving different constrained least-squares optimization problems.\n",
19-
"\n",
20-
"[1] https://www.colorado.edu/amath/sites/default/files/attached-files/mathcomp_88_fd_formulas.pdf\n",
21-
"\n",
22-
"[2] https://www.sciencedirect.com/science/article/pii/S0021999183711423\n",
23-
"\n",
24-
"[3] https://arxiv.org/pdf/2107.13525\n",
25-
"\n",
26-
"[4] https://www.sciencedirect.com/science/article/pii/S009830042100234X"
18+
"Both 3. and 4. above are achieved by computed by solving different constrained least-squares optimization problems."
2719
]
2820
},
2921
{
@@ -296,7 +288,7 @@
296288
"source": [
297289
"## 3. Investigation of the dispersion properties\n",
298290
"\n",
299-
"Chen, Peng and Li [4] provide the dispersion error in the form of a ratio\n",
291+
"Chen, Peng and Li [[4]](#References) provide the dispersion error in the form of a ratio\n",
300292
"$$\n",
301293
"\\frac{v_\\text{FD}}{v} = \\delta(r_a, r, \\alpha, \\beta),\n",
302294
"$$\n",
@@ -464,7 +456,7 @@
464456
"metadata": {},
465457
"source": [
466458
"## 3.1. Observing the effects on wavefields\n",
467-
"We now use Devito to solve the acoustic wave equation, propagating a Ricker wavelet source from the centre. We define a function that solves the acoustic wave equation in a 2D grid with reflecting boundary conditions. Recievers are also placed near the top of the domain allowing us to later generate a shot profile. The source in the centre of the 1km by 1km domain is simulated for at total of 0.6 seconds, long enough for the wave to pass the recievers and be reflected by the boundary\n",
459+
"We now use Devito to solve the acoustic wave equation, propagating a Ricker wavelet source from the centre. We define a function that solves the acoustic wave equation in a 2D grid with reflecting boundary conditions. Receivers are also placed near the top of the domain allowing us to later generate a shot profile. The source in the centre of the 1km by 1km domain is simulated for at total of 0.6 seconds, long enough for the wave to pass the receivers and be reflected by the boundary\n",
468460
"\n",
469461
"Note the use of the weights keyword argument passed to the derivative, `u.dx2(weights=weights)`, allowing the use of custom stencil weights."
470462
]
@@ -751,7 +743,7 @@
751743
"\n",
752744
"Note that in the unstable case there is a wavefront that is somehow ahead of the original wavelet! The large peak of the wavelet is also growing in amplitude and as a result the solution is not physically accurate. In the unstable case the first arrival time is also affected, and this becomes very obvious when looking at shot profiles.\n",
753745
"\n",
754-
"Again the plots from left to right go from the smallest Courtant number to largest. The correspinging wavelet is plotted along the green dashed line and the **true** first arrival time is clearly marked on all of the plots."
746+
"Again the plots from left to right go from the smallest Courant number to largest. The correspinging wavelet is plotted along the green dashed line and the **true** first arrival time is clearly marked on all of the plots."
755747
]
756748
},
757749
{
@@ -848,7 +840,7 @@
848840
"cell_type": "markdown",
849841
"metadata": {},
850842
"source": [
851-
"We are now free to choose an objective function for the stencil weights to optimise over. We follow the method of Tam and Webb [2] and executed in the Masters thesis of Caunt [3] (for second derivatives) and minimise the $L^2$ norm of the difference between the derivative and the stencil in Fourier space. We minimise $\\Phi(a_m)$, where\n",
843+
"We are now free to choose an objective function for the stencil weights to optimise over. We follow the method of Tam and Webb [2](#References) and executed in the Masters thesis of Caunt [3](#References) (for second derivatives) and minimise the $L^2$ norm of the difference between the derivative and the stencil in Fourier space. We minimise $\\Phi(a_m)$, where\n",
852844
"$$\n",
853845
"\\Phi(a_m):= \\int_{\\pi/2}^{\\pi/2} \\left| \\varphi^2 + \\sum_{m=-M}^M a_m e^{im\\varphi} \\right|^2 d\\varphi,\n",
854846
"$$\n",
@@ -981,7 +973,7 @@
981973
"source": [
982974
"Both sets of weight approximate the curve well, but it is clear from the right hand plot that the error is smaller for DRP stencil 1 over a larger range; less than $10^-4$ in the whole range $[0, \\pi/2]$.\n",
983975
"\n",
984-
"We can also look at the velocity error ratio for this set of coefficients, as outlined by Chen, Peng and Li [4] and we did above for the Fornberg weights."
976+
"We can also look at the velocity error ratio for this set of coefficients, as outlined by Chen, Peng and Li [4](#References) and we did above for the Fornberg weights."
985977
]
986978
},
987979
{
@@ -1310,7 +1302,7 @@
13101302
"name": "stderr",
13111303
"output_type": "stream",
13121304
"text": [
1313-
"/tmp/ipykernel_626021/237837473.py:12: RuntimeWarning: invalid value encountered in arccos\n",
1305+
"/tmp/ipykernel_904261/237837473.py:12: RuntimeWarning: invalid value encountered in arccos\n",
13141306
" diff = np.abs(np.acos(theta)/(k*dt) - v)\n"
13151307
]
13161308
},
@@ -1375,7 +1367,7 @@
13751367
"source": [
13761368
"It is, but only just! This should, however, come as no surprise: We used the value of $\\Delta t$ directly in the optimisation and we obtain a stencil that minimises dispersion as much as possible for that value.\n",
13771369
"\n",
1378-
"Next we plot the veloicty error ratio for DRP stencil 2:"
1370+
"Next we plot the velocity error ratio for DRP stencil 2:"
13791371
]
13801372
},
13811373
{
@@ -1450,7 +1442,7 @@
14501442
"cell_type": "markdown",
14511443
"metadata": {},
14521444
"source": [
1453-
"The two DRP stencils produce curves that are indistinguihable in the left hand plot, the difference is only really noticable in the logarithmic difference, where we se larger error for DRP stencil 2 than DRP stencil 1, while still remaining small.\n",
1445+
"The two DRP stencils produce curves that are indistinguishable in the left hand plot, the difference is only really noticable in the logarithmic difference, where we see larger error for DRP stencil 2 than DRP stencil 1, while still remaining small.\n",
14541446
"\n",
14551447
"We can also look at the isosurfaces of the velocity difference function $\\hat{\\delta}$, the expression that the objective function tried to minimise."
14561448
]
@@ -1598,6 +1590,44 @@
15981590
"source": [
15991591
"Again, the ripples following the wavelet have almost entirely disappeared. The unstable $\\Delta t$ still has an incorrect first arrival wave, but this is expected."
16001592
]
1593+
},
1594+
{
1595+
"cell_type": "markdown",
1596+
"metadata": {},
1597+
"source": [
1598+
"## References\n",
1599+
"[1] **Generation of Finite Difference Formulas on Arbitrarily Spaced Grids** (1988)\n",
1600+
"\n",
1601+
"_Bengt Fornberg_\n",
1602+
"\n",
1603+
"Mathematics of Computation, Vol. 51, No. 184\n",
1604+
"\n",
1605+
"http://dx.doi.org/10.1090/S0025-5718-1988-0935077-0, https://www.colorado.edu/amath/sites/default/files/attached-files/mathcomp_88_fd_formulas.pdf\n",
1606+
"\n",
1607+
"[2] **Dispersion-Relation-Preserving Finite Difference Schemes for Computational Acoustics** (1993)\n",
1608+
"\n",
1609+
"_Christopher K.W. Tam, Jay C. Webb_\n",
1610+
"\n",
1611+
"Journal of Computational Physics, Volume 107, Issue 2, Pages 262-281, ISSN 0021-9991\n",
1612+
"\n",
1613+
"https://doi.org/10.1006/jcph.1993.1142, https://www.sciencedirect.com/science/article/pii/S0021999183711423\n",
1614+
"\n",
1615+
"[3] **Spatially-optimized finite-difference schemes for numerical dispersion suppression in seismic applications** (2019)\n",
1616+
"\n",
1617+
"_Edward Caunt_\n",
1618+
"\n",
1619+
"Masters Thesis\n",
1620+
"\n",
1621+
"http://dx.doi.org/10.48550/arXiv.2107.13525, https://arxiv.org/pdf/2107.13525\n",
1622+
"\n",
1623+
"[4] **A framework for automatically choosing the optimal parameters of finite-difference scheme in the acoustic wave modeling** (2022)\n",
1624+
"\n",
1625+
"_Guiting Chen, Zhenming Peng, Yalin Li_\n",
1626+
"\n",
1627+
"Computers & Geosciences, Volume 159, 104948, ISSN 0098-3004,\n",
1628+
"\n",
1629+
"https://doi.org/10.1016/j.cageo.2021.104948, https://www.sciencedirect.com/science/article/pii/S009830042100234X"
1630+
]
16011631
}
16021632
],
16031633
"metadata": {

0 commit comments

Comments
 (0)