I would like to get a faster inversion of a local block-diagonal system, which I am using for time-stepping.
The bilinear form is
a += (grad(u)*grad(v)+lam*v+u*mu)*dx
where u, v belong to DG-Pk space, and lam and mu belong to DG-P0 space.
The matrix is symmetric and block diagonal. What is the fastest way to invert this matrix?
The following direct approach is a bit to slow
inva = a.mat.Inverse(fes.FreeDofs())
Attached is the code.
Here I found a solution, might be 3-4 times faster than sparse-direct solver:
First create a block Jacobi preconditioner (which is actually the exact solver), and then apply one-step smoothing. This seems to be the fastest way that I can think of…
Yes, I matrix-vector multiply with the block-Jacobi should be the fastest choice, it uses dense matrix-vector for the small blocks, and is parallel.
But double-check whether the geom_free for your bilinear-form is really doing what you want.
The geom_free flag is very mysterious to me. I may have to post another question on the subject. As for the current matrix problem, I think it is of no use, and should be turned off.
Well, the matrix-free approach looks quite promising for time domain DG simulation as was demonstrated in the i-tutorials 3.3 Discontinuous Galerkin Discretizations for linear transport — NGS-Py 6.2.2302 documentation
But, I need to restructure the code quite a bit to do geom_free…