Generation of a block system of equations

Question

Hi everyone, first post here. I have been looking for a while for a system in which I can solve a mixed formulation from the induced system without augmenting:

[A B^T]=[0]
[B 0 ]=[F]

But the problem is that there are no implementations of BlockMatrix whatsoever and cbc.block seems to be used in other contexts. My expected result is the following:

a = a(o, t)
bT = bT(t, u)
b = b(o,v)
f = f(v),

then assemble matrices:

A = assemble(a)
BT = assemble(bT)
B = assemble(b)
F = assemble(f)

big matrix = blockmatrix( [[A BT; B 0]])
big rhs = blockmatrix( [[0;F]] )

and finally solve it:

solve(big matrix, U, big rhs)

Thanks in advance

Answer 1

Hi, block structured systems can be conveniently assembled and solved with cbc.block

Comments

Thanks MiroK, I have tried using cbc.block, but I have had many issues with the installation in Debian Jessie regarding the connection between dolfin and petsc4py. Do you know if there is a way of implementing this being consistent with the boundary conditions? The problem is that using assembled matrices you get many lines of zeros. For example, if the spaces H and Q (mixed space HxQ) have discrete dimensions 80 and 50, then the assembled matrix is square of 130x130. This is not the expected result if we use a bilinear form including only functions from H:

a(sigma, tau) = integral dot (divergence(sigma) , divergence(tau)) * dx

and the assembled matrix should be of 80x80. Also, implementing the boundary conditions would nullify many rows and it is unclear how to reconstruct the solution afterwards.

Regards

---

Is there a reason why you don't implement your problem with MixedFunctionSpace, e.g. mixed Poisson?

---

Could you show us some code?

---

I have been reading about this, and after being able to partially execute the file mixedpoisson.py in the cbc.block demos, I reached the following difficulties that I haven't been able to overcome:

• Dolfin does not detect petsc4py, so not everything workd and I haven't been able to find a workaround for it

• I need to create some spaces that I do not know how to create, all from the PEERS elements. I am missing the extended RT with a bubble element, which I try to create through the following:

S = FunctionSpace(unitcirclemesh, 'RT', 1)
S = VectorFunctionSpace(same_mesh, 'B',1) + S

but throws an error (Exception: No point in creating empty RestrictedElement.). Same with the L2 base for skew symmetric tensors, if I just take the skew part, the system wont be invertible

• Because of the petsc4py issue, I cant use preconditioners and the routine I mentioned shows a warning of non convergence, and as saddle points are harder, even if I get to execute the code, it might not work.

I appreciate all help, as I have not been able these problems, and I have felt very comfortable with fenics in any other context, just the mixed approach through the [A B;BT 0] system seems to be tricky.

Thanks in advance

---

What does "Dolfin does not detect petsc4py" mean?

To enable petsc4py you would need to have it installed on your
system. You may have PETSc without petsc4py. If you have petsc4py
on your system the make sure it is also in PYTHONPATH
or else it might not be discovered.

In summary, does
from petsc4py import *
work ?