assembly of mixed spaces gives wrong results

[Doug Arnold]

Here is a simple piece of code high-lighting a bug which results in a stiffness matrix being completely wrong. It computes a bilinear form acting on a test function u from Lagrange P1 on the unit interval, and a pair of test functions (v, c) from a mixed function space consisting of piecewise constants for v and a real number for c. The bilinear form has two parts:

     a0 = u.dx(0) * v * dx # integral of u'*v over the interval
     a1 = u * c * Expression("1.-x[0]") * ds # value of u*c at the point x=0

If the mesh has n intervals, then the first bilinear form should assemble into the first n rows of the (n+1) x (n+1) stiffness matrix, and the second into the final row. The problem is that in fact the second bilinear form assembles into the wrong row. Therefore if the full matrix is assembled from a0 + a1 it has a row of all zeros at the end, so is of course singular.

I first reported this in FEniCS Q&A fenicsproject.org/qa/3218/apparent-bug-assembling-with-mixed-space-for-test-functions
for version 1.3.0. It persists in version 1.4.0 and the latest development version.

     from dolfin import *
     mesh = UnitIntervalMesh(4)
     V1 = FunctionSpace(mesh, 'CG', 1)
     V2 = FunctionSpace(mesh, 'DG', 0)
     R = FunctionSpace(mesh, 'R', 0)
     u = TrialFunction(V1)
     X = MixedFunctionSpace([V2, R])
     (v, c) = TestFunctions(X)
     # this bilinear form does not involve the real test function c,
     # so one row should be identically zero
     a0 = u.dx(0) * v * dx
     print assemble(a0).array()
     # this bilinear form only involves the real test function c
     # so the corresponding row should have a nonzero
     a1 = u * c * Expression("1 - x[0]") * ds
     print assemble(a1).array()

Finally, here is the output.

[[ 0. 0. 0. 1. -1.]
 [ 0. 1. -1. 0. 0.]
 [ 0. 0. 1. -1. 0.]
 [ 1. -1. 0. 0. 0.]
 [ 0. 0. 0. 0. 0.]]
[[ 0. 0. 0. 0. 0.]
 [ 0. 0. 0. 0. 0.]
 [ 0. 0. 0. 0. 1.]
 [ 0. 0. 0. 0. 0.]
 [ 0. 0. 0. 0. 0.]]

Note that the zero row of the first array is the 5th row, but the nonzero row of the
2nd array is the third, not the fifth row.

As an answer to the posting in Q&A MiroK added some diagnostic output that indicates the problems is in DofMapBuilder::reorder_local, and found that setting

parameters['reorder_dofs_serial'] = False

works around the problem.

[Johannes Ring (johannr)]

FEniCS no longer uses Launchpad for Questions & Answers. I have added your bug report on Bitbucket:

https://bitbucket.org/fenics-project/dolfin/issue/331/assembly-of-mixed-spaces-gives-wrong

[Johannes Ring (johannr)]

I meant: FEniCS no longer uses Launchpad for Questions & Answers or bug reports.