Hello again,
I have a function that is defined only on the surface of my mesh. In order to "extrude" the function into the interior of the mesh, so that all layers are equal, I solve the PDE
$$ \frac{\partial v}{\partial d} = 0, \hspace{5mm} v|_b = u(b), $$
where $d$ is the dimension in which we want to extrude the function $u$ and $b$ is the boundary of the mesh to where the function $u$ is defined.
The function which solves this problem works perfectly with one processor, but with two there is a zero boundary condition applied to the interface between processor sections. The following code run with mpirun will demonstrate this issue:
from dolfin import *
mesh = BoxMesh(-1, -1, 0, 1, 1, 2, 10, 10, 10)
ff = FacetFunction("size_t", mesh, 0)
# calculate boundaries :
for f in facets(mesh):
n = f.normal()
tol = 1e-3
# surface == 1 :
if n.z() >= tol and f.exterior():
ff[f] = 1
# base == 2 :
elif n.z() <= -tol and f.exterior():
ff[f] = 2
# define function space :
Q = FunctionSpace(mesh, "CG", 1)
u = Expression('sqrt(pow(x[0],2) + pow(x[1], 2) + pow(x[2], 2))')
def extrude(f, b, d, ff, Q):
# define test and trial based on function space :
phi = TestFunction(Q)
v = TrialFunction(Q)
# linear PDE :
a = v.dx(d) * phi * dx
L = DOLFIN_EPS * phi * dx # really close to zero to fool FFC
bc = DirichletBC(Q, f, ff, b)
# solve and return new Function
v = Function(Q)
solve(a == L, v, bc)
return v
def vert_integrate(u, ff, Q):
phi = TestFunction(Q)
v = TrialFunction(Q)
bc = DirichletBC(Q, 0, ff, 2)
a = v.dx(2) * phi * dx
L = u * phi * dx
v = Function(Q)
solve(a == L, v, bc)
return v
# plot solution :
test = vert_integrate(u, ff, Q) # vertical integral of u
ex = extrude(test, 1, 2, ff, Q) # extruded surface of test
File('u.pvd') << project(u, Q)
File('test.pvd') << test
File('ex.pvd') << ex
#plot(test)
#plot(ex)
#interactive()
The above code, when run with
mpirun -np 2 python code.py
provides this output:
