MueLu: block-iterative AMG for multiphysics

mayrmt commented 3 years ago

@trilinos/muelu @jhux2 @cgcgcg

Current state of multiphysics preconditioners in MueLu

The work by @tawiesn (now picked up by @pohm01) uses MueLu (and Xpetra) to create a family of "fully coupled AMG preconditioners for NxN multiphysics problems" (N = number of block rows in the matrix). "Fully coupled" refers to the fact, that the off-diagonal coupling blocks are transferred to the coarse levels and are then taken into account on each level by block smoothers implemented in MueLu and/or Xpetra.
RefMaxwell provides a specialized preconditioner tailored to a specific formulation of Maxwell's equations (2x2 block system)

Expectation

I want to perform a "block-iterative AMG for multiphysics". By that, I mean a preconditioner, that addresses the off-diagonal couplings on the fine level only, but allows to use different MueLu hierarchies to approximate the inverses of each main diagonal block (whenever required by the fine-level block method). This also includes the special case, where some (small) blocks are not tackled via multigrid, but just by any regular smoother or a direct solver.

Differences to existing implementations

The building blocks are exactly the same as i the "fully coupled AMG", however the nesting is reversed.
While RefMaxwell (as I understand it) follows the same idea, but it is specialized to 2x2 block systems arising from Maxwell's equations (cf. P. Bochev, J. Hu, C. Siefert, and R. Tuminaro. "An algebraic multigrid approach based on a compatible gauge reformulation of Maxwell's equations." SIAM Journal on Scientific Computing, 31(1), 557-583).

❗ Please let me know, if I have just overlooked an implementation in MueLu that already satisfies my needs.

Possible applications

I need this right now for a 2x2 system arising from embedding fibers into solid volumes (i.e. totally different physics than Maxwell), so I can't re-use the existing RefMaxwell implementation. Other possible applications for 2x2 systems can be found in mortar meshtying or contact mechanics (https://doi.org/10.1002/nme.6680). Of course, the idea is not limited to 2x2 systems, but has already been showcased in 3x3 systems (e.g. FSI, https://doi.org/10.1002/nme.3001).

Possible solutions

I'll outline my idea. However, that's the part where I hope for input from other @trilinos/muelu developers.

For the outer block iteration, use a 1-level hierarchy and use the desired block iteration method as coarse solver. Then, mask MueLu hierarchies as smoother objects and use them to approximate the required block inverses. Besides wrapping MueLu hierarchies as smoothers, this approach might not require much coding, but rather some advanced xml files.

Of course, there's the danger of circular dependencies. Not sure, if this is actually an issue and if and how it can be resolved.

cgcgcg commented 3 years ago

@mayrmt Do you want to coarsen the blocks together? Or would Teko work?

mayrmt commented 3 years ago

@cgcgcg I want to coarsen the blocks individually. The multigrid hierarchy will never see any coupling information.

Teko might work, but it seems to be tailored to 2x2 flow problems, right? I need a generic framework that allows for different numbers of blocks (2x2, 3x3, 4x4 for me). From a application point of view, it would be great if it lived just inside MueLu, so that one doesn't need to change the call to MueLu, just pass a different xml file.

cgcgcg commented 3 years ago

@mayrmt Teko should work independent of the number of blocks. It is possible that it has mostly been used for 2x2.

tawiesn commented 3 years ago

@cgcgcg @mayrmt Another issue of Teko is that it is based on Thyra. Matthias proposal would be based on the Xpetra stack that often is easier to use... With the RefMaxwell, it has already been demonstrated that one can have a special preconditioner for a specific application. One could have something similar with an outer BGS or SIMPLE (or whatever) iteration (based on Xpetra blocks) with MueLu subsmoothers.

cgcgcg commented 3 years ago

RefMaxwell is a bit of a special case. The way it is set up, we are preconditioning a single matrix, curlcurl+mass. When we are solving the Maxwell system in 2x2 block form, as we do in MiniEM, we build a Teko 2x2 preconditioner with RefMaxwell in one of the sub-blocks.

jhux2 commented 3 years ago

@mayrmt Just to recapitulate: what you are seeking is a way of solving blocked NxN linear systems, where N > 2, and currently you are solving for N=2. AMG would be used for inverting single blocks, and coupling between blocks is accounted for in some other fashion.

My preference would be see whether Teko is sufficient for your needs. If it isn't, I think it would be worthwhile to see where it's lacking. There are some internal applications which are beginning to use Teko, and it may be possible to leverage that work.

mayrmt commented 3 years ago

@jhux2 wrote:

Just to recapitulate: what you are seeking is a way of solving blocked NxN linear systems, where N > 2, and currently you are solving for N=2. AMG would be used for inverting single blocks, and coupling between blocks is accounted for in some other fashion.

Yes, that summarizes my intentions very well.

As of now, I'd like to pursue a pure "MueLu+Xpetra" implementation for the following reasons:

A pure "MueLu+Xpetra" implementation allows for a fair and direct comparison of the "Multigrid for Multiphysics" preconditioners developed by @tawiesn and picked-up by @pohm01 (i.e. with multiphysics coupling on all levels) vs. the proposed block-iterative approach. Our prior work indicates that the fully coupled variant is superior to the block-iterative approach in terms of iteration counts, but is often more expensive in terms of timings (especially during setup). To study this in detail and to enable a fair comparison, the "MueLu+Xpetra" route seems the best fit, as it actually uses exactly the same pieces of software -- just with a inverted order of nesting them into each other.
Doing block-iterative with "MueLu+Xpetra" would allow to use the same call to create the preconditioner, since all the magic would happen inside the MueLu xml-file.
Teko has required dependencies on Thyra (which we wouldn't need otherwise) and some more required dependencies on some of the Epetra-stack packages like ML/Aztec/Ifpack.
Teko uses Thyra-style numbering of GIDs, i.e. GID numbering starts from zero in each block. Some applications (e.g. ours) do not support that natively, but rather use Xpetra-style numbering, i.e. continuous GIDs across all blocks in a block system.
Our applications currently do not use Teko, just MueLu and Xpetra plus either Ifpack/Ifpack2 and Aztec/Belos. Using Teko would require major changes to our Trilinos configuration.

I have done some first prototypes locally to work towards the block-iterative approach in "MueLu+Xpetra". They should work without major changes to MueLu. (I have created some PRs, e.g. #9631, along the way, which are not particularly related to this issue, but just address some outstanding todos as a byproduct of my experiments.)

mayrmt commented 2 years ago

Pinging @maxfirmbach

github-actions[bot] commented 1 year ago

This issue has had no activity for 365 days and is marked for closure. It will be closed after an additional 30 days of inactivity. If you would like to keep this issue open please add a comment and/or remove the MARKED_FOR_CLOSURE label. If this issue should be kept open even with no activity beyond the time limits you can add the label DO_NOT_AUTOCLOSE. If it is ok for this issue to be closed, feel free to go ahead and close it. Please do not add any comments or change any labels or otherwise touch this issue unless your intention is to reset the inactivity counter for an additional year.

mayrmt commented 1 year ago

@maxfirmbach is working on this, so do not auto-close, yet.

github-actions[bot] commented 11 months ago

This issue has had no activity for 365 days and is marked for closure. It will be closed after an additional 30 days of inactivity. If you would like to keep this issue open please add a comment and/or remove the MARKED_FOR_CLOSURE label. If this issue should be kept open even with no activity beyond the time limits you can add the label DO_NOT_AUTOCLOSE. If it is ok for this issue to be closed, feel free to go ahead and close it. Please do not add any comments or change any labels or otherwise touch this issue unless your intention is to reset the inactivity counter for an additional year.

mayrmt commented 11 months ago

@maxfirmbach is working on this, so do not auto-close, yet.

maxfirmbach commented 8 months ago

@mayrmt and I decided to proceed our work regarding block-iterative preconditioners by using Teko and using MueLu for the inversion of sub-blocks (single/blocked). We are still exploring the functionalities and features of Teko and need to figure things out, but the 3x3 FSI case from https://doi.org/10.1002/nme.3001 is already working. Let's see if some questions arise in the near future :sweat_smile:

maxfirmbach commented 1 month ago

We now use Teko's block Gauss-Seidel for a handful of examples in our code, which works without any problems so far. Let's keep this open for a little bit longer to have a look at the long term behavior. We still need to test the saddle-point methods e.g. SIMPLE and the block LU.

trilinos / Trilinos