This page was generated from unit-6.1.2-surfacepde/surface_pdes.ipynb.

6.1.2 Surface PDEs

Surface Poisson equation

Homogeneous Dirichlet data

from ngsolve import *
from netgen.csg import *
from netgen.meshing import MeshingStep
from ngsolve.webgui import Draw

Given a half-sphere we want to solve the surface Poisson equation

\[\int_S\nabla_{\Gamma}u\cdot\nabla_{\Gamma}v\,ds = \int_S fv\,ds,\qquad u=0 \text{ on } \Gamma_D\]
geo          = CSGeometry()
sphere       = Sphere(Pnt(0,0,0), 1)
bot          = Plane(Pnt(0,0,0), Vec(0,0,-1))
finitesphere = sphere * bot

geo.AddSurface(sphere, finitesphere.bc("surface"))
geo.NameEdge(sphere,bot, "bottom")

mesh = Mesh(geo.GenerateMesh(maxh=0.3))

Therefore, we define the usual \(H^1\) finite element space and use the dirichlet_bbnd flag to indicate the BBoundary on which the Dirichlet condtions are prescribed. The test- and trial-functions are gien as usual.

fes = H1(mesh, order=2, dirichlet_bbnd="bottom")
u, v = fes.TnT()
0: 00000000000000000000011111111111111111111111111111
50: 11111111111111111111111111111111111111111111111110
100: 01101101101101101101101101101101101101101101101101
150: 10110110111111111111111111111111111111111111111111
200: 11111111111111111111111111111111111111111111111111
250: 11111111111111111111111111111111111111111111111111
300: 11111111111111111111111111111111111111111111111111
350: 1111111111111111111111

For the (bi-)linear form we have to take care that we have to define the according integrators on the boundary and that the Trace operator has to be used to obtain the tangential/surface derivative

a = BilinearForm(fes, symmetric=True)
a += grad(u).Trace()*grad(v).Trace()*ds

force = sin(x)*y*exp(z)

f = LinearForm(fes)
f += force*v*ds
<ngsolve.comp.LinearForm at 0x7f3e2c6418f0>

Solving is done as usual

gfu = GridFunction(fes) = a.mat.Inverse(fes.FreeDofs())*f.vec
Draw(gfu, mesh, "u")

Inhomogeneous Dirichlet data

To solve the same problem with non-homogenous Dirichlet data, \(u=u_D\) on \(\Gamma_D\) the same technique as in the volume case is used, where we have to set a function on the BBoundary instead on the boundary

gfu.Set(x, definedon=mesh.BBoundaries("bottom"))
r = f.vec.CreateVector() = f.vec - a.mat*gfu.vec += a.mat.Inverse(fes.FreeDofs())*r

Finite element spaces for surfaces

We differ between two types of finite element spaces for surfaces. The first class consists of spaces, where the restriction of a 3D element to the surface leads to a valid 2D element of the same type: H1, HCurl, HCurlCurl, NumberSpace.

These spaces can directly be used, one has to take care using the Trace operator. Otherwise an exception is thrown during assembling.

a += grad(u)*grad(v)*ds
NgException                               Traceback (most recent call last)
<ipython-input-1-afcf4da53a20> in <module>
----> 1 a += grad(u)*grad(v)*ds
      2 a.Assemble()

NgException: Trialfunction does not support BND-forms, maybe a Trace() operator is missing, type = grad

The NumberSpace is an exception as it represents only a number, where no Trace operator has to be used.

The second class is given by


Surface Space









Here, a 2D reference element is mapped directly onto the surface. To be consistent, also here the Trace operator has to be used.

[ ]: