What is a DuMu^x model

A DuMu^x model is an implementation of a discretized mathematical model, generally given by partial differential equations.

What is a DuMu^x model

Mathematical model (PDE): \(\begin{equation*} \small \frac{\partial S(u)}{\partial t} + \nabla \cdot \mathbf{F}(u) = q, \quad \forall (t,\mathbf{x}) \in (0,T] \times \Omega \end{equation*}\)

Discrete model, e.g. using finite volumes: \(\begin{equation*} \small |B| \frac{S_h(\mathbf{u}^{n+1}_h) - S_h(\mathbf{u}^{n}_h)}{\Delta t} + \sum_{\sigma \in \Sigma_B} F_{B,\sigma}(\mathbf{u}^{n+1}_h) = \int_{B} q^{n+1} \, dx, \quad \forall t_{n+1}\leq T, \; \forall B \end{equation*}\)

What is a DuMu^x model

Discrete model, e.g. using finite volumes: \(\begin{equation*} \small |B| \frac{S_h(\mathbf{u}^{n+1}_h) - S_h(\mathbf{u}^{n}_h)}{\Delta t} + \sum_{\sigma \in \Sigma_B} F_{B,\sigma}(\mathbf{u}^{n+1}_h) = \int_{B} q \, dx, \quad \forall t_{n+1}\leq T, \; \forall B \end{equation*}\)

• \(S_h:\) storage term
• \(F_{B,\sigma}:\) flux term over sub control volume face (scvf)
• \(q:\) source term

Where to implement these terms in DuMu^x?

LocalResidual

LocalResidual

Implements terms of a PDE in the functions

• computeStorage(...)
• computeFlux(...)
• computeSource(...)

Example: Diffusion equation

Mathematical model (PDE):

\(\begin{equation} \frac{\partial c}{\partial t} - \nabla \cdot (D \nabla c) = 0 \:\: \mathrm{in}\; \Omega \times (0,T] \\ \end{equation}\)

with

• \(c:\) concentration
• \(D:\) constant diffusion coefficient
• \(\Omega:\) spatial domain
• \(T:\) end time

Example: Diffusion equation

Discrete model using the Box discretization:

\(\begin{equation} \vert B \vert \frac{c_B^{n+1}-c_B^{n}}{\Delta t} - \sum_{\sigma \in \Sigma_{B}} \left[ D \nabla c_h^{n+1} \cdot \boldsymbol{n}_{B,\sigma} \vert \sigma \vert \right] = 0, \end{equation}\)

with

• \(c_B^n:\) concentration at time \(t_n\) and control volume \(B\)
• \(c^n_h:\) global discrete solution at time \(t_n\), interpolated using basis functions
• \(\mathbf{n}:\) unit outer normal vector
• \(\sigma:\) sub control volume face (scvf)

Example: Diffusion equation

Discrete model using the Box discretization:

\(\begin{equation} \vert B \vert \frac{c_B^{n+1}-c_B^{n}}{\Delta t} - \sum_{\sigma \in \Sigma_{B}} \left[ D \nabla c_h^{n+1} \cdot \boldsymbol{n}_{B,\sigma} \vert \sigma \vert \right] = 0, \end{equation}\)

LocalResidual

The local residual of the diffusion model:

template<class TypeTag>
class DiffusionModelLocalResidual
: public GetPropType<TypeTag, Properties::BaseLocalResidual>
{
    ...
}

Inherits from the BaseLocalResidual, which is chosen depending on the discretization scheme, here Box scheme.

Storage term

NumEqVector computeStorage(const Problem& problem,
                           const SubControlVolume& scv,
                           const VolumeVariables& volVars) const
 {
    NumEqVector storage;
    storage[Indices::massBalanceEqIdx]
        = volVars.priVar(Indices::concentrationIdx);
    return storage;
 }

Flux term

\(\begin{equation} F_{B,\sigma} = -D \nabla c_h^{n+1} \cdot \boldsymbol{n}_{B,\sigma} \vert \sigma \vert \end{equation}\)

with

• \(c^n_h:\) global discrete solution at time \(t_n\), interpolated using basis functions
• \(\mathbf{n}:\) unit outer normal vector
• \(\sigma:\) sub control volume face (scvf)

Flux term

NumEqVector computeFlux(const Problem& problem,
                        const Element& element,
                        const FVElementGeometry& fvGeometry,
                        const ElementVolumeVariables& elemVolVars,
                        const SubControlVolumeFace& scvf,
                        const ElementFluxVariablesCache& elemFluxVarsCache) const
{
    ...
}

Flux term

NumEqVector computeFlux(...) const
{
    // Compute ∇c
    const auto& fluxVarCache = elemFluxVarsCache[scvf];
    Dune::FieldVector<Scalar, dimWorld> gradConcentration(0.0);
    for (const auto& scv : scvs(fvGeometry))
    {
        const auto& volVars = elemVolVars[scv];
        gradConcentration.axpy(
            volVars.priVar(Indices::concentrationIdx),
            fluxVarCache.gradN(scv.indexInElement())
        );
    }
    ...
}

Flux term

NumEqVector computeFlux(...) const
{
    ...
    NumEqVector flux;

    // Compute the flux
    flux[Indices::massBalanceEqIdx] = -1.0*scvf.area()*vtmv(
        scvf.unitOuterNormal(),
        problem.diffusionCoefficient(),
        gradConcentration
    );

    return flux;
}

LocalResidual

A LocalResidual implements the discretized mathematical model.

For its implementation different model-specific properties have to be set

Model type tag

The property tag is an empty struct with the respective name, e.g. DiffusionModel

namespace Dumux::Properties::TTag {
//! The diffusion model tag that we can specialize properties for
struct DiffusionModel {};
} // end namespace Dumux::Properties::TTag

Model properties

We can set nodel properties for the DiffusionModel type tag

All properties are defined within the namespace Dumux::Properties

namespace Dumux::Properties {
    //define all properties
} // end namespace Dumux::Properties

Defining model properties

The type of the local residual is the class DiffusionModelLocalResidual defined from earlier

template<class TypeTag>
struct LocalResidual<TypeTag, TTag::DiffusionModel>
{ using type = DiffusionModelLocalResidual<TypeTag>; };

Defining model properties

The model traits specify information about the model:

template<class TypeTag>
struct ModelTraits<TypeTag, TTag::DiffusionModel>
{
    struct type
    {
        struct Indices
        {
            static constexpr int concentrationIdx = 0;
            static constexpr int massBalanceEqIdx = 0;
        };
        static constexpr int numEq() { return 1; }
    };
};

Defining model properties

Further model specific properties can be set accordingly by using the model property tag, i.e. TTag::DiffusionModel

Exercise: Model

Implementation of a nonlinear diffusion model for denoising of an MRI image

Tasks

• Implement local residual
• Set model properties
• Use model in test case
• Customize volume variables