Skip to content

Commit

Permalink
Answered more and add Brusselator page
Browse files Browse the repository at this point in the history
The Brusselator page as of this commit is fairly barebones and should be fleshed out with text. It can also be extended to showcase implicit solvers and boundary conditions.
  • Loading branch information
GeorgeR227 committed Dec 16, 2024
1 parent 6c51106 commit 3d755bb
Show file tree
Hide file tree
Showing 3 changed files with 132 additions and 0 deletions.
1 change: 1 addition & 0 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,7 @@ makedocs(
"Meshes" => "meshes/meshes.md",
"Vortices" => "navier_stokes/ns.md",
"Harmonics" => "harmonics/harmonics.md",
"Brusselator" => "brussel/brussel.md"
"Cahn-Hilliard" => "ch/cahn-hilliard.md",
"Klausmeier" => "klausmeier/klausmeier.md",
"CISM v2.1" => "cism/cism.md",
Expand Down
114 changes: 114 additions & 0 deletions docs/src/brussel/brussel.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,114 @@
# Brusselator

```@setup INFO
include(joinpath(Base.@__DIR__, ".." , "..", "docinfo.jl"))
info = DocInfo.Info()
```
## Dependencies
```@example DEC
using CairoMakie
import CairoMakie: wireframe, mesh, Figure, Axis
using CombinatorialSpaces
using ComponentArrays
using DiagrammaticEquations
using Decapodes
using LinearAlgebra
using MLStyle
using OrdinaryDiffEq
using GeometryBasics: Point2, Point3
Point2D = Point2{Float64}
Point3D = Point3{Float64}
```

## The Model
```@example DEC
Brusselator = @decapode begin
(U, V)::Form0
U2V::Form0
(U̇, V̇)::Form0
(α)::Constant
F::Parameter
U2V == (U .* U) .* V
U̇ == 1 + U2V - (4.4 * U) + (α * Δ(U)) + F
V̇ == (3.4 * U) - U2V + (α * Δ(V))
∂ₜ(U) == U̇
∂ₜ(V) == V̇
end
```

## The Mesh
```@example DEC
s = triangulated_grid(1,1,0.008,0.008,Point3D);
sd = EmbeddedDeltaDualComplex2D{Bool,Float64,Point2D}(s);
subdivide_duals!(sd, Circumcenter());
```

## Initial data
```@example DEC
U = map(sd[:point]) do (_,y)
22 * (y *(1-y))^(3/2)
end
V = map(sd[:point]) do (x,_)
27 * (x *(1-x))^(3/2)
end
F₁ = map(sd[:point]) do (x,y)
(x-0.3)^2 + (y-0.6)^2 ≤ (0.1)^2 ? 5.0 : 0.0
end
F₂ = zeros(nv(sd))
constants_and_parameters = (
α = 0.001,
F = t -> t ≥ 1.1 ? F₁ : F₂)
```

## Generate the Simulation
```@example DEC
sim = evalsim(Brusselator)
fₘ = sim(sd, nothing, DiagonalHodge())
u₀ = ComponentArray(U=U, V=V)
tₑ = 11.5
prob = ODEProblem(fₘ, u₀, (0, tₑ), constants_and_parameters)
soln = solve(prob, Tsit5())
```

## Visualize
```@example DEC
fig = Figure();
ax = Axis(fig[1,1])
mesh!(ax, s, color=soln(tₑ).U, colormap=:jet)
function save_dynamics(save_file_name)
time = Observable(0.0)
u = @lift(soln($time).U)
f = Figure()
ax = CairoMakie.Axis(f[1,1], title = @lift("Brusselator U Concentration at Time $($time)"))
gmsh = mesh!(ax, s, color=u, colormap=:jet,
colorrange=extrema(soln(tₑ).U))
Colorbar(f[1,2], gmsh)
timestamps = range(0, tₑ, step=1e-1)
record(f, save_file_name, timestamps; framerate = 15) do t
time[] = t
end
end
save_dynamics("brusselator.gif")
```

![Brusselator_results_flat](brusselator.gif)

```@example INFO
DocInfo.get_report(info) # hide
```

17 changes: 17 additions & 0 deletions docs/src/faq/faq.md
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,10 @@

## 1. How to incorporate scalar or vector field input data where you have a function of the embedded coordinates?

We can take a look at the [Brusselator page](../brussel/brussel.md#initial-data) which sets the values of each point on its mesh to a value as determined by some function. This can be done in a similar manner in 1D and 2D.

The Brusselator also demonstrates, with the variable `F`, how one can go about changing the function by which these values are set over time.

## 2. How to incorporate input data from a file with linear interpolation?

The Grigoriev Ice Cap model has a section where after the initial data is loaded from a TIF, the data is interpolated so that it may fit over a discrete mesh of our choosing. The link for that is [here](../grigoriev/grigoriev.md#loading-a-scientific-dataset).
Expand All @@ -10,6 +14,15 @@ The Grigoriev Ice Cap model has a section where after the initial data is loaded

## 4. How to plot a derived quantity?

Plotting in DECAPODES is commonly done with the [Makie](https://github.com/MakieOrg/Makie.jl) package in Julia. Makie allows for creating both still images, which are useful for visualizing the mesh itself and initial/final conditions, and videos, which can capture the full simulation from start to end.

- For [1D visualization](../ice_dynamics/ice_dynamics.md#visualize)

- For [2D visualization](../ice_dynamics/ice_dynamics.md#visualize-2d)

- For [3D visualization](../ice_dynamics/ice_dynamics.md#2-manifold-in-3d)


## 5. How to add artificial diffusion for 0- or 1-forms to improve stability?

## 6. How to use a Laplacian solver / multigrid?
Expand All @@ -27,3 +40,7 @@ To use multigrid methods in the Laplacian solver, you need to create a `PrimalGe
```

## 7. How to do a bunch of workflows?

A common workflow is to iterate through multiple different models as is done in the [Vorticity Model page](../navier_stokes/ns.md). A formulation is first done with a direct vorticity formulation but a quick run finds that this setup is unstable. A second formulation introduces a Laplacian solve which produces nice results.

Similar workflows may retain the same model but may iterate on the types of meshes/initial conditions used. An excellent example of this is found in the [Glacial Flow page](../ice_dynamics/ice_dynamics.md) where the model is first run in a [1D](../ice_dynamics/ice_dynamics.md#Define-a-mesh) setting and then quickly promoted to both [2D](../ice_dynamics/ice_dynamics.md#Define-our-mesh) and [3D](../ice_dynamics/ice_dynamics.md#2-Manifold-in-3D). This allows either running some dynamics in a more complicated setting, as just discussed, or allows for simplifying higher dimensional models by some sort of symmetry.

0 comments on commit 3d755bb

Please sign in to comment.