Cahn-Hilliard equation

Go back to Examples.

In this example we will show how to solve the Cahn-Hilliard equation in Medusa. The solution provides an implementation of periodic boundary conditions and shows how to combine implicit and explicit operators. This example uses a custom biharmonic operator, so we recommend firstly reading the tutorial on Customization.

Cahn-Hilliard equation

The Cahn-Hilliard equation describes the separation of binary fluid into pure domains. We define concentration, a scalar field $c(x, y) \in [-1, 1]$ to describe the fluid, with values $c \pm 1$ representing pure phases. In this case the equation can be written as

\[ \frac{\partial c}{\partial t} = D \nabla^2 (c^3 - c - \gamma \nabla^2 c), \] where $D$ is the diffusion constant and $\sqrt{\gamma}$ is the characteristic transition length between domains.

We want to avoid any boundary effects and solve the equation for infinite fluid which can be approximated by using periodic boundary conditions.

Periodic boundary conditions

The simplest method of implementing a periodic boundary is adding a band of periodic nodes that mirror values from nodes that lie close to the opposite boundary. We have to create periodic nodes based on locations of interior nodes and construct a mapping that will allow us to update their values based on their interior counterparts. Width of the periodic band should be chosen with support size in mind. Periodic band should be wide enough to include all nodes required to calculate interior values.

TODO: node visualization with same colors for periodic nodes, support in corner visualization

We start by creating a square domain centered on 0 and uniformly filling it with nodes. $border$ parameter determines the square size while $dx$, the typical distance between nodes, determines the node density.

mm::BoxShape<vec_t> box(-border, border);
mm::DomainDiscretization<vec_t> domain = box.discretizeBoundaryWithStep(dx);
mm::GeneralFill<vec_t> fill;
fill.seed(randomSeed);
domain.fill(fill, dx);

Solution

We need to handle two uncommon terms to solve the equation. Biharmonic ($\nabla^4$) operator can be implemented with custom operators as shown in Customization example. Nonlinear ($\nabla^2 c^3$) term proves to be more problematic, as it can not be solved implicitly. We either need to solve the entire equation with explicit time stepping or include the nonlinear term in the right hand side of our system of equations.

Mixed solution

Go back to Examples.