Momentum equation

momentum_equation

Derived terms and associated equations for the Stokes system.

All terms are considered as if they were on the right-hand side of the equation, leading to the following UFL expression returned by the residual method:

\[ (dq)/dt = sum "term.residual()" \]

This sign convention ensures compatibility with Thetis's time integrators. In general, however, we like to think about the terms as they are on the left-hand side. Therefore, in the residual methods below, we first sum the terms in the variable F as if they were on the left-hand side, i.e.

\[ (dq)/dt + F(q) = 0, \]

and then return -F.

Users should not interact with these classes; instead, please use the solver provided in the stokes_integrators module.

ViscosityTerm(test_space, trial_space, dx, ds, dS, **kwargs)

Bases: BaseTerm

Viscosity term \(-nabla * (mu nabla u)\) in the momentum equation.

Using the symmetric interior penalty method, the weak form becomes

\[ {:( -int_Omega nabla * (mu grad u) phi dx , = , int_Omega mu (grad phi) * (grad u) dx ), ( , - , int_(cc"I" uu cc"I"_v) "jump"(phi bb n) * "avg"(mu grad u) dS - int_(cc"I" uu cc"I"_v) "jump"(u bb n) * "avg"(mu grad phi) dS ), ( , + , int_(cc"I" uu cc"I"_v) sigma "avg"(mu) "jump"(u bb n) * "jump"(phi bb n) dS ) :} \]

where σ is a penalty parameter (see Epshteyn and Riviere, 2007).

Epshteyn, Y., & Rivière, B. (2007). Estimation of penalty parameters for symmetric interior penalty Galerkin methods. Journal of Computational and Applied Mathematics, 206(2), 843-872.

Source code in g-adopt/gadopt/equations.py

def __init__(
    self,
    test_space: firedrake.functionspaceimpl.WithGeometry,
    trial_space: firedrake.functionspaceimpl.WithGeometry,
    dx: firedrake.Measure,
    ds: firedrake.Measure,
    dS: firedrake.Measure,
    **kwargs,
):
    self.test_space = test_space
    self.trial_space = trial_space

    self.dx = dx
    self.ds = ds
    self.dS = dS

    self.mesh = test_space.mesh()
    self.dim = self.mesh.geometric_dimension()
    self.n = firedrake.FacetNormal(self.mesh)

    self.term_kwargs = kwargs

MomentumEquation(test_space, trial_space, quad_degree=None, **kwargs)

Bases: BaseEquation

Momentum equation with viscosity, pressure gradient, and source terms.

Source code in g-adopt/gadopt/equations.py

def __init__(
    self,
    test_space: firedrake.functionspaceimpl.WithGeometry,
    trial_space: firedrake.functionspaceimpl.WithGeometry,
    quad_degree: Optional[int] = None,
    **kwargs
):
    self.test_space = test_space
    self.trial_space = trial_space
    self.mesh = trial_space.mesh()

    p = trial_space.ufl_element().degree()
    if isinstance(p, int):  # isotropic element
        if quad_degree is None:
            quad_degree = 2*p + 1
    else:  # tensorproduct element
        p_h, p_v = p
        if quad_degree is None:
            quad_degree = 2*max(p_h, p_v) + 1

    if trial_space.extruded:
        # Create surface measures that treat the bottom and top boundaries similarly
        # to lateral boundaries. This way, integration using the ds and dS measures
        # occurs over both horizontal and vertical boundaries, and we can also use
        # "bottom" and "top" as surface identifiers, for example, ds("top").
        self.ds = CombinedSurfaceMeasure(self.mesh, quad_degree)
        self.dS = (
            firedrake.dS_v(domain=self.mesh, degree=quad_degree) +
            firedrake.dS_h(domain=self.mesh, degree=quad_degree)
        )
    else:
        self.ds = firedrake.ds(domain=self.mesh, degree=quad_degree)
        self.dS = firedrake.dS(domain=self.mesh, degree=quad_degree)

    self.dx = firedrake.dx(domain=self.mesh, degree=quad_degree)

    # self._terms stores the actual instances of the BaseTerm-classes in self.terms
    self._terms = []
    for TermClass in self.terms:
        self._terms.append(TermClass(test_space, trial_space, self.dx, self.ds, self.dS, **kwargs))

ContinuityEquation(test_space, trial_space, quad_degree=None, **kwargs)

Bases: BaseEquation

Mass continuity equation with a single divergence term.

Source code in g-adopt/gadopt/equations.py

def __init__(
    self,
    test_space: firedrake.functionspaceimpl.WithGeometry,
    trial_space: firedrake.functionspaceimpl.WithGeometry,
    quad_degree: Optional[int] = None,
    **kwargs
):
    self.test_space = test_space
    self.trial_space = trial_space
    self.mesh = trial_space.mesh()

    p = trial_space.ufl_element().degree()
    if isinstance(p, int):  # isotropic element
        if quad_degree is None:
            quad_degree = 2*p + 1
    else:  # tensorproduct element
        p_h, p_v = p
        if quad_degree is None:
            quad_degree = 2*max(p_h, p_v) + 1

    if trial_space.extruded:
        # Create surface measures that treat the bottom and top boundaries similarly
        # to lateral boundaries. This way, integration using the ds and dS measures
        # occurs over both horizontal and vertical boundaries, and we can also use
        # "bottom" and "top" as surface identifiers, for example, ds("top").
        self.ds = CombinedSurfaceMeasure(self.mesh, quad_degree)
        self.dS = (
            firedrake.dS_v(domain=self.mesh, degree=quad_degree) +
            firedrake.dS_h(domain=self.mesh, degree=quad_degree)
        )
    else:
        self.ds = firedrake.ds(domain=self.mesh, degree=quad_degree)
        self.dS = firedrake.dS(domain=self.mesh, degree=quad_degree)

    self.dx = firedrake.dx(domain=self.mesh, degree=quad_degree)

    # self._terms stores the actual instances of the BaseTerm-classes in self.terms
    self._terms = []
    for TermClass in self.terms:
        self._terms.append(TermClass(test_space, trial_space, self.dx, self.ds, self.dS, **kwargs))

StokesEquations(test_space, trial_space, quad_degree=None, **kwargs)

Stokes system involving the momentum and mass continuity equations.

Parameters:

Name	Type	Description	Default
test_space	WithGeometry	Firedrake function space of the test function	required
trial_space	WithGeometry	Firedrake function space of the trial function	required
quad_degree	Optional[int]	Quadrature degree. Default value is 2p + 1, where p is the polynomial degree of the trial space	None

Returns:

Type	Description
list[BaseEquation]	A list of equation instances for the Stokes system.

Source code in g-adopt/gadopt/momentum_equation.py

def StokesEquations(
    test_space: fd.functionspaceimpl.WithGeometry,
    trial_space: fd.functionspaceimpl.WithGeometry,
    quad_degree: Optional[int] = None,
    **kwargs,
) -> list[BaseEquation]:
    """Stokes system involving the momentum and mass continuity equations.

    Arguments:
      test_space: Firedrake function space of the test function
      trial_space: Firedrake function space of the trial function
      quad_degree: Quadrature degree. Default value is `2p + 1`, where
                     p is the polynomial degree of the trial space

    Returns:
      A list of equation instances for the Stokes system.

    """
    mom_eq = MomentumEquation(
        test_space.sub(0), trial_space.sub(0), quad_degree=quad_degree, **kwargs
    )
    cty_eq = ContinuityEquation(
        test_space.sub(1), trial_space.sub(1), quad_degree=quad_degree, **kwargs
    )
    return [mom_eq, cty_eq]