|
38 | 38 | "outputs": [], |
39 | 39 | "source": [ |
40 | 40 | "import numpy as np\n", |
| 41 | + "import numpy.testing as npt\n", |
| 42 | + "\n", |
41 | 43 | "import sksundae as sun\n", |
42 | 44 | "import matplotlib.pyplot as plt\n", |
43 | 45 | "\n", |
|
76 | 78 | "soln.y[:,1] *= 1e4 # scale the y1 values for plotting\n", |
77 | 79 | "\n", |
78 | 80 | "plt.semilogx(soln.t, soln.y)\n", |
79 | | - "plt.legend(['y0', 'y1', 'y2'])\n", |
80 | | - "plt.xlabel(r\"$t$\");\n", |
81 | | - "plt.ylabel(\"concentration\");" |
| 81 | + "plt.legend([\"y0\", \"y1\", \"y2\"])\n", |
| 82 | + "plt.xlabel(\"time, $t$\");\n", |
| 83 | + "plt.ylabel(\"concentration, $c$\");" |
82 | 84 | ] |
83 | 85 | }, |
84 | 86 | { |
|
114 | 116 | "soln.y[:,1] *= 1e4 # scale the y1 values for plotting\n", |
115 | 117 | "\n", |
116 | 118 | "plt.semilogx(soln.t, soln.y)\n", |
117 | | - "plt.legend(['y0', 'y1', 'y2'])\n", |
118 | | - "plt.xlabel(r\"$t$\");\n", |
119 | | - "plt.ylabel(\"concentration\");" |
| 119 | + "plt.legend([\"y0\", \"y1\", \"y2\"])\n", |
| 120 | + "plt.xlabel(\"time, $t$\");\n", |
| 121 | + "plt.ylabel(\"concentration, $c$\");" |
120 | 122 | ] |
121 | 123 | }, |
122 | 124 | { |
|
126 | 128 | "## Step-wise solutions\n", |
127 | 129 | "Solving step-by-step instead of across a full time span can be beneficial in some cases, especially for debugging. Therefore, a `step` method is also available in `IDA`. Before taking a step, the solver needs to know the initial conditions and time to determine the direction of integration for the following steps. Thus, before calling `step`, you should call `init_step`, as shown below. The initialization is handled automatically when using the `solve` method but must be done manually in a step-by-step approach. Since we are using the same solver instance that we initialized above, the `init_step` call will also provide a correction to `yp0` as specified by the `calc_initcond` option.\n", |
128 | 130 | "\n", |
129 | | - "Below, we run `init_step` and compare the solution `soln_0` to the initial values from the full solve, from above. Afterward, we take a step evaluated at `soln.t[10]` using the solution object from the full solve, allowing us to compare the step-by-step solution to a portion of the full solution. We only check that the solutions are within some tolerance (`1e-6`) because the solver's internal steps may differ from those in the full solution, meaning the values will be close but may not be exactly the same." |
| 131 | + "Below, we run `init_step` and compare the solution `soln_0` to the initial values from the full solve, from above. Afterward, we take a step evaluated at `soln.t[10]` using the solution object from the full solve, allowing us to compare the step-by-step solution to a portion of the full solution. We only check that the solutions are within some tolerance (`1e-7`) because the solver's internal steps may differ from those in the full solution, meaning the values will be close but may not be exactly the same." |
130 | 132 | ] |
131 | 133 | }, |
132 | 134 | { |
|
140 | 142 | "soln_0 = solver.init_step(tspan[0], y0, yp0)\n", |
141 | 143 | "print(soln_0)\n", |
142 | 144 | "\n", |
143 | | - "assert soln_0.t == soln.t[0]\n", |
144 | | - "assert np.all(soln_0.y - soln.y[0] < 1e-6)\n", |
| 145 | + "npt.assert_allclose(soln.t[0], soln_0.t)\n", |
| 146 | + "npt.assert_allclose(soln.y[0], soln_0.y, atol=1e-7)\n", |
145 | 147 | "\n", |
146 | 148 | "soln_1 = solver.step(soln.t[10])\n", |
147 | 149 | "print(soln_1)\n", |
148 | 150 | "\n", |
149 | | - "assert soln_1.t == soln.t[10]\n", |
150 | | - "assert np.all(soln_1.y - soln.y[10] < 1e-6)" |
| 151 | + "npt.assert_allclose(soln.t[10], soln_1.t)\n", |
| 152 | + "npt.assert_allclose(soln.y[10], soln_1.y, atol=1e-7)" |
151 | 153 | ] |
152 | 154 | }, |
153 | 155 | { |
|
190 | 192 | "soln.y_events[:,1] *= 1e4 # scale the y1 values for plotting\n", |
191 | 193 | "\n", |
192 | 194 | "plt.semilogx(soln.t, soln.y, '-', soln.t_events, soln.y_events, 'or')\n", |
193 | | - "plt.legend(['y0', 'y1', 'y2'])\n", |
194 | | - "plt.xlabel(r\"$t$\");\n", |
195 | | - "plt.ylabel(\"concentration\");\n", |
| 195 | + "plt.legend([\"y0\", \"y1\", \"y2\"])\n", |
| 196 | + "plt.xlabel(\"time, $t$\")\n", |
| 197 | + "plt.ylabel(\"concentration, $c$\")\n", |
196 | 198 | "\n", |
197 | 199 | "for t in soln.t_events:\n", |
198 | 200 | " plt.axvline(t, linestyle='--', color='k')" |
|
242 | 244 | "soln.y[:,1] *= 1e4 # scale the y1 values for plotting\n", |
243 | 245 | "\n", |
244 | 246 | "plt.semilogx(soln.t, soln.y)\n", |
245 | | - "plt.legend(['y0', 'y1', 'y2'])\n", |
246 | | - "plt.xlabel(r\"$t$\");\n", |
247 | | - "plt.ylabel(\"concentration\");" |
| 247 | + "plt.legend([\"y0\", \"y1\", \"y2\"])\n", |
| 248 | + "plt.xlabel(\"time, $t$\");\n", |
| 249 | + "plt.ylabel(\"concentration, $c$\");" |
248 | 250 | ] |
249 | 251 | }, |
250 | 252 | { |
|
273 | 275 | "name": "python", |
274 | 276 | "nbconvert_exporter": "python", |
275 | 277 | "pygments_lexer": "ipython3", |
276 | | - "version": "3.12.5" |
| 278 | + "version": "3.14.0" |
277 | 279 | } |
278 | 280 | }, |
279 | 281 | "nbformat": 4, |
|
0 commit comments