Trefftz + Tent pitching
Tent pitching techniques generate a space-time mesh, which complies with the causality properties of the hyperbolic PDE. The resulting mesh consists of tent shaped objects, the height of rach tent is complient with the local wavespeed. The PDE can be solved explicitly in each tent and independant tents can be solves in parallel.
[1]:
from ngstrefftz import *
from ngsolve.TensorProductTools import *
from netgen.occ import *
from ngsolve import *
from ngsolve.webgui import *
Constructing a tent-pitched mesh
We use ngstents to construct tent-pitched meshes.
[2]:
def Lmesh(maxh):
mu = 0.4
mp = MeshingParameters (maxh = maxh)
minh = math.pow(maxh,1/mu)
refpoints = int(1/minh)
for i in range(0, refpoints+1):
for j in range(0, refpoints+1):
xk = (i/refpoints-0.5)*2
yk = (j/refpoints-0.5)*2
r = sqrt(xk*xk+yk*yk)
mp.RestrictH (x=xk, y=yk, z=0, h=max(minh, maxh*math.pow(r,1-mu)) )
wp = WorkPlane()
wp.MoveTo(-1,-1).LineTo(0,-1).LineTo(0,0).LineTo(1,0).LineTo(1,1).LineTo(-1,1).Close()
geo = wp.Face()
geo = OCCGeometry(geo,dim=2)
return Mesh(geo.GenerateMesh(mp=mp))
maxh = 0.13
mesh = Lmesh(maxh)
Draw(mesh)
[2]:
BaseWebGuiScene
[3]:
dt = 0.05 # tent slab height (in time)
wavespeed = 1
ts = TentSlab(mesh, method="edge")
ts.SetMaxWavespeed(wavespeed)
ts.PitchTents(dt=dt, local_ct=True)
print("max slope", ts.MaxSlope())
print("n tents", ts.GetNTents())
max slope 1.0
n tents 780
Trefftz-DG for the acoustic wave equation + Tents
[4]:
for i in range(0,len(mesh.GetBoundaries())):
mesh.ngmesh.SetBCName(i,"neumann")
order=4
delta = 500
ee = exp(-(x+1)**2*delta )
initc=((ee,-2*(x+1)*delta * ee, 0, 0))
bdd = CoefficientFunction(( 0,0,0,0))
TT=TWave(order,ts,CoefficientFunction(wavespeed))
TT.SetInitial(initc)
TT.SetBoundaryCF(bdd)
V = L2(mesh,order=order)
u = GridFunction(V,"u")
gfshow = GridFunction(V, multidim=0)
t=0
SetNumThreads(4)
with TaskManager():
while t < 1.5:
TT.Propagate()
t += dt
TT.GetWave(u)
gfshow.AddMultiDimComponent(u.vec)
Draw (gfshow,mesh,interpolate_multidim=True,deformation=True,animate=True,autoscale=False,min=0,max=0.3)
[4]:
BaseWebGuiScene
Quasi-Trefftz + Tents
if a non-constant coefficient function as wavespeed is passed to \(\texttt{TWave}\) the quasi-Trefftz method will be used, see also the notebook on the quasi-Treffty method.
[5]:
order = 4
delta=500
dt=0.02
ee = exp(-(x)**2*delta )
initc=((ee,-2*(x)*delta * ee, 0, 0))
bdd = CoefficientFunction(( 0,0,0,0))
mesh = Mesh(unit_square.GenerateMesh(maxh = 0.05))
for i in range(0,len(mesh.GetBoundaries())):
mesh.ngmesh.SetBCName(i,"neumann")
ts = TentSlab(mesh, method="edge")
ts.SetMaxWavespeed(y+1)
ts.PitchTents(dt=dt, local_ct=True)
TT=TWave(order,ts,CoefficientFunction(y+1))
TT.SetInitial(initc)
TT.SetBoundaryCF(bdd)
V = L2(mesh,order=order)
u = GridFunction(V,"u")
gfshow = GridFunction(V, multidim=0)
t=0
SetNumThreads(4)
with TaskManager():
while t < 0.8:
TT.Propagate()
t += dt
TT.GetWave(u)
gfshow.AddMultiDimComponent(u.vec)
Draw (gfshow,mesh,interpolate_multidim=True,deformation=True,animate=True,autoscale=False,min=0,max=0.3)
[5]:
BaseWebGuiScene