Hi
I am trying to use DG method with FEnics. it doesnt work at all. I dont know if it is problem of periodic BC?
from dolfin import *
mesh = UnitInterval(16)
c0=Expression('exp(-pow((x[0]-0.5),2)/0.025)')
# Sub domain for Periodic boundary condition
class PeriodicBoundary(SubDomain):
def inside(self, x, on_boundary):
return bool(x[0] < DOLFIN_EPS and x[0] > -DOLFIN_EPS and on_boundary)
def map(self, x, y):
y[0] = x[0] - 1.0
# Create function space constrained by periodic BC
pbc = PeriodicBoundary()
V = FunctionSpace(mesh, 'DG', 1, constrained_domain=pbc)
c=TrialFunction(V)
v=TestFunction(V)
a_M = c*v*dx
M=assemble(a_M)
# Define normal component, mesh size and right-hand side
n = FacetNormal(mesh)
c=Function(V)
c_1=interpolate(c0,V)
c.assign(c_1)
k=[]
for i in range(4):
k.append(Function(V))
T = 0.002
t = 0
dt =0.001
alpha=0
co=(1.-alpha)/2.
while t<=T:
flux = avg(c)+co*jump(c)
a=v*c.dx(0)*dx - dot(n,flux*v)*dS
b=assemble(a)
solve(M,k[0].vector(),b)
c.vector()[:]=k[0].vector()
t += dt
plot(c,rescale=False)