Transport Local Operator for finite volume scheme
Description
The local operator must model the formulation of the transport equation in finite volumes. This local operator will be created taking into account the one created on !30 (closed), the one convectiondiffusionccfv.hh
implemented on PDELab, and the one implemented by @oklein in dune-modelling-example.
Tasks
-
Change the name of the local operator file for richards to richards_operator_DG.hh -
Implement the alpha_volume()
method. -
Implement the jacobian_volume()
method. -
Implement the alpha_skeleton()
method. -
Implement the jacobian_skeleton()
method. -
Implement the alpha_boundary()
method. -
Implement the jacobian_boundary()
method. -
Implement the lambda_volume()
method.
How to test the implementation?
Testing the local operator is hard without the external framework, so this task will be tied to those tasks implementing a TransportSimulation
.
Formulation
Warning: I have to check these equations again. It seems to be that there is a mistake. Source terms omitted by now.
The strong formulation for solute transport is
\begin{aligned}
\partial_t[\theta C_w] + \nabla\cdot [\textbf{j}_w C_w] - \nabla [\theta \mathsf{D}_{eff}\nabla C_w]=0 &\qquad \text{in } \Omega\\
C_w = g &\qquad \text{on } \Gamma_D \subseteq\partial\Omega\\
\nabla C_w \cdot \textbf{n} = \textbf{j}_{\scriptscriptstyle C_w}& \qquad \text{on } \Gamma_N =\partial\Omega \backslash \Gamma_D
\end{aligned}
with \textbf{j}_w = \theta \textbf{v}_w
and \mathsf{D}_{eff}(\textbf{j}_w,\theta)
. Now, following formulation in dune-pdelab-tutorial02
we use the ansatz function u=C_w
and the test function v
, we have that the weak formulation for the spatial part is
\int_\Omega \nabla\cdot [\textbf{j}_w u] v - \int_\Omega\nabla [\theta \mathsf{D}_{eff}\nabla u]v \qquad \forall v\in W_h
with W_h=\{w\in L^2(\Omega)\, : \, w|_T=\text{const for all } T\in\mathcal{T}_h\}
. Then, integrating by parts
\sum_{T\in\mathcal{T}_h}\left\{\int_{\partial T}(\textbf{j}_wu \cdot \textbf{n})v\,ds-\int_T(\textbf{j}_w\cdot\nabla v)u\,dx-\int_{\partial T}\theta\mathsf{D}_{eff}(\nabla u\cdot \textbf{n}) v\,ds+\int_T\theta\mathsf{D}_{eff}(\nabla u\cdot\nabla v)dx\right\}
Taking into account that \nabla u|_T=0
and approximating normal derivatives we have
\sum_{T\in\mathcal{T}_h}\left\{\int_{\partial T}(\textbf{j}_wu \cdot \textbf{n})v\,ds-\int_{\partial T}\theta\mathsf{D}_{eff}(\nabla u\cdot \textbf{n}) v\,ds\right\}
and rearranging
\begin{aligned}
\sum_{F\in\mathcal{F}_h^i}&\left\{\int_{\partial F}(\textbf{j}_w u\cdot\textbf{n}_F)[v(x_{T_F^-})-v(x_{T_F^+})]\right.\\
&-\left.\int_{\partial F}\mathsf{D}_{eff}(\nabla u\cdot \textbf{n}_F)[v(x_{T_F^-})\theta(x_{T_F^-})-v(x_{T_F^+})\theta(x_{T_F^+})] \,ds\right\} \\
+\sum_{F\in\mathcal{F}_h^{\partial\Omega}}&\left\{\int_{\partial F}\left[(\textbf{j}_wu \cdot \textbf{n}_F)-\mathsf{D}_{eff}(\nabla u\cdot \textbf{n}_F)\theta(x_{T_F^-})\right]v(x_{T_F^-})\,ds\right\},
\end{aligned}
\text{with}\quad\nabla u\cdot \textbf{n}_F=\frac{u_h(x_{T_F^+})-u_h(x_{T_F^-})}{||x_{T_F^+}-x_{T_F^-}||}+\text{error}
For the piecewise-constant elements in the test function, we have that the basis functions that generate the space W_h
are one on one element and zero on all others, i.e.
\phi_i(x) = \begin{cases}
1 &\text{if } x\in T_i \\
0 &\text{else}
\end{cases}
which is equivalent to v(x_{T_F^+})=0
on the equations above. Then,
\begin{aligned}
&\sum_{F\in\mathcal{F}_h^i}\int_{\partial F}\left[(\textbf{j}_w u \cdot\textbf{n}_F) - \mathsf{D}_{eff}(\nabla u\cdot \textbf{n}_F)\theta(x_{T_F^-})\right]\,ds \\
+&\sum_{F\in\mathcal{F}_h^{\partial\Omega}\cap\Gamma_D}\int_{\partial F}\left[(\textbf{j}_w g \cdot \textbf{n}_F)-\mathsf{D}_{eff}(\nabla^* u\cdot \textbf{n}_F)\theta(x_{T_F^-})\right]\,ds\\
-&\sum_{F\in\mathcal{F}_h^{\partial\Omega}\cap\Gamma_N}\int_{\partial F}\mathsf{D}_{eff}\textbf{j}_{\scriptscriptstyle C_w}\theta(x_{T_F^-})\,ds
\end{aligned}
where \nabla^* u\cdot \textbf{n}_F
is the finite difference between u
and the boundary condition g
in direction \textbf{n}_F
:
\nabla^* u\cdot \textbf{n}_F=\frac{u_h(x_F)-g(x_{T_F^-})}{||x_F-x_{T_F^-}||}+\text{error}
General Residual Form for Local Operators
\begin{aligned}
r(u,v)=&\sum_{T\in\mathcal{T}_h}\alpha_T^V(\mathcal{R}_Tu,\mathcal{R}_Tv)+\sum_{T\in\mathcal{T}_h}\lambda_T^V(\mathcal{R}_Tv)\\
+&\sum_{F\in\mathcal{F}_h^i}\alpha_F^S(\mathcal{R}_{T_F^-}u,\mathcal{R}_{T_F^+}u,\mathcal{R}_{T_F^-}v,\mathcal{R}_{T_F^+}v)\\
+&\sum_{F\in\mathcal{F}_h^{\partial\Omega}}\alpha_F^B(\mathcal{R}_{T_F^-}u,\mathcal{R}_{T_F^-}v)+\sum_{F\in\mathcal{F}_h^{\partial\Omega}}\alpha_F^B(\mathcal{R}_{T_F^-}v)
\end{aligned}