Skip to content
Snippets Groups Projects
Commit 8d127c89 authored by Andreas Dedner's avatar Andreas Dedner
Browse files

added a scalar test case (transport)

parent 826f0d47
No related branches found
No related tags found
1 merge request!4Latest features added to dune-fem-dg.
Pipeline #11248 failed
from ufl import *
from dune.ufl import Space
class Transport1D:
name = 'transport'
def F_c(U):
return [ [U[0], 0] ]
def outflowValue(u):
return u
boundaryValue = {1: outflowValue}
def maxLambda(U,n):
return abs(n[0])
def velocity(U):
return as_vector( [1,0] )
def physical(U):
return 1
def jump(U,V):
return U[0]-V[0]
space = Space(2,1)
def RiemanProblem(UL,UR):
x = SpatialCoordinate(space.cell())
return conditional(x[0]<-0.5,as_vector(UL),as_vector(UR))
def Constant():
return as_vector( [0.1] )
import time
from dune.grid import structuredGrid, cartesianDomain
from dune.fem import parameter
import dune.create as create
from dune.models.elliptic.formfiles import loadModels
from dune.femdg import createFemDGSolver
from scalar import RiemanProblem
from scalar import Transport1D as Model
parameter.append("parameter")
parameter.append({"fem.verboserank": -1})
x0,x1,N = [-1, 0], [1, 0.1], [20, 5]
grid = structuredGrid(x0,x1,N)
# grid = create.grid("ALUSimplex", cartesianDomain(x0,x1,N))
dimR = 1
t = 0
saveStep = 0.01
saveTime = saveStep
count = 0
endTime = 1.
def useODESolver():
global count, t, saveTime
initial = RiemanProblem( [1], [0] )
spaceName = "dgonb"
polOrder = 2
space = create.space(spaceName, grid, order=polOrder, dimrange=dimR)
u_h = space.interpolate(initial, name='u_h')
u_h_n = u_h.copy(name="previous")
operator = createFemDGSolver( Model, space )
start = time.time()
grid.writeVTK(Model.name, celldata=[u_h], number=count)
while t < endTime:
u_h_n.assign(u_h)
operator(u_h_n, u_h)
dt = operator.deltaT()
t += dt
if t > saveTime:
print('dt = ', dt, 'time = ',t, 'count = ',count )
count += 1
grid.writeVTK(Model.name, celldata=[u_h], number=count)
saveTime += saveStep
grid.writeVTK(Model.name, celldata=[u_h], number=count)
print("time loop:",time.time()-start)
useODESolver()
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment