Free surface calculations

2.12 Free surface calculations

In reality the boundary conditions of a convecting Earth are not no-slip or free slip (i.e., no normal velocity). Instead, we expect that a free surface is a more realistic approximation, since air and water should not prevent the flow of rock upward or downward. This means that we require zero stress on the boundary, or $σ \cdot n = 0$ , where $σ = 2 η ε (u)$ . In general there will be flow across the boundary with this boundary condition. To conserve mass we must then advect the boundary of the domain in the direction of fluid flow. Thus, using a free surface necessitates that the mesh be dynamically deformable.

2.12.1 Arbitrary Lagrangian-Eulerian implementation

The question of how to handle the motion of the mesh with a free surface is challenging. Eulerian meshes are well behaved, but they do not move with the fluid motions, which makes them difficult for use with free surfaces. Lagrangian meshes do move with the fluid, but they quickly become so distorted that remeshing is required. ASPECT implements an Arbitrary Lagrangian-Eulerian (ALE) framework for handling motion of the mesh. The ALE approach tries to retain the benefits of both the Lagrangian and the Eulerian approaches by allowing the mesh motion $u_{m}$ to be largely independent of the fluid. The mass conservation condition requires that $u_{m} \cdot n = u \cdot n$ on the free surface, but otherwise the mesh motion is unconstrained, and should be chosen to keep the mesh as well behaved as possible.

ASPECT uses a Laplacian scheme for calculating the mesh velocity. The mesh velocity is calculated by solving

\begin{align} - Δ u_{m} & = 0 & in Ω, & (36) \\ u_{m} & = (u \cdot n) n & on \partial Ω_{free surface}, & (37) \\ u_{m} \cdot n & = 0 & on \partial Ω_{free slip}, & (38) \\ u_{m} & = 0 & on \partial Ω_{Dirichlet} . & (39) \end{align}

After this mesh velocity is calculated, the mesh vertices are time-stepped explicitly. This scheme has the effect of choosing a minimally distorting perturbation to the mesh. Because the mesh velocity is no longer zero in the ALE approach, we must then correct the Eulerian advection terms in the advection system with the mesh velocity (see, e.g. [DHPRF04]). For instance, the temperature equation (31) becomes

ρ C_{p} (\frac{\partial T}{\partial t} + (u - u_{m}) \cdot \nabla T) - \nabla \cdot k \nabla T = ρ H in Ω .

2.12.2 Free surface stabilization

Small disequilibria in the location of a free surface can cause instabilities in the surface position and result in a “sloshing” instability. This may be countered with a quasi-implicit free surface integration scheme described in [KMM10]. This scheme enters the governing equations as a small stabilizing surface traction that prevents the free surface advection from overshooting its true position at the next time step. ASPECT implements this stabilization, the details of which may be found in [KMM10].

An example of a simple model which uses a free surface may be found in Section 5.2.6.

[next] [prev] [prev-tail] [front] [up]