Periodicity¶

Netgen is able to create periodic meshes, with them NGSolve can create periodic finite element spaces by mapping the degrees of freedom on the slave boundary to the ones on the master boundary. To use this functionality create a periodic mesh first:

Periodic Meshes¶

3D¶

If a geometry object has a periodic identification, Netgen creates a periodic mesh for the specified surfaces. In 3 dimensions, this identification is created using CSGeometry.PeriodicSurfaces(master,slave). Here we define a cube with periodicity in x-direction and create a periodic NGSolve mesh for it:

from netgen.csg import *

geo = CSGeometry()
left = Plane(Pnt(0,0,0),Vec(-1,0,0))
right = Plane(Pnt(1,0,0),Vec(1,0,0))
brick = OrthoBrick(Pnt(-1,0,0),Pnt(2,1,1)).bc("outer")
cube = brick * left * right
geo.Add(cube)
geo.PeriodicSurfaces(left,right)

from ngsolve import *
mesh = Mesh(geo.GenerateMesh(maxh=0.3))

2D¶

In two dimensions you can create periodic geometries with the copy parameter of the SplineGeometry.Append method. This defines the appended spline as the slave of the one given in copy:

Caution

The slave boundaries must be defined in the same direction as the master ones!

from netgen.geom2d import *

periodic = SplineGeometry()
pnts = [ (0,0), (1,0), (1,1), (0,1) ]
pnums = [periodic.AppendPoint(*p) for p in pnts]

periodic.Append ( ["line", pnums[0], pnums[1]],bc="outer")
# This should be our master edge so we need to save its number.
lright = periodic.Append ( ["line", pnums[1], pnums[2]], bc="periodic")
periodic.Append ( ["line", pnums[2], pnums[3]], bc="outer")
# Minion boundaries must be defined in the same direction as master ones,
# this is why the the point numbers of this spline are defined in the reverse direction,
# leftdomain and rightdomain must therefore be switched as well!
# We use the master number as the copy argument to create a slave edge.
periodic.Append ( ["line", pnums[0], pnums[3]], leftdomain=0, rightdomain=1, copy=lright, bc="periodic")

from ngsolve import *
mesh = Mesh(periodic.GenerateMesh(maxh=0.2))

1D¶

In one dimension you can add periodic identifications to points (This is not necessary for mesh creation of course, but lets NGSolve find the periodic vertices the same way as it finds the ones in higher dimensions). The added point numbers must be the ones returned from mesh.Add(MeshPoint). This is a manually created mesh with periodic boundaries:

from netgen.meshing import *

mesh = Mesh(dim=1)
pids = []
n = 50
for i in range(n+1):
    pids.append (mesh.Add (MeshPoint(Pnt(i/n, 0, 0))))
for i in range(n):
    mesh.Add(Element1D([pids[i],pids[i+1]],index=1))
mesh.Add (Element0D( pids[0], index=1))
mesh.Add (Element0D( pids[n], index=2))
mesh.AddPointIdentification(pids[0],pids[n],1,2)

Periodic FESpace¶

The periodic finite element space is created using the generator function Periodic. A quasi-periodic space can be created the same way if the optional parameter 'phase' is set.

Caution

The phase shifts must be given in the same order as the definition of the periodic boundaries in the geometry!

We can solve a periodic problem using the above defined mesh:

fes = Periodic(H1(mesh,order=3,dirichlet="outer"))

u,v = fes.TrialFunction(), fes.TestFunction()

a = BilinearForm(fes,symmetric=True)
a += SymbolicBFI(grad(u) * grad(v))

f = LinearForm(fes)
f += SymbolicLFI(x*v)

u = GridFunction(fes,"u")

with TaskManager():
    a.Assemble()
    f.Assemble()
    u.vec.data = a.mat.Inverse(fes.FreeDofs()) * f.vec

Draw(u)

The solution of this problem plotted in 2D and in 3D:

Caution

Since the dofs are mapped to the other side, periodic spaces are not possible on meshes with only one element in the direction of the periodic boundaries (then the elementmatrix would need to be of different size and NGSolve cannot handle that). So make sure that the meshsize is less than half your domain width!

You can download these examples here: 1D, 2D, 3D.

Function Documentation¶

ngsolve.Periodic()¶

Periodic or quasi-periodic Finite Element Spaces. The periodic fespace is a wrapper around a standard fespace with an additional dof mapping for the periodic degrees of freedom. All dofs on minion boundaries are mapped to their master dofs. Because of this, the mesh needs to be periodic. Low order fespaces are currently not supported, so methods using them will not work.

Parameters:

fespacengsolve.comp.FESpace: finite element space
phaselist of Complex = None: phase shift for quasi-periodic finite element space. The basis functions on the minion boundary are multiplied by the factor given in this list. If None (default) is given, a periodic fespace is created. The order of the list must match the order of the definition of the periodic boundaries in the mesh.
used_idnrslist of int = None: identification numbers to be made periodic if you don't want to use all periodic identifications defined in the mesh, if None (default) all available periodic identifications are used.

Keyword arguments can be:
order: int = 1: order of finite element space
complex: bool = False: Set if FESpace should be complex
dirichlet: regexpr: Regular expression string defining the dirichlet boundary. More than one boundary can be combined by the | operator, i.e.: dirichlet = 'top|right'
dirichlet_bbnd: regexpr: Regular expression string defining the dirichlet bboundary, i.e. points in 2D and edges in 3D. More than one boundary can be combined by the | operator, i.e.: dirichlet_bbnd = 'top|right'
dirichlet_bbbnd: regexpr: Regular expression string defining the dirichlet bbboundary, i.e. points in 3D. More than one boundary can be combined by the | operator, i.e.: dirichlet_bbbnd = 'top|right'
definedon: Region or regexpr: FESpace is only defined on specific Region, created with mesh.Materials('regexpr') or mesh.Boundaries('regexpr'). If given a regexpr, the region is assumed to be mesh.Materials('regexpr').
dim: int = 1: Create multi dimensional FESpace (i.e. [H1]^3)
dgjumps: bool = False: Enable discontinuous space for DG methods, this flag is needed for DG methods, since the dofs have a different coupling then and this changes the sparsity pattern of matrices.
autoupdate: bool = False: Automatically update on a change to the mesh.
low_order_space: bool = True: Generate a lowest order space together with the high-order space, needed for some preconditioners.
hoprolongation: bool = False: Create high order prolongation operators, only available for H1 and L2 on simplicial meshes
order_policy: ORDER_POLICY = ORDER_POLICY.OLDSTYLE: CONSTANT .. use the same fixed order for all elements, NODAL ..... use the same order for nodes of same shape, VARIABLE ... use an individual order for each edge, face and cell, OLDSTYLE .. as it used to be for the last decade
print: bool = False: (historic) print some output into file set by 'SetTestoutFile'