From cc85ab7d9821ec651bd588f80f7128a307dad3b6 Mon Sep 17 00:00:00 2001
From: ckhurana <ckhurana@lap-khurana.psfc.mit.edu>
Date: Mon, 20 Oct 2025 11:15:43 -0400
Subject: [PATCH 01/11] added settings chapter

---
 book/_toc.yml                     |   4 ++
 book/content/settings/settings.md | 101 ++++++++++++++++++++++++++++++
 2 files changed, 105 insertions(+)
 create mode 100644 book/content/settings/settings.md

diff --git a/book/_toc.yml b/book/_toc.yml
index 6cac4fa6..1c10d913 100644
--- a/book/_toc.yml
+++ b/book/_toc.yml
@@ -27,6 +27,10 @@ parts:
     - file: content/species_reactions/reactions
     - file: content/species_reactions/surface_reactions
 
+  - caption: Settings
+    chapters:
+    - file: content/settings/settings
+
   - caption: Applications
     chapters:
     - file: content/applications/task02
diff --git a/book/content/settings/settings.md b/book/content/settings/settings.md
new file mode 100644
index 00000000..44fa6e59
--- /dev/null
+++ b/book/content/settings/settings.md
@@ -0,0 +1,101 @@
+---
+jupytext:
+  formats: ipynb,md:myst
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.17.3
+kernelspec:
+  display_name: festim-workshop
+  language: python
+  name: python3
+---
+
+# Settings #
+
++++
+
+The settings of a FESTIM simulation are defined with a `festim.Settings` object. This tutorial provides information for defining required and optional settings for users to customize their simulations.
+
+Objectives:
+* Defining tolerances and solver settings
+* Setting up transient or steady-state simulations
+
++++
+
+## Defining tolerances and solver settings ##
+
+The required settings for any FESTIM simulation are the absolute and relative tolerances, while users can optionally specify the maximum number of iterations for the solver, degree order for finite element, and whether to use residual or incremental convergence criterion (for Newton solvers).
+
+We can define the tolerances using `atol` (absolute) and `rtol` (relative):
+
+```{code-cell} ipython3
+import festim as F
+
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10
+)
+```
+
+To specify the maximum number of iterations (which defaults to 30), we can use `max_iterations`:
+
+```{code-cell} ipython3
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    max_iterations=50
+)
+```
+
+To specify the degree order of the finite element (which defaults to 1), we can use `element_degree`:
+
+```{code-cell} ipython3
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    element_degree=2
+    )
+```
+
+To specify the convergence criterion, we can use `convergence_criterion` and strings for `residual` and `incremental`. For a residual-based convergence:
+
+```{code-cell} ipython3
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    convergence_criterion='residual'
+)
+```
+
+## Setting up transient or steady-state simulations ##
+
+For transient simulations, we need to define `final_time` and `stepsize`, while for steady-state problems, we simply need to set `transient` to `False`.
+
+For example, if we have an absolute and relative tolerance of `1e10` and `1e-10`, respectively, we can define the steady-state settings as:
+
+```{code-cell} ipython3
+import festim as F
+
+my_settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    transient=False,
+)
+```
+
+For a transient simulation with a run-time of 10 seconds and stepsize of 2 seconds:
+
+```{code-cell} ipython3
+my_settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    final_time=10,
+    stepsize=2
+)
+```
+
+```{note}
+FESTIM defaults the `transient` setting to `True`, while the stepsize and final time defaults to `None`.
+```

From 62954a7de4b2d9ed3fbe659e8a844f97c1ae8813 Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 14:23:10 -0500
Subject: [PATCH 02/11] added section on adaptive time stepping

---
 book/_toc.yml                     |   8 +-
 book/content/misc/settings.md     | 216 ++++++++++++++++++++++++++++++
 book/content/settings/settings.md | 101 --------------
 3 files changed, 220 insertions(+), 105 deletions(-)
 create mode 100644 book/content/misc/settings.md
 delete mode 100644 book/content/settings/settings.md

diff --git a/book/_toc.yml b/book/_toc.yml
index 8add9dc8..c1f31f6a 100644
--- a/book/_toc.yml
+++ b/book/_toc.yml
@@ -32,10 +32,6 @@ parts:
     - file: content/temperatures/temperatures_basic
     - file: content/temperatures/temperatures_advanced
     
-  - caption: Settings
-    chapters:
-    - file: content/settings/settings
-
   - caption: Post-processing
     chapters:
     - file: content/post_process/intro
@@ -44,6 +40,10 @@ parts:
     - file: content/post_process/paraview
     - file: content/post_process/pyvista
 
+  - caption: Misc
+    chapters:
+    - file: content/misc/settings
+
   - caption: Applications
     chapters:
     - file: content/applications/task02
diff --git a/book/content/misc/settings.md b/book/content/misc/settings.md
new file mode 100644
index 00000000..51eccba4
--- /dev/null
+++ b/book/content/misc/settings.md
@@ -0,0 +1,216 @@
+---
+jupytext:
+  formats: md:myst,ipynb
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.20.0
+kernelspec:
+  display_name: festim-workshop
+  language: python
+  name: python3
+---
+
+# Settings #
+
++++
+
+The settings of a FESTIM simulation are defined with a `festim.Settings` object. This tutorial provides information for defining required and optional settings for users to customize their simulations.
+
+Objectives:
+* Defining tolerances and solver settings
+* Setting up transient or steady-state simulations
+
++++
+
+## Defining tolerances and solver settings ##
+
+The required settings for any FESTIM simulation are the absolute and relative tolerances, while users can optionally specify the maximum number of iterations for the solver, degree order for finite element, and whether to use residual or incremental convergence criterion (for Newton solvers).
+
+We can define the tolerances using `atol` (absolute) and `rtol` (relative):
+
+```{code-cell} ipython3
+import festim as F
+
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10
+)
+```
+
+To specify the maximum number of iterations (which defaults to 30), we can use `max_iterations`:
+
+```{code-cell} ipython3
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    max_iterations=50
+)
+```
+
+To specify the degree order of the finite element (which defaults to 1), we can use `element_degree`:
+
+```{code-cell} ipython3
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    element_degree=2
+    )
+```
+
+To specify the convergence criterion, we can use `convergence_criterion` and strings for `residual` and `incremental`. For a residual-based convergence:
+
+```{code-cell} ipython3
+settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    convergence_criterion='residual'
+)
+```
+
+## Setting up transient or steady-state simulations ##
+
+For transient simulations, we need to define `final_time` and `stepsize`, while for steady-state problems, we simply need to set `transient` to `False`.
+
+For example, if we have an absolute and relative tolerance of `1e10` and `1e-10`, respectively, we can define the steady-state settings as:
+
+```{code-cell} ipython3
+import festim as F
+
+my_settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    transient=False,
+)
+```
+
+For a transient simulation with a run-time of 10 seconds and stepsize of 2 seconds:
+
+```{code-cell} ipython3
+my_settings = F.Settings(
+    atol=1e10,
+    rtol=1e-10,
+    final_time=10,
+    stepsize=2
+)
+```
+
+```{note}
+FESTIM defaults the `transient` setting to `True`, while the stepsize and final time defaults to `None`.
+```
+
+## Adaptive stepsize ##
+
+It is often useful to have an adaptive stepsize that grows or shrink based on the difficulty of the solution.
+
+FESTIM does that by allowing users to define their own `F.Stepsize` object.
+
+The `target_nb_iterations` sets the optimal number of Newton iterations. If more iterations are required in order to converge, this might suggest that the problem is hard to solve and a smaller stepsize is required. On the other hand, when the solver converges very quickly, this may be possible to have larger stepsizes.
+
+The parameter `growth_factor` defines by how much the stepsize is increased, and `cutback_factor` defines by how much it is shrunk.
+
+Let's demonstrate this with a simple example. Here we create an "empty" transient problem with no BCs, no source terms, nothing! The solution is $c=0$ everywhere. We do this so that the number of iterations required to "converge" is always below `target_nb_iterations`, and the stepsize is increased everytime.
+
+```{code-cell} ipython3
+import festim as F
+from dolfinx.mesh import create_unit_square
+from mpi4py import MPI
+
+fenics_mesh = create_unit_square(MPI.COMM_WORLD, 10, 10)
+festim_mesh = F.Mesh(fenics_mesh)
+
+my_model = F.HydrogenTransportProblem()
+
+material_top = F.Material(D_0=1, E_D=0)
+
+vol = F.VolumeSubdomain(id=1, material=material_top)
+
+my_model.mesh = festim_mesh
+my_model.subdomains = [vol]
+
+H = F.Species("H")
+my_model.species = [H]
+
+my_model.settings = F.Settings(atol=1e-8, rtol=1e-8, final_time=1000)
+
+my_model.temperature = 300
+
+my_model.exports = [F.TotalVolume(field=H, volume=vol)]
+```
+
+We define a `F.Stepsize` object with an initial value of 10 and some typical control parameters:
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=10,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+)
+```
+
+Let's run the adaptive time stepping model:
+
+```{code-cell} ipython3
+print("Running with adaptive time stepping...")
+my_model.initialise()
+my_model.run()
+times_fast = my_model.exports[0].t
+```
+
+Now, let's replace the stepsize by a fixed stepsize and see how they compare:
+
+```{code-cell} ipython3
+print("Running with fixed time stepping...")
+my_model.settings.stepsize = 10
+my_model.initialise()
+my_model.run()
+times_slow = my_model.exports[0].t
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+import matplotlib.pyplot as plt
+import numpy as np
+
+fig, axs = plt.subplots(2, 1, sharex=True)
+
+
+def plot(times, label=None):
+    steps = np.arange(len(times))
+    dt = np.diff(times)
+
+    axs[0].plot(steps[1:], dt, label=label)
+    (l,) = axs[1].plot(steps, times, label=label)
+    axs[1].scatter(steps[-1], times[-1], color=l.get_color())
+    axs[1].annotate(
+        "Final time",
+        (steps[-1], times[-1]),
+        textcoords="offset points",
+        xytext=(-10, -5),
+        ha="right",
+        color=l.get_color(),
+    )
+
+
+plot(times_fast, label="Adaptive time stepping")
+plot(times_slow, label="Fixed time stepping")
+
+axs[0].set_ylabel(r"Step size $\Delta t$")
+axs[0].set_ylim(bottom=0)
+axs[1].set_ylabel("Time $t$")
+axs[0].legend(loc="upper right")
+plt.xlabel("Timestep")
+
+plt.show()
+```
+
+As expected, the stepsize is growing at each time step, meaning the final time is reached in just above 20 timesteps. Whereas for the fixed time stepping, it takes 100 iterations.
+
++++
+
+
+
+## Custom PETSC parameters ##
diff --git a/book/content/settings/settings.md b/book/content/settings/settings.md
deleted file mode 100644
index 44fa6e59..00000000
--- a/book/content/settings/settings.md
+++ /dev/null
@@ -1,101 +0,0 @@
----
-jupytext:
-  formats: ipynb,md:myst
-  text_representation:
-    extension: .md
-    format_name: myst
-    format_version: 0.13
-    jupytext_version: 1.17.3
-kernelspec:
-  display_name: festim-workshop
-  language: python
-  name: python3
----
-
-# Settings #
-
-+++
-
-The settings of a FESTIM simulation are defined with a `festim.Settings` object. This tutorial provides information for defining required and optional settings for users to customize their simulations.
-
-Objectives:
-* Defining tolerances and solver settings
-* Setting up transient or steady-state simulations
-
-+++
-
-## Defining tolerances and solver settings ##
-
-The required settings for any FESTIM simulation are the absolute and relative tolerances, while users can optionally specify the maximum number of iterations for the solver, degree order for finite element, and whether to use residual or incremental convergence criterion (for Newton solvers).
-
-We can define the tolerances using `atol` (absolute) and `rtol` (relative):
-
-```{code-cell} ipython3
-import festim as F
-
-settings = F.Settings(
-    atol=1e10,
-    rtol=1e-10
-)
-```
-
-To specify the maximum number of iterations (which defaults to 30), we can use `max_iterations`:
-
-```{code-cell} ipython3
-settings = F.Settings(
-    atol=1e10,
-    rtol=1e-10,
-    max_iterations=50
-)
-```
-
-To specify the degree order of the finite element (which defaults to 1), we can use `element_degree`:
-
-```{code-cell} ipython3
-settings = F.Settings(
-    atol=1e10,
-    rtol=1e-10,
-    element_degree=2
-    )
-```
-
-To specify the convergence criterion, we can use `convergence_criterion` and strings for `residual` and `incremental`. For a residual-based convergence:
-
-```{code-cell} ipython3
-settings = F.Settings(
-    atol=1e10,
-    rtol=1e-10,
-    convergence_criterion='residual'
-)
-```
-
-## Setting up transient or steady-state simulations ##
-
-For transient simulations, we need to define `final_time` and `stepsize`, while for steady-state problems, we simply need to set `transient` to `False`.
-
-For example, if we have an absolute and relative tolerance of `1e10` and `1e-10`, respectively, we can define the steady-state settings as:
-
-```{code-cell} ipython3
-import festim as F
-
-my_settings = F.Settings(
-    atol=1e10,
-    rtol=1e-10,
-    transient=False,
-)
-```
-
-For a transient simulation with a run-time of 10 seconds and stepsize of 2 seconds:
-
-```{code-cell} ipython3
-my_settings = F.Settings(
-    atol=1e10,
-    rtol=1e-10,
-    final_time=10,
-    stepsize=2
-)
-```
-
-```{note}
-FESTIM defaults the `transient` setting to `True`, while the stepsize and final time defaults to `None`.
-```

From d5df5b24fd09edae94ea17a0cd9be9ec6904d2bd Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 14:28:13 -0500
Subject: [PATCH 03/11] capped stepsize

---
 book/content/misc/settings.md | 53 +++++++++++++++++++++++++++++------
 1 file changed, 44 insertions(+), 9 deletions(-)

diff --git a/book/content/misc/settings.md b/book/content/misc/settings.md
index 51eccba4..eaaeec90 100644
--- a/book/content/misc/settings.md
+++ b/book/content/misc/settings.md
@@ -178,21 +178,22 @@ import numpy as np
 fig, axs = plt.subplots(2, 1, sharex=True)
 
 
-def plot(times, label=None):
+def plot(times, label=None, annotate=True):
     steps = np.arange(len(times))
     dt = np.diff(times)
 
     axs[0].plot(steps[1:], dt, label=label)
     (l,) = axs[1].plot(steps, times, label=label)
     axs[1].scatter(steps[-1], times[-1], color=l.get_color())
-    axs[1].annotate(
-        "Final time",
-        (steps[-1], times[-1]),
-        textcoords="offset points",
-        xytext=(-10, -5),
-        ha="right",
-        color=l.get_color(),
-    )
+    if annotate:
+        axs[1].annotate(
+            "Final time",
+            (steps[-1], times[-1]),
+            textcoords="offset points",
+            xytext=(-10, -5),
+            ha="right",
+            color=l.get_color(),
+        )
 
 
 plot(times_fast, label="Adaptive time stepping")
@@ -211,6 +212,40 @@ As expected, the stepsize is growing at each time step, meaning the final time i
 
 +++
 
+Stepsize can be capped by setting the parameter `max_stepsize`:
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=10,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+    max_stepsize=40, # maximum step size
+)
+
+my_model.initialise()
+my_model.run()
+
+capped_times = my_model.exports[0].t
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+fig, axs = plt.subplots(2, 1, sharex=True)
+plot(times_fast, label="Adaptive time stepping", annotate=False)
+plot(times_slow, label="Fixed time stepping", annotate=False)
+plot(capped_times, label="Adaptive with max", annotate=False)
+
+axs[0].set_ylabel(r"Step size $\Delta t$")
+axs[0].set_ylim(bottom=0)
+axs[1].set_ylabel("Time $t$")
+axs[0].legend(loc="upper right")
+plt.xlabel("Timestep")
+
+plt.show()
+```
+
 
 
 ## Custom PETSC parameters ##

From 15db7f40e8f8aab6a55e862fcd44e339bec1e98d Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 14:29:41 -0500
Subject: [PATCH 04/11] separate section for stepsize

---
 book/_toc.yml                 |   1 +
 book/content/misc/settings.md | 148 -----------------------------
 book/content/misc/stepsize.md | 172 ++++++++++++++++++++++++++++++++++
 3 files changed, 173 insertions(+), 148 deletions(-)
 create mode 100644 book/content/misc/stepsize.md

diff --git a/book/_toc.yml b/book/_toc.yml
index c1f31f6a..37e5a02d 100644
--- a/book/_toc.yml
+++ b/book/_toc.yml
@@ -43,6 +43,7 @@ parts:
   - caption: Misc
     chapters:
     - file: content/misc/settings
+    - file: content/misc/stepsize
 
   - caption: Applications
     chapters:
diff --git a/book/content/misc/settings.md b/book/content/misc/settings.md
index eaaeec90..c8cc5948 100644
--- a/book/content/misc/settings.md
+++ b/book/content/misc/settings.md
@@ -100,152 +100,4 @@ my_settings = F.Settings(
 FESTIM defaults the `transient` setting to `True`, while the stepsize and final time defaults to `None`.
 ```
 
-## Adaptive stepsize ##
-
-It is often useful to have an adaptive stepsize that grows or shrink based on the difficulty of the solution.
-
-FESTIM does that by allowing users to define their own `F.Stepsize` object.
-
-The `target_nb_iterations` sets the optimal number of Newton iterations. If more iterations are required in order to converge, this might suggest that the problem is hard to solve and a smaller stepsize is required. On the other hand, when the solver converges very quickly, this may be possible to have larger stepsizes.
-
-The parameter `growth_factor` defines by how much the stepsize is increased, and `cutback_factor` defines by how much it is shrunk.
-
-Let's demonstrate this with a simple example. Here we create an "empty" transient problem with no BCs, no source terms, nothing! The solution is $c=0$ everywhere. We do this so that the number of iterations required to "converge" is always below `target_nb_iterations`, and the stepsize is increased everytime.
-
-```{code-cell} ipython3
-import festim as F
-from dolfinx.mesh import create_unit_square
-from mpi4py import MPI
-
-fenics_mesh = create_unit_square(MPI.COMM_WORLD, 10, 10)
-festim_mesh = F.Mesh(fenics_mesh)
-
-my_model = F.HydrogenTransportProblem()
-
-material_top = F.Material(D_0=1, E_D=0)
-
-vol = F.VolumeSubdomain(id=1, material=material_top)
-
-my_model.mesh = festim_mesh
-my_model.subdomains = [vol]
-
-H = F.Species("H")
-my_model.species = [H]
-
-my_model.settings = F.Settings(atol=1e-8, rtol=1e-8, final_time=1000)
-
-my_model.temperature = 300
-
-my_model.exports = [F.TotalVolume(field=H, volume=vol)]
-```
-
-We define a `F.Stepsize` object with an initial value of 10 and some typical control parameters:
-
-```{code-cell} ipython3
-my_model.settings.stepsize = F.Stepsize(
-    initial_value=10,
-    growth_factor=1.1,  # grow by 10%
-    cutback_factor=0.9,  # shrink by 10%
-    target_nb_iterations=4,  # target number of iterations per time step
-)
-```
-
-Let's run the adaptive time stepping model:
-
-```{code-cell} ipython3
-print("Running with adaptive time stepping...")
-my_model.initialise()
-my_model.run()
-times_fast = my_model.exports[0].t
-```
-
-Now, let's replace the stepsize by a fixed stepsize and see how they compare:
-
-```{code-cell} ipython3
-print("Running with fixed time stepping...")
-my_model.settings.stepsize = 10
-my_model.initialise()
-my_model.run()
-times_slow = my_model.exports[0].t
-```
-
-```{code-cell} ipython3
-:tags: [hide-input]
-
-import matplotlib.pyplot as plt
-import numpy as np
-
-fig, axs = plt.subplots(2, 1, sharex=True)
-
-
-def plot(times, label=None, annotate=True):
-    steps = np.arange(len(times))
-    dt = np.diff(times)
-
-    axs[0].plot(steps[1:], dt, label=label)
-    (l,) = axs[1].plot(steps, times, label=label)
-    axs[1].scatter(steps[-1], times[-1], color=l.get_color())
-    if annotate:
-        axs[1].annotate(
-            "Final time",
-            (steps[-1], times[-1]),
-            textcoords="offset points",
-            xytext=(-10, -5),
-            ha="right",
-            color=l.get_color(),
-        )
-
-
-plot(times_fast, label="Adaptive time stepping")
-plot(times_slow, label="Fixed time stepping")
-
-axs[0].set_ylabel(r"Step size $\Delta t$")
-axs[0].set_ylim(bottom=0)
-axs[1].set_ylabel("Time $t$")
-axs[0].legend(loc="upper right")
-plt.xlabel("Timestep")
-
-plt.show()
-```
-
-As expected, the stepsize is growing at each time step, meaning the final time is reached in just above 20 timesteps. Whereas for the fixed time stepping, it takes 100 iterations.
-
-+++
-
-Stepsize can be capped by setting the parameter `max_stepsize`:
-
-```{code-cell} ipython3
-my_model.settings.stepsize = F.Stepsize(
-    initial_value=10,
-    growth_factor=1.1,  # grow by 10%
-    cutback_factor=0.9,  # shrink by 10%
-    target_nb_iterations=4,  # target number of iterations per time step
-    max_stepsize=40, # maximum step size
-)
-
-my_model.initialise()
-my_model.run()
-
-capped_times = my_model.exports[0].t
-```
-
-```{code-cell} ipython3
-:tags: [hide-input]
-
-fig, axs = plt.subplots(2, 1, sharex=True)
-plot(times_fast, label="Adaptive time stepping", annotate=False)
-plot(times_slow, label="Fixed time stepping", annotate=False)
-plot(capped_times, label="Adaptive with max", annotate=False)
-
-axs[0].set_ylabel(r"Step size $\Delta t$")
-axs[0].set_ylim(bottom=0)
-axs[1].set_ylabel("Time $t$")
-axs[0].legend(loc="upper right")
-plt.xlabel("Timestep")
-
-plt.show()
-```
-
-
-
 ## Custom PETSC parameters ##
diff --git a/book/content/misc/stepsize.md b/book/content/misc/stepsize.md
new file mode 100644
index 00000000..02613842
--- /dev/null
+++ b/book/content/misc/stepsize.md
@@ -0,0 +1,172 @@
+---
+jupytext:
+  formats: md:myst,ipynb
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.20.0
+kernelspec:
+  display_name: festim-workshop
+  language: python
+  name: python3
+---
+
+# Stepsize #
+
++++
+
+The settings of a FESTIM simulation are defined with a `festim.Settings` object. This tutorial provides information for defining required and optional settings for users to customize their simulations.
+
+Objectives:
+* Defining tolerances and solver settings
+* Setting up transient or steady-state simulations
+
++++
+
+## Adaptive stepsize ##
+
+It is often useful to have an adaptive stepsize that grows or shrink based on the difficulty of the solution.
+
+FESTIM does that by allowing users to define their own `F.Stepsize` object.
+
+The `target_nb_iterations` sets the optimal number of Newton iterations. If more iterations are required in order to converge, this might suggest that the problem is hard to solve and a smaller stepsize is required. On the other hand, when the solver converges very quickly, this may be possible to have larger stepsizes.
+
+The parameter `growth_factor` defines by how much the stepsize is increased, and `cutback_factor` defines by how much it is shrunk.
+
+Let's demonstrate this with a simple example. Here we create an "empty" transient problem with no BCs, no source terms, nothing! The solution is $c=0$ everywhere. We do this so that the number of iterations required to "converge" is always below `target_nb_iterations`, and the stepsize is increased everytime.
+
+```{code-cell} ipython3
+import festim as F
+from dolfinx.mesh import create_unit_square
+from mpi4py import MPI
+
+fenics_mesh = create_unit_square(MPI.COMM_WORLD, 10, 10)
+festim_mesh = F.Mesh(fenics_mesh)
+
+my_model = F.HydrogenTransportProblem()
+
+material_top = F.Material(D_0=1, E_D=0)
+
+vol = F.VolumeSubdomain(id=1, material=material_top)
+
+my_model.mesh = festim_mesh
+my_model.subdomains = [vol]
+
+H = F.Species("H")
+my_model.species = [H]
+
+my_model.settings = F.Settings(atol=1e-8, rtol=1e-8, final_time=1000)
+
+my_model.temperature = 300
+
+my_model.exports = [F.TotalVolume(field=H, volume=vol)]
+```
+
+We define a `F.Stepsize` object with an initial value of 10 and some typical control parameters:
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=10,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+)
+```
+
+Let's run the adaptive time stepping model:
+
+```{code-cell} ipython3
+print("Running with adaptive time stepping...")
+my_model.initialise()
+my_model.run()
+times_fast = my_model.exports[0].t
+```
+
+Now, let's replace the stepsize by a fixed stepsize and see how they compare:
+
+```{code-cell} ipython3
+print("Running with fixed time stepping...")
+my_model.settings.stepsize = 10
+my_model.initialise()
+my_model.run()
+times_slow = my_model.exports[0].t
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+import matplotlib.pyplot as plt
+import numpy as np
+
+fig, axs = plt.subplots(2, 1, sharex=True)
+
+
+def plot(times, label=None, annotate=True):
+    steps = np.arange(len(times))
+    dt = np.diff(times)
+
+    axs[0].plot(steps[1:], dt, label=label)
+    (l,) = axs[1].plot(steps, times, label=label)
+    axs[1].scatter(steps[-1], times[-1], color=l.get_color())
+    if annotate:
+        axs[1].annotate(
+            "Final time",
+            (steps[-1], times[-1]),
+            textcoords="offset points",
+            xytext=(-10, -5),
+            ha="right",
+            color=l.get_color(),
+        )
+
+
+plot(times_fast, label="Adaptive time stepping")
+plot(times_slow, label="Fixed time stepping")
+
+axs[0].set_ylabel(r"Step size $\Delta t$")
+axs[0].set_ylim(bottom=0)
+axs[1].set_ylabel("Time $t$")
+axs[0].legend(loc="upper right")
+plt.xlabel("Timestep")
+
+plt.show()
+```
+
+As expected, the stepsize is growing at each time step, meaning the final time is reached in just above 20 timesteps. Whereas for the fixed time stepping, it takes 100 iterations.
+
++++
+
+Stepsize can be capped by setting the parameter `max_stepsize`:
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=10,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+    max_stepsize=40, # maximum step size
+)
+
+my_model.initialise()
+my_model.run()
+
+capped_times = my_model.exports[0].t
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+fig, axs = plt.subplots(2, 1, sharex=True)
+plot(times_fast, label="Adaptive time stepping", annotate=False)
+plot(times_slow, label="Fixed time stepping", annotate=False)
+plot(capped_times, label="Adaptive with max", annotate=False)
+
+axs[0].set_ylabel(r"Step size $\Delta t$")
+axs[0].set_ylim(bottom=0)
+axs[1].set_ylabel("Time $t$")
+axs[0].legend(loc="upper right")
+plt.xlabel("Timestep")
+
+plt.show()
+```
+

From e1bfc613ed7b601df0f640cb8fd94889fd9c2ae8 Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 14:32:43 -0500
Subject: [PATCH 05/11] changed intro

---
 book/content/misc/stepsize.md | 11 +++++++----
 1 file changed, 7 insertions(+), 4 deletions(-)

diff --git a/book/content/misc/stepsize.md b/book/content/misc/stepsize.md
index 02613842..ac92612f 100644
--- a/book/content/misc/stepsize.md
+++ b/book/content/misc/stepsize.md
@@ -16,14 +16,17 @@ kernelspec:
 
 +++
 
-The settings of a FESTIM simulation are defined with a `festim.Settings` object. This tutorial provides information for defining required and optional settings for users to customize their simulations.
-
 Objectives:
-* Defining tolerances and solver settings
-* Setting up transient or steady-state simulations
+* Ensure the simulation hits certain time milestones
+* Set a stepsize for your simulation 
+* Accelerate your simulation with adaptive time stepping
 
 +++
 
+## Milestones ##
+
+
+
 ## Adaptive stepsize ##
 
 It is often useful to have an adaptive stepsize that grows or shrink based on the difficulty of the solution.

From 01a1093d963d8b3087699e9eec62379322f7cd55 Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 15:13:27 -0500
Subject: [PATCH 06/11] added milestone section

---
 book/content/misc/stepsize.md | 165 ++++++++++++++++++++++++++++++++--
 1 file changed, 160 insertions(+), 5 deletions(-)

diff --git a/book/content/misc/stepsize.md b/book/content/misc/stepsize.md
index ac92612f..75f288dd 100644
--- a/book/content/misc/stepsize.md
+++ b/book/content/misc/stepsize.md
@@ -17,16 +17,12 @@ kernelspec:
 +++
 
 Objectives:
-* Ensure the simulation hits certain time milestones
 * Set a stepsize for your simulation 
 * Accelerate your simulation with adaptive time stepping
+* Ensure the simulation hits certain time milestones
 
 +++
 
-## Milestones ##
-
-
-
 ## Adaptive stepsize ##
 
 It is often useful to have an adaptive stepsize that grows or shrink based on the difficulty of the solution.
@@ -173,3 +169,162 @@ plt.xlabel("Timestep")
 plt.show()
 ```
 
+## Milestones ##
+
+Now that we know how to grow the stepsize to accelerate simulations, we have a new problem: what if we want the simulation to pass by a specific point in time but the timestep could be so big it completely misses it?
+
+Let's illustrate this by setting up a problem with a particle source only turning on only between 100 and 105 s, and $c=0$ on the boundary.
+
+```{code-cell} ipython3
+import festim as F
+from dolfinx.mesh import create_unit_square
+from mpi4py import MPI
+
+fenics_mesh = create_unit_square(MPI.COMM_WORLD, 10, 10)
+festim_mesh = F.Mesh(fenics_mesh)
+
+my_model = F.HydrogenTransportProblem()
+
+material_top = F.Material(D_0=0.1, E_D=0)
+
+vol = F.VolumeSubdomain(id=1, material=material_top)
+boundary = F.SurfaceSubdomain(id=2, locator=lambda x: np.full_like(x[0], True))
+
+my_model.mesh = festim_mesh
+my_model.subdomains = [vol, boundary]
+
+H = F.Species("H")
+my_model.species = [H]
+
+my_model.settings = F.Settings(atol=1e-8, rtol=1e-8, final_time=200)
+
+my_model.temperature = 300
+
+my_model.boundary_conditions = [
+    F.FixedConcentrationBC(value=0, species=H, subdomain=boundary)
+]
+```
+
+We set a time-dependent particle source term:
+
+```{code-cell} ipython3
+
+source_start = 100
+source_end = 105
+my_model.sources = [
+    F.ParticleSource(value=lambda t: 1 if t <= source_end and t >= source_start else 0, species=H, volume=vol)
+]
+```
+
+We track the total quantity of H by adding a `TotalVolume` derived quantity:
+
+```{code-cell} ipython3
+my_model.exports = [F.TotalVolume(field=H, volume=vol)]
+```
+
+Then we set an adaptive timestep with a fairly large initial value:
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=20,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+)
+
+my_model.initialise()
+my_model.run()
+```
+
+By plotting the values of `TotalVolume` we see that:
+* the timesteps kind of jump over the time period of interest
+* the value is zero all the time, which is WRONG!
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+import matplotlib.pyplot as plt
+
+plt.plot(my_model.exports[0].t, my_model.exports[0].data, marker="o")
+plt.axvline(
+    source_start, color="tab:green", alpha=0.5, linestyle="--", label="Source start"
+)
+plt.axvline(source_end, color="tab:red", alpha=0.5, linestyle="--", label="Source end")
+plt.legend()
+plt.xlabel("Time $t$")
+plt.ylabel("Total amount of H")
+plt.ylim(bottom=-0.001)
+plt.show()
+```
+
+We can try and solve it first using `milestones`. Here we set a list of two milestones at the beginning and at the end of the source period:
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=20,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+    milestones=[source_start, source_end]
+)
+
+my_model.initialise()
+my_model.run()
+```
+
+The result is already better, the solution is not zero all the time. But still, the solution looks a bit whacky... This is because, while the stepsize is modified (truncated) to hit the milestones, it is still fairly large during the time period of interest.
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+plt.plot(my_model.exports[0].t, my_model.exports[0].data, marker="o")
+plt.axvline(
+    source_start, color="tab:green", alpha=0.5, linestyle="--", label="Source start"
+)
+plt.axvline(source_end, color="tab:red", alpha=0.5, linestyle="--", label="Source end")
+plt.legend()
+plt.xlabel("Time $t$")
+plt.ylabel("Total amount of H")
+plt.ylim(bottom=-0.001)
+plt.show()
+```
+
+To improve it, let's set the `max_stepsize` argument. Here we want the stepsize to be capped at `0.5` when the time is bewteen `source_start - 5` and `source_end + 5`, otherwise, no limit (`None`).
+
+We also modify the first milestone to hit just before the source is turned on.
+
+We pass a `lambda` funtion to `max_stepsize` which is a function of `t`.
+
+```{code-cell} ipython3
+my_model.settings.stepsize = F.Stepsize(
+    initial_value=20,
+    growth_factor=1.1,  # grow by 10%
+    cutback_factor=0.9,  # shrink by 10%
+    target_nb_iterations=4,  # target number of iterations per time step
+    milestones=[source_start - 5, source_end],
+    max_stepsize=lambda t: 0.5 if source_start - 5 <= t <= source_end + 5 else None
+)
+
+my_model.initialise()
+my_model.run()
+```
+
+Now the solution looks much smoother! 🎉
+
+We can also see that the stepsize starts increasing again after ~115 s.
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+plt.plot(my_model.exports[0].t, my_model.exports[0].data, marker="o")
+plt.axvline(
+    source_start, color="tab:green", alpha=0.5, linestyle="--", label="Source start"
+)
+plt.axvline(source_end, color="tab:red", alpha=0.5, linestyle="--", label="Source end")
+plt.legend()
+plt.xlabel("Time $t$")
+plt.ylabel("Total amount of H")
+plt.ylim(bottom=-0.001)
+plt.xlim(source_start - 10, source_end + 20)
+plt.show()
+```

From 78ae91551a4ba5f53c9e4510ada279b91a03d1d5 Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 15:32:28 -0500
Subject: [PATCH 07/11] added from scratch

---
 book/content/misc/festim_from_scratch.md | 24 ++++++++++++++++++++++++
 1 file changed, 24 insertions(+)
 create mode 100644 book/content/misc/festim_from_scratch.md

diff --git a/book/content/misc/festim_from_scratch.md b/book/content/misc/festim_from_scratch.md
new file mode 100644
index 00000000..672d7328
--- /dev/null
+++ b/book/content/misc/festim_from_scratch.md
@@ -0,0 +1,24 @@
+---
+jupytext:
+  formats: md:myst,ipynb
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.20.0
+kernelspec:
+  display_name: festim-workshop
+  language: python
+  name: python3
+---
+
+# Build FESTIM from scratch #
+
++++
+
+Objectives:
+* Understand the FESTIM foundations
+* Be better armed for contributing to FESTIM
+
++++
+

From 620478d425647ab8c5d3ddff711851ef4b5b713e Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 15:35:13 -0500
Subject: [PATCH 08/11] from scratch in toc

---
 book/_toc.yml | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/book/_toc.yml b/book/_toc.yml
index 37e5a02d..4a38f0ed 100644
--- a/book/_toc.yml
+++ b/book/_toc.yml
@@ -44,6 +44,8 @@ parts:
     chapters:
     - file: content/misc/settings
     - file: content/misc/stepsize
+    - file: content/misc/festim_from_scratch
+
 
   - caption: Applications
     chapters:

From e805e17b26b564fdae020bd7340616d3168a0a1d Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 17:22:17 -0500
Subject: [PATCH 09/11] plot source

---
 book/content/misc/stepsize.md | 32 ++++++++++++++++++++++++++++++--
 1 file changed, 30 insertions(+), 2 deletions(-)

diff --git a/book/content/misc/stepsize.md b/book/content/misc/stepsize.md
index 75f288dd..15aec137 100644
--- a/book/content/misc/stepsize.md
+++ b/book/content/misc/stepsize.md
@@ -207,15 +207,43 @@ my_model.boundary_conditions = [
 
 We set a time-dependent particle source term:
 
-```{code-cell} ipython3
+```{math}
+\begin{cases}
+1, & \text{for } 100\leq t\leq 105\\
+0, & \text{otherwise }
+\end{case}
+```
 
+```{code-cell} ipython3
 source_start = 100
 source_end = 105
+
+def source_value(t):
+    if t <= source_end and t >= source_start:
+        return 1
+    else:
+        return 0
+
 my_model.sources = [
-    F.ParticleSource(value=lambda t: 1 if t <= source_end and t >= source_start else 0, species=H, volume=vol)
+    F.ParticleSource(value=source_value, species=H, volume=vol)
 ]
 ```
 
+```{code-cell} ipython3
+:tags: [hide-input]
+
+import matplotlib.pyplot as plt
+
+times = np.linspace(0, 200, 1000)
+source_values = [source_value(t) for t in times]
+plt.plot(times, source_values)
+plt.fill_between(times, source_values, alpha=0.3)
+
+plt.xlabel("Time")
+plt.ylabel("Source value")
+plt.show()
+```
+
 We track the total quantity of H by adding a `TotalVolume` derived quantity:
 
 ```{code-cell} ipython3

From 4a7aef3f400558f112ddbd4a270e3e4749380adf Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 21:44:58 -0500
Subject: [PATCH 10/11] double dollar instead

---
 book/content/misc/stepsize.md | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/book/content/misc/stepsize.md b/book/content/misc/stepsize.md
index 15aec137..628da802 100644
--- a/book/content/misc/stepsize.md
+++ b/book/content/misc/stepsize.md
@@ -207,12 +207,12 @@ my_model.boundary_conditions = [
 
 We set a time-dependent particle source term:
 
-```{math}
+$$
 \begin{cases}
 1, & \text{for } 100\leq t\leq 105\\
 0, & \text{otherwise }
 \end{case}
-```
+$$
 
 ```{code-cell} ipython3
 source_start = 100

From f2c697162daeed55fef14b48672cfb067a7ec665 Mon Sep 17 00:00:00 2001
From: RemDelaporteMathurin <rdelaportemathurin@gmail.com>
Date: Mon, 23 Feb 2026 21:48:53 -0500
Subject: [PATCH 11/11] fixed typo

---
 book/content/misc/stepsize.md | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/book/content/misc/stepsize.md b/book/content/misc/stepsize.md
index 628da802..469c4d8c 100644
--- a/book/content/misc/stepsize.md
+++ b/book/content/misc/stepsize.md
@@ -205,13 +205,13 @@ my_model.boundary_conditions = [
 ]
 ```
 
-We set a time-dependent particle source term:
+We set a time-dependent particle source term $S$:
 
 $$
-\begin{cases}
+S = \begin{cases}
 1, & \text{for } 100\leq t\leq 105\\
 0, & \text{otherwise }
-\end{case}
+\end{cases}
 $$
 
 ```{code-cell} ipython3