Slope limiting for DG

variational formulation
1

Hi,

has anyone attempted to use slope limiters for DG schemes in fenics? Right now I’m solving a 2D euler-like system.

Is there any good references on slope limiters for DG methods or projects that use them?

There’s a description of a slope limiting procedure in this reference but right now I wouldn’t know how to start coding it.

Thanks in advance.

2

In my experience, “good” slope limiters are a real pain to implement in dolfin since you need to increase the compact support of bases to extend further than the local neighbours of cells.

Perhaps consider the work in Ocellaris, e.g. here.

1 Like
3

Thanks, I will check it.

At first glance the code looks quite complex indeed.

Do you know of any other options to tame problems arising from negative density solutions due to oscillations in higher order discretisations?

4

The (little) work I’ve done on approximations of non-smooth compressible flow solutions was not with FEniCS. I got around it with adaptive hp-refinement as developed by my PhD advisor. I’d resolve regions with large transitions or shocks with low order approximations.

With the new developments in basix I really hope hp-refinement could be in dolfinx’s future. But this requires a research question to make it worth the time to implement.

You could also consider artificial viscosity, e.g. eq (14) of here.

1 Like
5

Great, thanks for all the information.

6

Hi, @nate, I’m having some problems trying to implement the artificial viscosity solution. Namely, when I calculate any differential operator on the flux tensor I get a Zero object, for example:

import dolfin as df
import ufl as ufl
import numpy as np

mesh = df.UnitSquareMesh(32, 32)

V = df.VectorFunctionSpace(mesh, 'DG', 0, 3)

u = df.Expression(('sin(pi*x[0]) * sin(pi*x[1])',
                        '0.',
                        '0.'), pi = np.pi, degree=2)

u = df.project(u, V)

def F(u):
    return df.as_tensor([[u[1], u[2]],
                         [u[1]*u[1]/u[0], u[1]*u[2]/u[0]],
                         [u[1]*u[2]/u[0], u[2]*u[2]/u[0]]])

F = F(u)

div_F = ufl.div(F)

print(div_F)

outputs:

0 (shape (3,))

This doesn’t happen when i calculate derivatives of test functions in the convective part of the flux:

v = df.TestFunction(V)
print(ufl.grad(v))

Output:

grad(v_0)

Am I doing something wrong in the calculation of the diffusive flux?

7

Looks like your function space V_\mathrm{DG}^{h,p} is degree p=0, so \mathrm{div}(F) = 0 \quad \forall F \in V_\mathrm{DG}^{h,0}. A first order accurate DG method should be highly dissipative as-is and not require stabilisation.

1 Like
8

True, my bad. I was testing everything with the order zero discretization which indeed works well and I did not realize the gradients vanish there. Thanks for the quick reply.