|
Basilisk CFD
Adaptive Cartesian mesh PDE framework
|
#include "poisson.h"Go to the source code of this file.
Functions | |
| void | vertical_diffusion (Point point, scalar h, scalar s, double dt, double D, double dst, double s_b, double lambda_b) |
| void | event_viscous_term (void) |
| In the layered solver, vertical viscosity is applied to the velocity field just after advection, but before the pressure gradient/acceleration term is applied. | |
| void | horizontal_diffusion (scalar *list, double D, double dt) |
| trace mgstats | implicit_horizontal_diffusion (scalar *list, double D, double dt) |
Variables | |
| double | nu = 0. |
| const vector | lambda_b [] = {0,0,0} |
| const vector | dut [] = {0,0,0} |
| const vector | u_b [] = {0,0,0} |
In the layered solver, vertical viscosity is applied to the velocity field just after advection, but before the pressure gradient/acceleration term is applied.
To take the pressure gradient into account, we first apply the acceleration of the previous timestep, apply vertical viscosity and then substract the previous acceleration.
Event: viscous_term (i++,last)
Definition at line 147 of file diffusion.h.
References _i, dimension, dry, dt, dut, foreach_layer, h, ha, hf, lambda_b, nu, point, u, u_b, vertical_diffusion(), vector::x, and x.
This approximates
\[ h \partial_t s = D \nabla \cdot (h \nabla s) \]
with \(D\) the diffusion coefficient. Note that metric terms linked to the slope of the layers are not taken into account. Note also that the time discretisation is explicit so that the timestep must be limited (manually) by \(\min(\Delta^2/D)\).
Definition at line 176 of file diffusion.h.
References _i, a, cm, D, dimension, dry, dt, fm, foreach_layer, free(), h, hf, list, list_clone(), s, sq(), vector::x, and x.
When the timestep constraint due to explicit time integration is too strong, one can use the time-implicit version below. The diffusion equation
\[ h \partial_t s = D \nabla \cdot (h \nabla s) \]
is discretised implicitly as
\[ h \frac{s^{n + 1} - s^n}{\Delta t} = \nabla \cdot (Dh \nabla s^{n + 1}) \]
which can be reordered as
\[ \nabla \cdot (\Delta tDh \nabla s^{n + 1}) - hs^{n + 1} = - hs^n \]
i.e. a Poisson–Helmoltz equation of the form
\[ \nabla \cdot (\alpha \nabla s^{n + 1}) + \lambda s^{n + 1} = b \]
which is solved using the multigrid Poisson–Helmholz solver.
The function returns convergence statistics for the multigrid solver.
Definition at line 230 of file diffusion.h.
References _i, alpha, b, cm, D, dt, fm, foreach_layer, h, hf, lambda, poisson, s, vector::x, and x.
| void vertical_diffusion | ( | Point | point, |
| scalar | h, | ||
| scalar | s, | ||
| double | dt, | ||
| double | D, | ||
| double | dst, | ||
| double | s_b, | ||
| double | lambda_b | ||
| ) |
We consider the vertical diffusion of a tracer \(s\) with a diffusion coefficient \(D\) for the multilayer solver.
For stability, we discretise the vertical diffusion equation implicitly as
\[ \frac{(hs_l)^{n + 1} - (hs_l)^{\star}}{\Delta t} = D \left( \frac{s_{l + 1} - s_l}{h_{l + 1 / 2}} - \frac{s_l - s_{l - 1}}{h_{l - 1 / 2}} \right)^{n + 1} \]
which can be expressed as the linear system
\[ \mathbf{Ms}^{n + 1} = \mathrm{rhs} \]
where \(\mathbf{M}\) is a tridiagonal matrix. The lower, principal and upper diagonals are a, b and c respectively.
Boundary conditions on the top and bottom layers need to be added to close the system. We chose to impose a Neumann condition on the free-surface i.e.
\[ \partial_z s |_t = \dot{s}_t \]
and a Navier slip condition on the bottom i.e.
\[ s|_b = s_b + \lambda_b \partial_z s|_b \]
The rhs of the tridiagonal system is \(h_l s_l\).
The lower, principal and upper diagonals \(a\), \(b\) and \(c\) are given by
\[ a_{l > 0} = - \left( \frac{D \Delta t}{h_{l - 1 / 2}} \right)^{n + 1} \]
\[ c_{l < \mathrm{nl} - 1} = - \left( \frac{D \Delta t}{h_{l + 1 / 2}} \right)^{n + 1} \]
\[ b_{0 < l < \mathrm{nl} - 1} = h_l^{n + 1} - a_l - c_l \]
For the top layer the boundary conditions give the (ghost) boundary value
\[ s_{\mathrm{nl}} = s_{\mathrm{nl} - 1} + \dot{s}_t h_{\mathrm{nl} - 1}, \]
which gives the diagonal coefficient and right-hand-side
\[ b_{\mathrm{nl} - 1} = h_{\mathrm{nl} - 1}^{n + 1} - a_{\mathrm{nl} - 1} \]
\[ \mathrm{rhs}_{\mathrm{nl} - 1} = (hs)_{\mathrm{nl} - 1}^{\star} + D \Delta t \dot{s}_t \]
For the bottom layer a third-order discretisation of the Navier slip condition gives
\[ \begin{aligned} b_0 & = h_0 + 2 \Delta t D \left( \frac{1}{h_0 + h_1} + \frac{h^2_1 + 3 h_0 h_1 + 3 h^2_0}{\det} \right),\\ c_0 & = - 2 \Delta t D \left( \frac{1}{h_0 + h_1} + \frac{h^2_0}{\det} \right),\\ \text{rhs}_0 & = (hs_0)^{\star} + 2 \Delta t D s_b \frac{h^2_1 + 3 h_0 h_1 + 2 h^2_0}{\det},\\ \det & = h_0 (h_0 + h_1)^2 + 2\lambda (3\,h_0 h_1 + 2\,h_0^2 + h_1^2), \end{aligned} \]
We can now solve the tridiagonal system using the Thomas algorithm.
Definition at line 33 of file diffusion.h.
References a, b, c, D, dst, dt, foreach_layer, h, l, lambda_b, nl, point, s, sq(), and x.
Referenced by event_viscous_term().
Definition at line 136 of file diffusion.h.
Referenced by event_viscous_term().
Definition at line 136 of file diffusion.h.
Referenced by event_viscous_term(), remap_c(), remap_robin(), right_value(), and vertical_diffusion().
| double nu = 0. |
By default the viscosity is zero and we impose free-slip on the free-surface and no-slip on the bottom boundary i.e. \(\dot{\mathbf{u}}_t = 0\), \(\mathbf{\lambda}_b = 0\), \(\mathbf{u}_b = 0\).
Definition at line 135 of file diffusion.h.
Referenced by event_tracer_advection(), event_viscous_term(), and fenep().
Definition at line 136 of file diffusion.h.
Referenced by event_viscous_term().