Skip to content

Commit

Permalink
DOC: add harmonics explainer (#264)
Browse files Browse the repository at this point in the history
* DOC: add harmonics explainer

* Use Cairo, link to eigen method

* Temp unlink gif for test build

* add copy to harmonics page

* DOC: remove draft status

* Finishing touches

* Add back paragraph

---------

Co-authored-by: Luke Morris <[email protected]>
  • Loading branch information
jpfairbanks and lukem12345 authored Sep 4, 2024
1 parent c19742c commit e9779ad
Show file tree
Hide file tree
Showing 4 changed files with 103 additions and 1 deletion.
1 change: 1 addition & 0 deletions docs/Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0"
Catlab = "134e5e36-593f-5add-ad60-77f754baafbe"
CombinatorialSpaces = "b1c52339-7909-45ad-8b6a-6e388f7c67f2"
ComponentArrays = "b0b7db55-cfe3-40fc-9ded-d10e2dbeff66"
CoordRefSystems = "b46f11dc-f210-4604-bfba-323c1ec968cb"
Decapodes = "679ab3ea-c928-4fe6-8d59-fd451142d391"
DiagrammaticEquations = "6f00c28b-6bed-4403-80fa-30e0dc12f317"
Distributed = "8ba89e20-285c-5b6f-9357-94700520ee1b"
Expand Down
2 changes: 2 additions & 0 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ makedocs(
sitename = "Decapodes.jl",
doctest = false,
checkdocs = :none,
draft = false;
pagesonly = true,
linkcheck = true,
linkcheck_ignore = [r"agupubs\.onlinelibrary\.wiley\.com", # This gives a 403 Forbidden
Expand All @@ -49,6 +50,7 @@ makedocs(
"Overview" => "overview/overview.md",
"Equations" => "equations/equations.md",
"Vortices" => "navier_stokes/ns.md",
"Harmonics" => "harmonics/harmonics.md",
"Cahn-Hilliard" => "ch/cahn-hilliard.md",
"Klausmeier" => "klausmeier/klausmeier.md",
"CISM v2.1" => "cism/cism.md",
Expand Down
100 changes: 100 additions & 0 deletions docs/src/harmonics/harmonics.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,100 @@
# Harmonics of the Sphere

This page shows how to use Decapodes tooling to explore the harmonics of a discrete manifold. This isn't using any Decapodes specific code, but it
is emblematic of a more advanced analysis you might want to do on your Decapode.

In this case we are trying to visualize the roots of the Laplacian on a discrete manifold.

Load the dependencies

```@example Harmonics
# Meshing:
using CombinatorialSpaces
using CoordRefSystems
using GeometryBasics: Point3
const Point3D = Point3{Float64};
# Visualization:
using CairoMakie
# Simulation:
using LinearAlgebra
```

Load the mesh

```@example Harmonics
const RADIUS = 1.0
s = loadmesh(Icosphere(3, RADIUS));
sd = EmbeddedDeltaDualComplex2D{Bool, Float64, Point3D}(s);
subdivide_duals!(sd, Barycenter());
```

Compute the Laplacian eigenvectors using [LinearAlgebra.eigen](https://docs.julialang.org/en/v1/stdlib/LinearAlgebra/#LinearAlgebra.eigen). This requires making the sparse Laplacian matrix dense with `collect`. Alternatively, use [Arpack.jl](https://arpack.julialinearalgebra.org/stable/).

```@example Harmonics
Δ0 = -Δ(0,sd)
λ = eigen(collect(Δ0))
```

Let's check that our eigenvalues satisfy the right equation.
The first eigenvector should be the kernel of the laplacian.
So the following norm should be close to 0.

```@example Harmonics
q1 = λ.vectors[:,1]
norm(Δ0 *q1)
```

The first eigenvector is boring to visualize, because it is constant.
So we will make some plots of the second eigenvector.
If you run this on the desktop, you can use GLMakie and get an interactive plot to explore. We will just draw two angles.

```@example Harmonics
q = λ.vectors[:,2]
fig = Figure()
Label(fig[1, 1, Top()], "Default Angle", padding = (0, 0, 5, 0))
ax = LScene(fig[1,1], scenekw=(lights=[],))
msh = CairoMakie.mesh!(ax, s, color=q)
Colorbar(fig[1,2], msh, size=32)
# Second Angle
Label(fig[2, 1, Top()], "Bottom Angle", padding = (0, 0, 5, 0))
ax = LScene(fig[2,1], scenekw=(lights=[],))
update_cam!(ax.scene, Vec3f(-1/2,-1/2,1.0/2), Vec3f(1,1,1), Vec3f(0, 0, 1))
msh = CairoMakie.mesh!(ax, s, color=q)
Colorbar(fig[2,2], msh, size=32)
fig
```

```@example Harmonics
q = λ.vectors[:,12]
fig = Figure()
Label(fig[1, 1, Top()], "Default Angle", padding = (0, 0, 5, 0))
ax = LScene(fig[1,1], scenekw=(lights=[],))
msh = CairoMakie.mesh!(ax, s, color=q)
Colorbar(fig[1,2], msh, size=32)
# Second Angle
Label(fig[2, 1, Top()], "Bottom Angle", padding = (0, 0, 5, 0))
ax = LScene(fig[2,1], scenekw=(lights=[],))
update_cam!(ax.scene, Vec3f(-1/2,-1/2,1.0/2), Vec3f(1,1,1), Vec3f(0, 0, 1))
msh = CairoMakie.mesh!(ax, s, color=q)
Colorbar(fig[2,2], msh, size=32)
fig
```

# Exploring solutions with Krylov methods

We can also use the information about the eigenvectors for spectral techniques in solving the equations. Krylov methods are a bridge between
linear solvers and spectral information.

```julia
using Krylov
b = zeros(nv(sd))
b[1] = 1
b[end] = -1
x, stats = Krylov.gmres(Δ0, b, randn(nv(sd)), restart=true, memory=20, atol = 1e-10, rtol=1e-8, history=true, itmax=10000)
= x .- sum(x)./length(x)
norm(x̂)
stats
norm(Δ0*(x) - b)
```
1 change: 0 additions & 1 deletion docs/src/index.md
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,6 @@ Decapodes.jl is the synthesis of Applied Category Theory (ACT) techniques for fo

By combining the power of ACT and the DEC, we seek to improve the scientific computing workflow. Decapodes simulations are [hierarchically composable](https://algebraicjulia.github.io/Decapodes.jl/dev/budyko_sellers_halfar/), generalize over [any type of manifold](https://algebraicjulia.github.io/Decapodes.jl/dev/ice_dynamics/), and are [performant and accurate](https://www.cise.ufl.edu/~luke.morris/cism.html) with a declarative domain specific language (DSL) that is [human-readable](https://algebraicjulia.github.io/Decapodes.jl/dev/klausmeier/#Model-Representation).

![Grigoriev Ice Cap Dynamics](grigoriev/grigoriev.gif)

## Getting started

Expand Down

0 comments on commit e9779ad

Please sign in to comment.