Mass conservation law

Question

I am trying to implement a simple mass conservation equation:

$\frac{dh}{dt} + \frac{dq}{dx}= a$

on [0,1], where q=u*h. I set h=1 on the first half of the domain (x<0.5) and zero on the second half. Setting u=1 and a=1, I expect this shock to simply move to the right. I have implemented this UFL file:

cell = interval

element = FiniteElement("Discontinuous Lagrange", cell, 2)

h = TrialFunction(element)
phi = TestFunction(element)
ac = Coefficient(element)
u = Coefficient(element)
h0 = Coefficient(element)

dt = Constant(cell)

flux = h*u

a = phi*h/dt*dx - phi.dx(0)*flux*dx + phi("+")*flux("+")*dS - phi("-")*flux("-")*dS
L = phi*h0/dt*dx + ac*phi*dx

but when I implement the test case described above the solution blows up at the shock.

I don't think I've implemented the weak form correctly --- how should I represent the flux between elements?

Thanks!

Comments

What do you see if you use a piecewise constant approximation? High order DG methods are not usually stable at shocks. Do you have a link to a document which details the DG FEM formulation you're trying to implement?

---

I don't have a great link to a document. I'm not 100% confident that I'm deriving the DG FEM formulation correctly --- I'm quite inexperienced with the theory. Starting from

du/dt - b*du/dx = f

with f = 0, b = 1, and u=1 at x=0 and u=1 if x<0.5; u=0 else. I have implemented the (naive) ufl file

cell = interval 

scalar = FiniteElement("Lagrange", cell, 2)
vector = VectorElement("Lagrange", cell, 2)

u = TrialFunction(scalar)
v = TestFunction(scalar)
u0 = Coefficient(scalar)
b = Coefficient(vector)
f = Coefficient(scalar)

dt = Constant(cell)

a = u*v*dx + 0.5*dt*dot(b, grad(u)*v)*dx
L = u0*v*dx - 0.5*dt*dot(b, grad(u0)*v)*dx + dt*f*v*dx

which behaves exactly as you would expect --- the shock propagates across the domain but a bunch of oscillations appear. I assume this is because of the "hyperbolic-ness" of the PDE. I switched to DG hoping to deal with this but very quickly my inexperience with DG became a problem.