Create Block Preconditioner that Requires Matrix Multiplication

Hi.
The attached file solves the [color=blue]mixed formulation of the diffusion equation[/color] as a [color=blue]block system[/color]. I would like to create a specific preconditioner:

[center]BlockMatrix( [ [D, None], [None, M] ] )[/center]

in which D = diag(a.mat) and M is the symmetric Gauss-Seidel (GS) matrix for

[center]Y = b.mat @ D[sup]-1[/sup] @ b.mat.T[/center]

where @ means matrix multiplication. I know one can create [color=blue]Jacobi and GS preconditioners, but I first need to compute matrix Y. Is this possible in NGSolve, so that Y can be a BaseMatrix for GS?

I would appreciate any advice about best practices for doing linear algebra in NGSolve.

Thank you,
Barry

Attachment: block_2020-01-21.py

If you want to create a Jacobi or GS smoother, you explicitely need Y as a sparse matrix. You can create D^-1 as a Sparse matrix with SparseMatrixd.CreateFromCOO. The @ operator can perform sparse matrix multiplication.

Note that the L2 space does not have any DOFs associated with boundaries, so setting Dirichlet values for an L2 space does not do anything.

Best,
Lukas

Attachment: block_2020-01-21_2020-01-21.py