Add SolverDirect for Tpetra Wrappers

[jpthiele]

This will add support for (sparse) direct solvers implemented
or interfaced in Amesos2, which is part of the Tpetra stack.

Status:

Implementation itself (with A as LA::TpetraWrappers::SparseMatrix)

• Interface for initialize(A)
• Interface for solve(x,b) with LA::TpetraWrappers::Vector
• Interface for solve(A,x,b) with LA::TpetraWrappers::Vector
• Interface for setting Teuchos::ParameterList

[bangerth]

Might help with #10414.

[bangerth]

You'll eventually also have to add the .cc and/or .templates.h file.

[jpthiele]

You'll eventually also have to add the .cc and/or .templates.h file.

That file is already in the works.

I also looked at the Epetra tests for the direct solvers but they CG solver of AztecOO as a baseline.
I would suggest using KLU2 and compare it to the solutions of all other solvers that are configured?
Maybe we will also have to configure a test setup for the currently uninstalled solvers.

[jpthiele]

Might help with #10414.

Yes definitely, there are multiple solvers that were at least not mentioned in the Epetra wrappers but partially also not in Amesos(1).

SuperLU_MT is also wrapped in Amesos2.
As is MKL_Pardiso which has an OpenMP version.
And ShyLU Basker and Tacho which are a shared-memory LU and Cholesky factorizations.

But looking at the documentation of the solver wrapper classes in Amesos2 I see that we should also give the option
to pass a Teuchos::Parameter list, so users can fine tune the solvers should they want to.

[jpthiele]

Related to #16614 the question remains what the best way to convert a sequential Vector and a distributed Vector from deal.II to Tpetra actually is.
As it is also needed here a specialized function would be great @kinnewig.

[jpthiele]

Working on the parameters a question would be whether we want some Utility classes or structs
to handle creation of a valid parameter list for a solver?

Especially for MUMPS a user would have to dig deep to find out what ICNTL(1) or ICNTL(9) would actually set.
Having a simple struct maybe with enums to set things like solving Ax=b vs. A^Tx=b vs. A^*x=b or whether or not a matrix is symmetric would be a lot easier then.

[bangerth]

We usually do this via AdditionalData structs. That could be used here as well, but I will be ok with just a patch that provides the bare minimum to make things work, and then come back to making it more convenient later on.

I think it would be great if we could make the majority of the tests in #16611 work by the end of April or so, so we can be sure that this can all be shipped with the next release. There will be things we want to improve later on, and making things convenient may fall in this category.

[jpthiele]

The additional data would work if there is a nonempty intersection of all possible parameters for the solvers.

Another option would be to write derived classes for each of the supported solvers with their data sets.
But that seems like doubling work as that is already what Amesos2 does basically.

[jpthiele]

In terms of staying in one ecosystem, i.e. all TpetraWrappers the SolverDirect wrapper is now feature complete.
All we need is a mapping from dealii::Vector and dealii::LinearAlgebra::distributed::Vector but that could also be done
in a seperate PR.

Considering testing we might think about also testing the float, std::complex<float> and std::complex<double> configurations.

[jpthiele]

I just had a look at the Ifpack1 preconditioners in deal.II and noticed that we basically do
exactly what I proposed, i.e. a PreconditionBase that could be used with any ifpack preconditioner given the right parameter set
(and we should probably expose that functionality somewhere) for Ifpack2
and we have derived classes to set the right set of parameters.

So I would propose to actually use the chance of reworking the interface to have the same structure for both.

The structure would be as follows:

• PreconditionBase and SolverDirectBase
• derived classes for specific solvers, e.g. PreconditionsSSOR, SolverDirectKLU(2), SolverDirectMUMPS.
• and a derived class for full customization using the Teuchos::ParameterList, SolverDirect and Precondition

What do you think @bangerth, @kinnewig, @kronbichler, @masterleinad ?

[kronbichler]

I would be very happy if we could have a class that can be set with parameter files alongside some generic classes that implement often-used options, exactly as you suggest. I think it often is desirable for more complex codes to switch between different options for solvers and preconditioners via a parameter file. Of course, one can address that by adding the necessary switches in the C++ code (via a base class or similar mechanisms) for cases where the code needs different implementations. But since Trilinos/Ifpack2/Amesos2 already does control via a single solver instance that for us, it would be a pity not to use it. So 👍 from my end to this suggestion.

[jpthiele]

I have reworked the class layout as discussed.
So now we have a SolverDirectBase that can't be instantiated on its own (protected constructor)
a SolverDirect that works basically like the intermediate implementation by supplying a solver_name
in the AdditionalData field and then optionally passing a Teuchos::ParameterList to fine tune.

Additionally, we have SolverDirektKLU2 which has all the AdditionalData fields that the parameter set provides.

We also have a test 01 to ensure setting a parameter sets works (here with KLU2 again) and to ensure that
the exception is raised when an unsupported solver is about to be created.
And a test 02 to see that KLU2 works.

I have also finished a SolverDirectLAPACK and SolverDirectTachoderived class,
but both fail with some Lapack error

** On entry to DPOTRF parameter number  4 had an illegal value (Tacho)
** On entry to DGEEQU parameter number  4 had an illegal value (Lapack)

but that might possibly stem from some 64bit vs. 32bit clash.
I'd have to recompile Trilinos with 32bit to find out.
If anyone has a setup like that already running you can try to see if you get these errors by passing "Tacho" or "LAPACK" to the SolverDirect class as solver_name

[jpthiele]

I have undrafted it because I feel that any further work like supporting specific solvers could be done in separate PRs.

[jpthiele]

Having the same problem with Ifpack2, I just noticed that we should check whether Amesos2 is actually part of Trilinos and added that to the CMake files.
This works, but changing the filenames of the direct tests and inserting trilinos_with_amesos2=on
resulted in the tests not being picked up by the setup_tests_trilinos_tpetra target.
So is there anything more to do to let the test setup know to pickup that new CMake variable?

[bangerth]

/rebuild

[jpthiele]

Okay figures out how to ensure amesos2 is set up.
The test filename needed an additional 'with', i.e. with_trilinos_with_amesos2=true.

[bangerth]

Looks good, just a few small issues.

[bangerth]

Each of the main classes in this file needs a brief documentation describing what it does.

[bangerth]

Same here and elsewhere.

[jpthiele]

I am not sure if this relates to the three lines rule or if you meant anything else?

[bangerth]

This related to a comment that must have gone lost. I don't think the verbose_cout thing is very useful. Let's use the opportunity to get rid of it here and in all of the other places where we use it in these direct solvers, along with the flag that controls whether to output.

[bangerth]

Our usual convention is 3 empty lines between functions. This whole file doesn't satisfy this.

[jpthiele]

Done

[bangerth]

This is purposefully a non-existent solver, right? Perhaps add a comment in front of it?

[jpthiele]

Done.

[jpthiele]

If have added the missing documentation for the classes and structs.
If any formulation is difficult to understand or could be improved I'm happy to fix it.

[bangerth]

This looks good with the exception of the output thing. Let's merge since it is ready; if you feel like addressing the output thing in a separate PR, go ahead -- or we can leave that as a follow-up and I can just open a github issue.

[bangerth]

This related to a comment that must have gone lost. I don't think the verbose_cout thing is very useful. Let's use the opportunity to get rid of it here and in all of the other places where we use it in these direct solvers, along with the flag that controls whether to output.

[jpthiele]

Ahh yes that must have gone missing.
I've removed it in a new branch and when it is done compiling and testing I will open a PR.