Hi,
I am trying to use the local solver to perform a projection similar to BDM on a velocity field. The equation system is not square globally, but it is square locally on each element. For the reduced test case below the equation system is 6x6 locally.
The documentation of LocalSolver is rather sparse and on my real code I get NaN answers, while for the simplified code below I get segfaults (with FEniCS 1.6). I may be using the LocalSolver in a wrong manner (I do not quite understand the difference between solve_local_rhs and solve_global_rhs, but with the code below both gives segfaults).
I know the matrix of the real problem is well defined locally. I have my own local solver and assembler in Python which work, but are very slow even with pypy. Have I misunderstood how to use the LocalSolver or is there a bug that should be reported?
from dolfin import *
mesh = UnitTriangleMesh()
#mesh = UnitSquareMesh(1, 1)
Vu = VectorFunctionSpace(mesh, 'DG', 1)
Vv = FunctionSpace(mesh, 'DGT', 1)
u = TrialFunction(Vu)
v = TestFunction(Vv)
n = FacetNormal(mesh)
w = Constant([1, 1])
a = dot(u, n)*v*ds
L = dot(w, n)*v*ds
for R in '+-':
a += dot(u(R), n(R))*v(R)*dS
L += dot(w(R), n(R))*v(R)*dS
A = assemble(a)
b = assemble(L)
U = Function(Vu)
if True:
ls = LocalSolver(a, L)
ls.solve_local_rhs(U)
else:
solve(A, U.vector(), b)
print U.vector().array()
