- marked as major
- edited description
system_assembler(): Bilinear and linear forms do not have same exterior facet subdomains in SystemAssembler.
Issue #257
duplicate
I'm getting warnings of the form
*** Warning: Bilinear and linear forms do not have same exterior facet subdomains in SystemAssembler. Taking subdomains from bilinear form
which I had overlooked in the deluge that is the application output. At other points in FEniCS, such inconsistencies are treated as errors (e.g., inconsistencies in dx
). Shouldn't the same be the case here? (At least I would have profited from it.)
Here is an MWE that highlights how separate assembling of lhs
and rhs
works correctly, assemble_system()
does not:
from dolfin import *
mesh = RectangleMesh(0.0, 0.0, 1.0, 1.0, 100, 100)
class Obstacle(SubDomain):
def inside(self, x, on_boundary):
return (between(x[1], (0.5, 0.7)) and between(x[0], (0.5, 1.0)))
obstacle = Obstacle()
subdomains = CellFunction('size_t', mesh)
subdomains.set_all(0)
obstacle.mark(subdomains, 1)
V = FunctionSpace(mesh, 'CG', 1)
f = project(Constant(1.0), V)
submesh = SubMesh(mesh, subdomains, 1)
Vsub = FunctionSpace(submesh, 'CG', 1)
v = TestFunction(Vsub)
class Left(SubDomain):
def inside(self, x, on_boundary):
return on_boundary \
and x[0] < 0.5 + DOLFIN_EPS
left = Left()
class Right(SubDomain):
def inside(self, x, on_boundary):
return on_boundary \
and x[0] > 1.0 - DOLFIN_EPS
right = Right()
submesh_boundaries = FacetFunction('size_t', submesh)
submesh_boundaries.set_all(0)
right.mark(submesh_boundaries, 1)
ds_subdomain = Measure('ds')[submesh_boundaries]
#ds_subdomain = Measure('ds', domain=submesh)[submesh_boundaries]
#ds_subdomain = Measure('ds', domain=submesh, domain_data=submesh_boundaries)
u = TrialFunction(Vsub)
v = TestFunction(Vsub)
n = FacetNormal(submesh)
alpha = 100.0
beta = 100.0
F = dot(grad(u), grad(v)) * dx \
- Constant(alpha) * v * ds_subdomain(0) \
- Constant(beta) * v * ds_subdomain(1) \
- f * v * dx
bc = DirichletBC(Vsub, 0.0, left)
# works:
#A = assemble(lhs(F))
#b = assemble(rhs(F))
# yields warning and produces the wrong result:
A, b = assemble_system(lhs(F), rhs(F))
bc.apply(A, b)
sol = Function(Vsub)
solver = KrylovSolver('cg', 'amg')
solver.parameters['relative_tolerance'] = 1.0e-12
solver.parameters['absolute_tolerance'] = 0.0
solver.solve(A, sol.vector(), b)
plot(sol)
interactive()
Comments (4)
-
reporter -
reporter - edited description
-
reporter - changed title to system_assembler(): Bilinear and linear forms do not have same exterior facet subdomains in SystemAssembler.
-
- changed status to duplicate
Duplicate of
#78. - Log in to comment