Skip to content

examples: Refresh several notebooks in examples/ #2542

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
georgebisbas wants to merge 12 commits into
base: main
Choose a base branch
from
Open

examples: Refresh several notebooks in examples/ #2542

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
12 commits
Select commit Hold shift + click to select a range
504aa89
examples: Refresh 08_snapshotting
georgebisbas Feb 19, 2025
4bc5610
examples: Better parametrization and error messages
georgebisbas Feb 19, 2025
f2cb28b
examples: Update /examples/seismic/tutorials/08_snapshotting.ipynb
georgebisbas Feb 21, 2025
394994b
examples: Minor edits
georgebisbas Feb 25, 2025
b85b2e3
examples: Update 13/15/16
georgebisbas Mar 6, 2025
08c02b8
examples: Update bs_ivbp
georgebisbas Mar 6, 2025
b7b3eb5
examples: Rename black-scholes tutorial
georgebisbas Mar 13, 2025
000aa9f
Revert "examples: Rename black-scholes tutorial"
georgebisbas Mar 13, 2025
fc1a2c9
examples: Update by comments
georgebisbas Mar 13, 2025
74a8cd7
examples: Relax tolerance
georgebisbas Mar 13, 2025
1ff03de
examples: Further reviews and updates
georgebisbas Mar 21, 2025
016b17c
examples: More updates and fixes
georgebisbas Mar 25, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
62 changes: 25 additions & 37 deletions examples/finance/bs_ivbp.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,23 +27,22 @@
"name": "stdout",
"output_type": "stream",
"text": [
"dt,tmax,nt; 0.0005 1.0 2001\n",
"shape; (81,)\n",
"origin; (60,)\n",
"spacing; (1.0,)\n",
"extent; 80\n"
"dt,tmax,nt: 0.0005, 1.0, 2001\n",
"shape: (81,)\n",
"origin: (60,)\n",
"spacing: (1.0,)\n",
"extent: 80\n"
]
}
],
Copy link
Contributor

@EdCaunt EdCaunt Mar 11, 2025 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should this notebook be renamed "01_black_scholes" for consistency with other tutorials?


Reply via ReviewNB

Copy link
Contributor Author

@georgebisbas georgebisbas Mar 13, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would be happy with that

Copy link
Contributor

@EdCaunt EdCaunt Mar 11, 2025 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more. 

Line #19.    plt.plot(t_range, np.array(vals))

Needs a plt.show()


Reply via ReviewNB

Copy link
Contributor Author

@georgebisbas georgebisbas Mar 13, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

why?

"source": [
"from devito import (Eq, Grid, TimeFunction, Operator, solve, Constant, \n",
" SpaceDimension, configuration, SubDomain, centered)\n",
" SpaceDimension, configuration, centered)\n",
"\n",
"from mpl_toolkits.mplot3d import Axes3D\n",
"from mpl_toolkits.mplot3d.axis3d import Axis\n",
"import matplotlib.pyplot as plt\n",
"import matplotlib as mpl\n",
"from matplotlib import cm\n",
"\n",
"from sympy.stats import Normal, cdf\n",
"import numpy as np\n",
Expand Down Expand Up @@ -101,11 +100,11 @@
"spacing = (ds0, )\n",
"extent = int(ds0 * (ns - 1))\n",
"\n",
"print(\"dt,tmax,nt;\", dt0,tmax,nt)\n",
"print(\"shape; \", shape)\n",
"print(\"origin; \", origin)\n",
"print(\"spacing; \", spacing)\n",
"print(\"extent; \", extent)"
"print(f\"dt,tmax,nt: {dt0}, {tmax}, {nt}\")\n",
"print(f\"shape: {shape}\")\n",
"print(f\"origin: {origin}\")\n",
"print(f\"spacing: {spacing}\")\n",
"print(f\"extent: {extent}\")"
]
},
{
Expand Down Expand Up @@ -242,8 +241,8 @@
"name": "stderr",
"output_type": "stream",
"text": [
"Operator `Kernel` run in 0.01 s\n",
"Operator `Kernel` run in 0.01 s\n"
"Operator `Kernel` ran in 0.01 s\n",
"Operator `Kernel` ran in 0.01 s\n"
]
},
{
Expand All @@ -263,7 +262,7 @@
"\n",
"# Run our operators\n",
"startDevito = timer.time()\n",
" \n",
"\n",
"# Apply operator\n",
"op.apply(dt=dt0)\n",
"\n",
Expand Down Expand Up @@ -293,7 +292,7 @@
"#NBVAL_IGNORE_OUTPUT\n",
"\n",
"# Get an appropriate ylimit\n",
"slice_smax = v.data[:,int(smax-smin-padding)]\n",
"slice_smax = v.data[:, int(smax-smin-padding)]\n",
"ymax = max(slice_smax) + 2\n",
"\n",
"# Plot\n",
Expand Down Expand Up @@ -477,8 +476,8 @@
"name": "stdout",
"output_type": "stream",
"text": [
"devito pde timesteps: 2000, 0.205647s runtime\n",
"call_value_bs timesteps: 5, 3.813932s runtime\n"
"devito pde timesteps: 2000, 0.019737s runtime\n",
"call_value_bs timesteps: 5, 2.596800s runtime\n"
]
},
{
Expand All @@ -501,10 +500,10 @@
"# https://aaronschlegel.me/black-scholes-formula-python.html\n",
"def call_value_bs(S, K, T, r, sigma):\n",
" N = Normal('x', 0.0, 1.0)\n",
" \n",
"\n",
" d1 = (np.log(S / K) + (r + 0.5 * sigma ** 2) * T) / (sigma * np.sqrt(T))\n",
" d2 = (np.log(S / K) + (r - 0.5 * sigma ** 2) * T) / (sigma * np.sqrt(T))\n",
" \n",
"\n",
" call = (S * cdf(N)(d1) - K * np.exp(-r * T) * cdf(N)(d2))\n",
" return call\n",
"\n",
Expand Down Expand Up @@ -563,7 +562,7 @@
{
"data": {
"text/plain": [
"[<matplotlib.lines.Line2D at 0x7f5495584640>]"
"[<matplotlib.lines.Line2D at 0x7fe2669fda50>]"
]
},
"execution_count": 9,
Expand All @@ -587,8 +586,8 @@
"#NBVAL_IGNORE_OUTPUT\n",
"\n",
"# Plot the l2 norm of the formula and our solution over time\n",
"t_range = np.linspace(dt0,1.0,50)\n",
"x_range = range(padding, smax-smin-padding*2, 1)\n",
"t_range = np.linspace(dt0, 1.0, 50)\n",
"x_range = range(padding, smax-smin-padding*2, 1)\n",
"vals = []\n",
"\n",
"for t in t_range:\n",
Expand All @@ -600,7 +599,7 @@
"\n",
" rms = np.sqrt(np.float64(l2 / len(x_range)))\n",
" vals.append(rms)\n",
" \n",
"\n",
"plt.figure(figsize=(12,10))\n",
"plt.plot(t_range, np.array(vals))"
]
Expand All @@ -609,22 +608,11 @@
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.00581731890853893"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"outputs": [],
"source": [
"#NBVAL_IGNORE_OUTPUT\n",
"\n",
"np.mean(vals)"
"assert np.isclose(np.mean(vals), 0.005885, rtol=1e-2)\n"
]
},
{
Expand Down
35 changes: 17 additions & 18 deletions examples/seismic/tutorials/03_fwi.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -173,7 +173,7 @@
"from examples.seismic import AcquisitionGeometry\n",
"\n",
"t0 = 0.\n",
"tn = 1000. \n",
"tn = 1000.\n",
"f0 = 0.010\n",
"# First, position source centrally in all dimensions, then set depth\n",
"src_coordinates = np.empty((1, 2))\n",
Expand Down Expand Up @@ -419,11 +419,10 @@
"outputs": [],
"source": [
"# Create FWI gradient kernel \n",
"from devito import Function, TimeFunction, norm\n",
"from devito import Function, norm\n",
"from examples.seismic import Receiver\n",
"\n",
"import scipy\n",
"def fwi_gradient(vp_in): \n",
"def fwi_gradient(vp_in):\n",
" # Create symbols to hold the gradient\n",
" grad = Function(name=\"grad\", grid=model.grid)\n",
" # Create placeholders for the data residual and data\n",
Expand All @@ -440,19 +439,19 @@
" for i in range(nshots):\n",
" # Update source location\n",
" geometry.src_positions[0, :] = source_locations[i, :]\n",
" \n",
"\n",
" # Generate synthetic data from true model\n",
" _, _, _ = solver.forward(vp=model.vp, rec=d_obs)\n",
" \n",
"\n",
" # Compute smooth data and full forward wavefield u0\n",
" _, u0, _ = solver.forward(vp=vp_in, save=True, rec=d_syn)\n",
" \n",
"\n",
" # Compute gradient from data residual and update objective function \n",
" compute_residual(residual, d_obs, d_syn)\n",
" \n",
"\n",
" objective += .5*norm(residual)**2\n",
" solver.gradient(rec=residual, u=u0, vp=vp_in, grad=grad)\n",
" \n",
"\n",
" return objective, grad"
]
},
Expand Down Expand Up @@ -571,11 +570,11 @@
"name": "stdout",
"output_type": "stream",
"text": [
"Objective value is 39292.605833 at iteration 1\n",
"Objective value is 24506.628229 at iteration 2\n",
"Objective value is 14386.573665 at iteration 3\n",
"Objective value is 7907.616850 at iteration 4\n",
"Objective value is 3960.106497 at iteration 5\n"
"Objective value is 39293.304688 at iteration 1\n",
"Objective value is 24506.697266 at iteration 2\n",
"Objective value is 14386.468750 at iteration 3\n",
"Objective value is 7907.354492 at iteration 4\n",
"Objective value is 3959.922607 at iteration 5\n"
]
}
],
Expand All @@ -590,19 +589,19 @@
" # Compute the functional value and gradient for the current\n",
" # model estimate\n",
" phi, direction = fwi_gradient(model0.vp)\n",
" \n",
"\n",
" # Store the history of the functional values\n",
" history[i] = phi\n",
" \n",
"\n",
" # Artificial Step length for gradient descent\n",
" # In practice this would be replaced by a Linesearch (Wolfe, ...)\n",
" # that would guarantee functional decrease Phi(m-alpha g) <= epsilon Phi(m)\n",
" # where epsilon is a minimum decrease constant\n",
" alpha = .05 / mmax(direction)\n",
" \n",
"\n",
" # Update the model estimate and enforce minimum/maximum values\n",
" update_with_box(model0.vp , alpha , direction)\n",
" \n",
"\n",
" # Log the progress made\n",
" print('Objective value is %f at iteration %d' % (phi, i+1))"
]
Expand Down
8 changes: 5 additions & 3 deletions examples/seismic/tutorials/05_staggered_acoustic.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,12 +15,14 @@
},
"outputs": [],
"source": [
"from devito import *\n",
"from devito import (SpaceDimension, Grid, TimeFunction, Constant, solve,\n",
" Operator, Eq, VectorTimeFunction, norm)\n",
"from devito.types import NODE\n",
"from examples.seismic.source import DGaussSource, TimeAxis\n",
"from examples.seismic import plot_image\n",
"import numpy as np\n",
"\n",
"from sympy import init_printing, latex\n",
"from sympy import init_printing\n",
"init_printing(use_latex='mathjax')"
]
},
Expand All @@ -32,7 +34,7 @@
},
"outputs": [],
"source": [
"# Initial grid: 1km x 1km, with spacing 100m\n",
"# Set up initial grid: 1km x 1km, with spacing 100m\n",
"extent = (2000., 2000.)\n",
"shape = (81, 81)\n",
"x = SpaceDimension(name='x', spacing=Constant(name='h_x', value=extent[0]/(shape[0]-1)))\n",
Expand Down
10,323 changes: 3,337 additions & 6,986 deletions examples/seismic/tutorials/08_snapshotting.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

Loading
Loading