Instationary transport equation

23.2. Instationary transport equation#

Let \(w \in L_{\infty}^d\) be a given vectorfield, called the wind. We search for a concentration \(u : \Omega \times [0,T] \) such that

\[ \frac{\partial u}{\partial t} + \operatorname{div} (w u) = f \]

We split the boundary into inflow and outflow boundaries:

\[\begin{align*} \Gamma_{in} &= \{ x \in \partial \Omega : w \cdot n < 0 \} \\ \Gamma_{out} &= \{ x \in \partial \Omega : w \cdot n \geq 0 \} \\ \end{align*}\]

On the inflow, we can specify the concentration as a Dirichlet boundary condition:

\[ u(x,t) = u_{in} \qquad \forall \, x \in \Gamma_{in} \, \forall \, t \in [0,T] \]
from ngsolve import *
from ngsolve.webgui import Draw
mesh = Mesh(unit_square.GenerateMesh(maxh=0.05))

fes = H1(mesh, order=4, dirichlet="left|bottom")
u,v = fes.TnT()
gfu = GridFunction(fes)

wind = CF( (1, 0.2) )
gfu.Set (exp(-100*( (x-0.2)**2 + (y-0.5)**2) ))

n = specialcf.normal(2)
Atrans = BilinearForm(-u*wind*grad(v)*dx + wind*n*u*v*ds).Assemble()
f = LinearForm(fes)
mass = BilinearForm(u*v*dx).Assemble()
massinv = mass.mat.Inverse(freedofs=fes.FreeDofs(), inverse="sparsecholesky")

scene = Draw(gfu, min=0, max=1)

dt = 0.0002
tend = 0.8
t = 0
while t < tend:
    gfu.vec.data += dt * massinv * (f.vec-Atrans.mat*gfu.vec)
    t += dt
    scene.Redraw()

23.2.1. Exercises:#

  • try different wind fields, e.g. \(w(x,y) = (1, 0.1 \sin (30 x))\)

  • improve timestepping to \(2^{nd}\) order Runge Kutta

  • add a diffusion term to the transport equation:

    \[ \frac{\partial u}{\partial t} - \varepsilon \Delta u + \operatorname{div} (w u) = f \]

  • add a source term \(f\)

  • discretization methods tailored for transport equations are Discontinuous Galerkin methods, https://jschoeberl.github.io/iFEM/DG/stationary.html

23.2.2. Reynolds’ transport theorem#

The transport equation arises from conservation equations on moving control volumina.

Let

\[ \Phi : \Omega \times [0,T] \rightarrow {\mathbb R}^n \]

describe the deformation of a reference domain \(\Omega\), e.g. the particles of a continuum.

Call \(V_0 \subset \Omega\) a control volume, e.g. all particles colored with red paint.

Define

\[ V(t) = \{ \Phi(x,t) : x \in V_0 \} \]

Theorem: For proper functions \(u(x,t)\) there holds

\[ \frac{d}{dt} \int_{V(t)} u(x,t) \, dx = \int_{V(t)} \frac{\partial u}{\partial t} \, dx + \int_{\partial V(t)} w \cdot n \, u \, ds \]

with velocity \(w(x,t) = \frac{\partial \Phi}{\partial t}(x,t)\).

Next, apply Gauss’ theorem to obtain:

\[ \frac{d}{dt} \int_{V(t)} u(x,t) \, dx = \int_{V(t)} \frac{\partial u}{\partial t} + \operatorname{div} (w u) \]

If some quantity \(\int_V u \, dx\) is preserved over time for all moving control volumina, then it satisfies the transport equation

\[ \frac{\partial u}{\partial t} + \operatorname{div} (w u) = 0 \]

Examples of preserved quantities in fluid dynamics are mass \(\int \rho \, dx\), and momentum \(\int \rho u dx\)