[petsc-users] Direct solver for 3D-FEM

[Jed Brown]

On Fri, May 20, 2011 at 15:16, Thomas Witkowski <
thomas.witkowski at tu-dresden.de> wrote:

> Thanks for the information. One related question: The system that is solved
> here is a 3x3 block matrix: Two of these blocks are empty, two are the
> discretized (FEM) identity, three are the discretized Laplace operator


What equation? Have you assumed special boundary conditions? Problems like
elasticity and Stokes normally need a symmetric stress tensor. Only for
constant coefficients and certain boundary conditions, is it okay to use the
"block Laplacian" form.


> and one block is the discretization of a complex sum of operators. At all,
> there are only three different blocks. Can I make somehow advantages of this
> structure? I could create the three inverses of the three different blocks,
> but exists than a cheap/simple way to produce the action of the inverse of
> the whole 3x3 block matrix?


For the special case of a Kronecker product, you can solve using solves with
the blocks, but not in the general case. Consider the structured grid 2D
Laplacian [A, -I; -I, A] where A is the tridiagonal problem for a 1D
Laplacian plus 2*I on the diagonal. If you could build a direct solver for
the 2D using direct solvers for the 1D problem, then you would have an O(n)
algorithm for a direct solve of a 2D Laplacian.

You can, however, build an iterative method from the blocks. If each block
is built separately, you can (with petsc-dev) put them into MatNest and use
PCFieldSplit with direct solves in splits.

An alternative is to interlace the degrees of freedom and use the (S)BAIJ
matrix format. This will reduce the amount of memory needed for your matrix
(by using a more efficient storage format) and will improve the memory
performance. Unfortunately, there are no third-party direct solvers for
BAIJ. Note that your assembly code can always call MatSetValuesBlocked() (or
..BlockedLocal()) regardless of the matrix format, so you can switch between
AIJ and BAIJ at run-time.

Why bother with a direct solver if your problem has constant coefficients?
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20110520/7a01f285/attachment-0001.htm>