# Solve Surface PDE on stl geometry with boundary

Hi all i like to do the surface pde example but with a geometry that comes from some stl file.

is this possible ?

i would be happy if some one could tell me I have now combined the surface pde example with the, manual mesh generation example, i think thats a good startingpoint but i do not know how to set dirichlet bnd conditions, would be great if some one can help me :).

I get in the console “used dof inconsistency” … if i run the following code …

``````from ngsolve import *
from netgen.csg import *
from netgen.geom2d import unit_square
from ngsolve import Draw

ngmesh = netgen.meshing.Mesh(dim=2)

N=5
pnums = []
for i in range(N + 1):
for j in range(N + 1):
pnums.append(ngmesh.Add(netgen.meshing.MeshPoint(Pnt(i / N, j / N, 0))))

for j in range(N):
for i in range(N):
ngmesh.Add(netgen.meshing.Element2D(idx_dom, [pnums[i + j * (N + 1)],
pnums[i + (j + 1) * (N + 1)],
pnums[i + 1 + (j + 1) * (N + 1)],
pnums[i + 1 + j * (N + 1)]]))

# horizontal boundaries
for i in range(N):
ngmesh.Add(netgen.meshing.Element1D([pnums[N + i * (N + 1)],
pnums[N + (i + 1) * (N + 1)]], index=1))
ngmesh.Add(netgen.meshing.Element1D([pnums[0 + i * (N + 1)], pnums[0 + (i + 1) * (N + 1)]], index=1))

# vertical boundaries
for i in range(N):
ngmesh.Add(netgen.meshing.Element1D([pnums[i + N * (N + 1)], pnums[i + 1 + N * (N + 1)]], index=2))

mesh = Mesh(ngmesh)
Draw(mesh)

# geo          = CSGeometry()
# sphere       = Sphere(Pnt(0,0,0), 1)
# bot          = Plane(Pnt(0,0,0), Vec(0,0,-1))
# finitesphere = sphere * bot
#
# geo.NameEdge(sphere,bot, "bottom")
#
# mesh = Mesh(geo.GenerateMesh(maxh=0.3))
# mesh.Curve(2)
# Draw(mesh)
#
fes = H1(mesh, order=2, dirichlet=[1, 2])
u, v = fes.TnT()
print(fes.FreeDofs())

a = BilinearForm(fes, symmetric=True)
a.Assemble()

force = 1.0

f = LinearForm(fes)
f += force*v*ds
f.Assemble()

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

Hi Kai,

with ngmesh = netgen.meshing.Mesh(dim=2) you are creating a 2D mesh, not a surface mesh in 3D.
Thus, with a += grad(u).Trace()*grad(v).Trace()*ds you are integrating only over the four edges. This gives the used_dof_inconsistence warning.

I changed/added the following lines, compare attached file, to get the desired result. As stated in the surface pde tutorial you need the keyword dirichlet_bbnd instead of dirichlet to set Dirichlet boundary conditions for surface PDEs.

[code]ngmesh = netgen.meshing.Mesh(dim=3)

# to access edges easier

ngmesh.SetCD2Name(1, “top_bot”)
ngmesh.SetCD2Name(2, “left_right”)

fes = H1(mesh, order=2, dirichlet_bbnd=“top_bot|left_right”)[/code]
https://ngsolve.org/media/kunena/attachments/889/surface_pde.py

Best,
Michael

Attachment: 70741d8fdcc071269c78e9fa9db83d65