Differences between parts of mixed FEM and single FEM

Asked by Jan on 2013-02-05

I am wondering, why I get different sizes and entries in the system matrices when I assemble them
1.) as a submatrix of a mixed FEM ansatz
2.) directly using the same single test/trial spaces

Or do I something fundamentally wrong.

Please consider the following example in python.
For example, the length of the right hand side vector differs when changing the parameter mixed = True, False

from dolfin import *
from scipy.sparse import csr_matrix
import numpy as np

parameters.linear_algebra_backend = "uBLAS"

mixed = False

mesh = UnitSquareMesh(2, 2)

if mixed is True:
 # Define mixed FEM function spaces
 V = VectorFunctionSpace(mesh, "CG", 2)
 Q = FunctionSpace(mesh, "CG", 1)
 W = V * Q

 (u, p) = TrialFunctions(W)
 (v, q) = TestFunctions(W)

else:
 # Define FEM function spaces
 V = VectorFunctionSpace(mesh, "CG", 2)

 u = TrialFunction(V)
 v = TestFunction(V)

# define Form
aa = inner(grad(u), grad(v))*dx

# Assemble system
A = assemble(aa)

# righthandside
fvhomo = Constant((1,2))
Lvh = inner(fvhomo,v)*dx
b = assemble(Lvh)

## Convert DOLFIN representation to numpy arrays
rows, cols, values = A.data()
Aa = csr_matrix((values, cols, rows))
ba = b.array()
ba = ba.reshape(len(ba), 1)

print ba.shape

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
Jan Blechta
Solved:
2013-02-07
Last query:
2013-02-07
Last reply:
2013-02-06
Best Jan Blechta (blechta) said : #1

V.dim()=50
Q.dim()=9
W.dim()=59

If you assemble forms on space V you got matrix shaped (50, 50) and vector shaped 50. Conversely if you assemble on space W you got (59, 59) and 59 respectively. Even if you use test functions of only first subspace. So the entries in matrix and vector corresponding to V.sub(1) are zero making it singular.

If your intent is assembling the problem blockwise consider using cbc.block package.

Jan

Jan Blechta (blechta) said : #2

Typo fix: there should be W.sub(1) instead of V.sub(1).

Jan (jnhlnd) said : #3

Thanks Jan Blechta, that solved my question.