diff --git a/Project.toml b/Project.toml index 4a6bec5..5120d0e 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "FiniteVolumeMethod" uuid = "d4f04ab7-4f65-4d72-8a28-7087bc7f46f4" authors = ["Daniel VandenHeuvel "] -version = "0.4.0" +version = "0.4.1" [deps] ChunkSplitters = "ae650224-84b6-46f8-82ea-d812ca08434e" diff --git a/docs/src/annulus.md b/docs/src/annulus.md index c9766bf..e9ccf1e 100644 --- a/docs/src/annulus.md +++ b/docs/src/annulus.md @@ -2,13 +2,13 @@ Please ensure you have read the List of Examples section before proceeding. -We now consider another diffusion problem, except now we demonstrate how we can solve a problem on a multiply-connected domain. The problem we consider comes from http://onelab.info/wiki/Tutorial/Heat_equation_with_Dirichlet_boundary_control, namely +We now consider another diffusion problem, except now we demonstrate how we can solve a problem on a multiply-connected domain. The problem we consider comes from [http://onelab.info/wiki/Tutorial/Heat_equation_with_Dirichlet_boundary_control](http://onelab.info/wiki/Tutorial/Heat_equation_with_Dirichlet_boundary_control), namely ```math \begin{equation*} \begin{array}{rcll} \dfrac{\partial u}{\partial t} & = & \boldsymbol{\nabla}^2 u & \boldsymbol{x} \in \Omega, t>0, \\ -\boldsymbol{\nabla} u \boldsymbol{\cdot} \hat{\boldsymbol n} = 0 & \boldsymbol{x} \in \mathcal D(0, 1), t>0, \\ +\boldsymbol{\nabla} u \boldsymbol{\cdot} \hat{\boldsymbol n} &=& 0 & \boldsymbol{x} \in \mathcal D(0, 1), t>0, \\ u(x, y, t) & = & c(t) & \boldsymbol{x} \in \mathcal D(0, 0.2), t>0, \\ u(x, y, 0) & = & u_0(x, y) & \boldsymbol{x} \in \Omega, \end{array} @@ -21,7 +21,7 @@ where $\mathcal D(0, r)$ is a circle of radius $r$ centred at the origin, $\Omeg u_0(x) =10\mathrm{e}^{-25\left(\left(x + \frac12\right)^2 + \left(y + \frac12\right)^2\right)} - 10\mathrm{e}^{-45\left(\left(x - \frac12\right)^2 + \left(y - \frac12\right)^2\right)} - 5\mathrm{e}^{-50\left(\left(x + \frac{3}{10}\right)^2 + \left(y + \frac12\right)^2\right)}. ``` -To define this problem, we define the problem as we have been doing, but now we take special care to define the multiply-connected domain. In particular, we define the boundary nodes according to the specification in DelaunayTriangulation.jl (see the boundary nodes discussion here https://danielvandh.github.io/DelaunayTriangulation.jl/stable/interface/interface/). The complete code is below, where we generate the mesh, visualise the solution, and also animate it. +To define this problem, we define the problem as we have been doing, but now we take special care to define the multiply-connected domain. In particular, we define the boundary nodes according to the specification in DelaunayTriangulation.jl (see the boundary nodes discussion here [https://danielvandh.github.io/DelaunayTriangulation.jl/stable/interface/interface/](https://danielvandh.github.io/DelaunayTriangulation.jl/stable/interface/interface/)). The complete code is below, where we generate the mesh, and then visualise the solution. ```julia ## Generate the mesh. @@ -75,19 +75,6 @@ mesh!(ax, pt_mat, T_mat, color=sol.u[6], colorrange=(-10, 20), colormap=:viridis ax = Axis(fig[1, 3], width=600, height=600) mesh!(ax, pt_mat, T_mat, color=sol.u[11], colorrange=(-10, 20), colormap=:viridis) save("figures/annulus_test.png", fig) - -fig = Figure() -t_rng = LinRange(0, 2, 361) -j = Observable(1) -ax = Axis(fig[1, 1], xlabel=L"x", ylabel=L"y", - title=Makie.lift(_j -> L"t = %$(rpad(round(t_rng[_j], digits = 5), 7, '0'))", j), - titlealign=:left) -mesh!(ax, pt_mat, T_mat, color=@lift(sol(t_rng[$j])), colorrange=(-10, 20), colormap=:viridis) -record(fig, "figures/annulus_test.mp4", eachindex(t_rng)) do _j - j[] = _j -end ``` -![Annulus solution](https://github.com/DanielVandH/FiniteVolumeMethod.jl/blob/main/test/figures/annulus_test.png?raw=true) - -![Annulus animation](https://github.com/DanielVandH/FiniteVolumeMethod.jl/blob/main/test/figures/annulus_test.mp4?raw=true) \ No newline at end of file +![Annulus solution](https://github.com/DanielVandH/FiniteVolumeMethod.jl/blob/main/test/figures/annulus_test.png?raw=true) \ No newline at end of file diff --git a/docs/src/example_list.md b/docs/src/example_list.md index 6f12919..467c6d1 100644 --- a/docs/src/example_list.md +++ b/docs/src/example_list.md @@ -123,7 +123,7 @@ The purpose of this example is to show how we can solve problems on a multiply-c \begin{equation*} \begin{array}{rcll} \dfrac{\partial u}{\partial t} & = & \boldsymbol{\nabla}^2 u & \boldsymbol{x} \in \Omega, t>0, \\ -\boldsymbol{\nabla} u \boldsymbol{\cdot} \hat{\boldsymbol n} = 0 & \boldsymbol{x} \in \mathcal D(0, 1), t>0, \\ +\boldsymbol{\nabla} u \boldsymbol{\cdot} \hat{\boldsymbol n} & = & 0 & \boldsymbol{x} \in \mathcal D(0, 1), t>0, \\ u(x, y, t) & = & c(t) & \boldsymbol{x} \in \mathcal D(0, 0.2), t>0, \\ u(x, y, 0) & = & u_0(x, y) & \boldsymbol{x} \in \Omega, \end{array} @@ -144,8 +144,8 @@ The purpose of this example is to show how we can solve a steady problem. We sol \begin{equation*} \begin{array}{rcll} \boldsymbol{\nabla}^2 u(x, y) & = & 0 & 0 < x, y < \pi, \\ -u(x, 0) & = & \sinh(x) & 0 < x < \pi, -u(x, \pi) & = & -\sinh(x) & 0 < x < \pi, +u(x, 0) & = & \sinh(x) & 0 < x < \pi, \\ +u(x, \pi) & = & -\sinh(x) & 0 < x < \pi, \\ u(0, y) & = & 0 & 0 < y < \pi, \\ u(\pi, y) & = & \sinh(\pi)\cos(y) & 0 < x < \pi, \end{array} @@ -166,4 +166,4 @@ D\boldsymbol{\nabla}^2 T & = & -1, \end{equation*} ``` -where the boundary conditions could be either absorbing, reflecting, or a combination of the two, and $T$ is the mean exit time of a particle with diffusivity $D$ out of the given domain. We cover many examples, reproducing some of the work in my papers at https://iopscience.iop.org/article/10.1088/1367-2630/abe60d and https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d. \ No newline at end of file +where the boundary conditions could be either absorbing, reflecting, or a combination of the two, and $T$ is the mean exit time of a particle with diffusivity $D$ out of the given domain. We cover many examples, reproducing some of the work in my papers at [https://iopscience.iop.org/article/10.1088/1367-2630/abe60d](https://iopscience.iop.org/article/10.1088/1367-2630/abe60d) and [https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d](https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d). \ No newline at end of file diff --git a/docs/src/laplace.md b/docs/src/laplace.md index a5fa4ac..bd2664b 100644 --- a/docs/src/laplace.md +++ b/docs/src/laplace.md @@ -6,8 +6,8 @@ Now we consider Laplace's equaton. This is a steady problem, whereas all previou \begin{equation*} \begin{array}{rcll} \boldsymbol{\nabla}^2 u(x, y) & = & 0 & 0 < x, y < \pi, \\ -u(x, 0) & = & \sinh(x) & 0 < x < \pi, -u(x, \pi) & = & -\sinh(x) & 0 < x < \pi, +u(x, 0) & = & \sinh(x) & 0 < x < \pi,\\ +u(x, \pi) & = & -\sinh(x) & 0 < x < \pi, \\ u(0, y) & = & 0 & 0 < y < \pi, \\ u(\pi, y) & = & \sinh(\pi)\cos(y) & 0 < x < \pi, \end{array} @@ -32,7 +32,7 @@ we are solving problems with $\partial u/\partial t = 0$ so that \end{equation*} ``` -For the reaction-diffusion formula, this is instead given by $0 = \boldsymbol\nabla[D\boldsymbol\nabla u] + R$. Moreover, the initial condition that we provide is now actually used as an initial estimate in the nonlinear solver that computes the steady state (defined via NonlinearSolve.jl). Lastly, the final time that we provide is replaced with infinity. With this in mind, the code for solving the above problem is given below. +For the reaction-diffusion formula, this is instead given by $0 = \boldsymbol\nabla\boldsymbol\cdot[D\boldsymbol\nabla u] + R$. Moreover, the initial condition that we provide is now actually used as an initial estimate in the nonlinear solver that computes the steady state (defined via NonlinearSolve.jl). Lastly, the final time that we provide is replaced with infinity. With this in mind, the code for solving the above problem is given below. ```julia using FiniteVolumeMethod diff --git a/docs/src/mean_exit_time.md b/docs/src/mean_exit_time.md index ba84dd8..34d38d2 100644 --- a/docs/src/mean_exit_time.md +++ b/docs/src/mean_exit_time.md @@ -1,6 +1,6 @@ # Example IX: Mean exit time problems -In this example, we will investigate some simple mean exit time problems. In particular, the aim is to reproduce some of the figures from the numerical solutions in my other work on mean exit time at https://iopscience.iop.org/article/10.1088/1367-2630/abe60d and https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d. +In this example, we will investigate some simple mean exit time problems. In particular, the aim is to reproduce some of the figures from the numerical solutions in my other work on mean exit time at [https://iopscience.iop.org/article/10.1088/1367-2630/abe60d](https://iopscience.iop.org/article/10.1088/1367-2630/abe60d) and [https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d](https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d). Fistly, mean exit time problems with linear diffusion can be defined according to @@ -254,7 +254,7 @@ resize_to_layout!(fig) ## Annulus -Now we consider mean exit time on an annulus. The paper https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d also considers heterogeneous annuli, but we will not do that here. We will consider a reflecting inner boundary condition with an absorbing outer boundary condition, so that $\boldsymbol{\nabla}^2 T = -1/D$ in the annulus, $\boldsymbol{\nabla} T \vdot \hat{\boldsymbol n} = 0$ on the inner ring, and $T = 0$ on the outer ring. The annulus we consider is $R_1 < r < R_2$ with $R_1 = 2$ and $R_2 = 3$, which corresponds to an exact solution $T(x, y) = (R_2^2 - r^2)/(4D) + R_1^2\log(r/R_2)/(2D)$, where $r = \sqrt{x^2 + y^2}$. +Now we consider mean exit time on an annulus. The paper [https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d](https://iopscience.iop.org/article/10.1088/1751-8121/ac4a1d) also considers heterogeneous annuli, but we will not do that here. We will consider a reflecting inner boundary condition with an absorbing outer boundary condition, so that $\boldsymbol{\nabla}^2 T = -1/D$ in the annulus, $\boldsymbol{\nabla} T \boldsymbol{\cdot} \hat{\boldsymbol n} = 0$ on the inner ring, and $T = 0$ on the outer ring. The annulus we consider is $R_1 < r < R_2$ with $R_1 = 2$ and $R_2 = 3$, which corresponds to an exact solution $T(x, y) = (R_2^2 - r^2)/(4D) + R_1^2\log(r/R_2)/(2D)$, where $r = \sqrt{x^2 + y^2}$. ```julia using FiniteVolumeMethod @@ -370,4 +370,4 @@ Colorbar(fig[1, 2], msh) resize_to_layout!(fig) ``` -![Perturbed annulus mean exit time](https://github.com/DanielVandH/FiniteVolumeMethod.jl/blob/main/test/figures/annulus_mean_exit_time.png?raw=true) \ No newline at end of file +![Perturbed annulus mean exit time](https://github.com/DanielVandH/FiniteVolumeMethod.jl/blob/main/test/figures/perturbed_annulus_mean_exit_time.png?raw=true) \ No newline at end of file diff --git a/test/annulus_example.jl b/test/annulus_example.jl index cb309f9..948c0c4 100644 --- a/test/annulus_example.jl +++ b/test/annulus_example.jl @@ -55,17 +55,4 @@ ax = Axis(fig[1, 2], width=600, height=600) mesh!(ax, pt_mat, T_mat, color=sol.u[6], colorrange=(-10, 20), colormap=:viridis) ax = Axis(fig[1, 3], width=600, height=600) mesh!(ax, pt_mat, T_mat, color=sol.u[11], colorrange=(-10, 20), colormap=:viridis) -SAVE_FIGURE && save("figures/annulus_test.png", fig) - -if SAVE_FIGURE - fig = Figure() - t_rng = LinRange(0, 2, 361) - j = Observable(1) - ax = Axis(fig[1, 1], xlabel=L"x", ylabel=L"y", - title=Makie.lift(_j -> L"t = %$(rpad(round(t_rng[_j], digits = 5), 7, '0'))", j), - titlealign=:left) - mesh!(ax, pt_mat, T_mat, color=@lift(sol(t_rng[$j])), colorrange=(-10, 20), colormap=:viridis) - record(fig, "figures/annulus_test.mp4", eachindex(t_rng)) do _j - j[] = _j - end -end \ No newline at end of file +SAVE_FIGURE && save("figures/annulus_test.png", fig) \ No newline at end of file diff --git a/test/figures/annulus_test.mp4 b/test/figures/annulus_test.mp4 deleted file mode 100644 index 0250378..0000000 Binary files a/test/figures/annulus_test.mp4 and /dev/null differ