[Trilinos-Users] Unexpected Error In Amesos2 After Upgrading Trilinos
Jonathan Hu
jhu at sandia.gov
Wed Dec 19 15:34:34 EST 2018
Corey,
The root issue is that the example is trying to do multigrid
preconditioning on an identity matrix. For the two smaller matrix
sizes, MueLu will not do any coarsening. Instead, it will create a
"one-level" multigrid method, which is a (trivial) direct solve. For
the largest matrix, MueLu coarsens once, throwing away all the DOFs
because they are identified as boundary DOFs. The resulting coarse grid
system is zero-by-zero, which the direct solver does not like. To see
this behavior, change the line
mueLuParameterList.set("verbosity", "none");
to
mueLuParameterList.set("verbosity", "high");
You've basically hit a weird corner case, not a bug in the direct solver.
Regards,
Jonathan
Corey A. Henderson wrote on 12/15/2018 11:42 AM:
> Any chance I could get some feedback on this error?
>
> On Sun, Dec 2, 2018 at 11:59 AM Corey A. Henderson
> <cahenderson at wisc.edu <mailto:cahenderson at wisc.edu>> wrote:
>
> Hello,
>
> I am struggling to debug an error being thrown by Amesos2 in
> v12.12.1 (Release), when using Belos & MueLu for a GMRES solve.
> The identical code run against v12.10.1 (Debug) does not throw an
> error. Furthermore, the error does not happen for small matrix
> sizes. It works for N=1000 and N=2000 but fails for N=3000.
>
> The error I get is:
> MueLu::TpetraOperator::apply : detected an exception
> /dev/shm/Trilinos-trilinos-release-12-12-1/packages/amesos2/src/Amesos2_KLU2_def.hpp:144:
> Throw number = 39
> Throw test that evaluated to true: sp_colind == nullptr
> Amesos2 Runtime Error: sp_colind returned null
> terminate called after throwing an instance of 'std::runtime_error'
> what():
> /dev/shm/Trilinos-trilinos-release-12-12-1/packages/amesos2/src/Amesos2_KLU2_def.hpp:144:
> Throw number = 39
> Throw test that evaluated to true: sp_colind == nullptr
> Amesos2 Runtime Error: sp_colind returned null
> Aborted
>
> In my simple test code I believe I am correctly setting column
> indices when I set up my matrix (A should be the identity matrix).
> I've pasted code below that reproduces the error. I'd appreciate
> help identifying where my mistake is, please.
>
> --
> Corey A. Henderson
> PhD Candidate and NSF Graduate Fellow
> Dept. of Engineering Physics
> Univ. of Wisconsin - Madison
>
> ********
>
> #include <iostream>
> #include "BelosTpetraAdapter.hpp"
> #include "Tpetra_CrsMatrix.hpp"
> #include "MueLu_CreateTpetraPreconditioner.hpp"
> #include "BelosSolverFactory.hpp"
> const bool testIt(const ulong numEntries) {
> bool success = false;
> Teuchos::RCP<const Teuchos::Comm<int> > localComm(new
> Teuchos::SerialComm<int>());
> Teuchos::RCP<Tpetra::Map<int, long long> > testMap(
> new Tpetra::Map<int, long long>(static_cast<long
> long>(numEntries), 0LL, localComm, Tpetra::GloballyDistributed));
> Teuchos::RCP<Tpetra::CrsMatrix<double, int, long long> > testA(
> new Tpetra::CrsMatrix<double, int, long long>(testMap, 1UL,
> Tpetra::StaticProfile));
> Teuchos::Array<double> rowDataArray(1UL, 1.0);
> Teuchos::Array<long long> columnIndicesArray(1UL, 1LL);
> for (long long row=0; row < static_cast<long long>(numEntries);
> ++row) {
> columnIndicesArray[0] = row;
> testA->insertGlobalValues(row, columnIndicesArray(),
> rowDataArray());
> }
> testA->fillComplete();
> Teuchos::RCP<Tpetra::MultiVector<double, int, long long> > testB(
> new Tpetra::MultiVector<double, int, long long>(testMap, 1,
> false));
> for (long long row=0; row < static_cast<long long>(numEntries);
> ++row)
> testB->replaceGlobalValue(row, 0, static_cast<double>(row));
> Teuchos::RCP<Tpetra::MultiVector<double, int, long long> > testX(
> new Tpetra::MultiVector<double, int, long long>(testMap, 1,
> false));
> testX->putScalar(0.0);
> Teuchos::RCP<Tpetra::Operator<double, int, long long> >
> testAOperator =
> Teuchos::rcp_dynamic_cast<Tpetra::Operator<double, int, long
> long> >(testA, true);
> Teuchos::ParameterList mueLuParameterList;
> mueLuParameterList.set("verbosity", "none");
> mueLuParameterList.set("problem: type", "Poisson-2D");
> Teuchos::RCP<MueLu::TpetraOperator<double, int, long long> >
> mueLuPreconditioner =
> MueLu::CreateTpetraPreconditioner(testAOperator,
> mueLuParameterList);
> Belos::SolverFactory<double, Tpetra::MultiVector<double, int,
> long long>, Tpetra::Operator<double, int, long long> > factory;
> Teuchos::RCP<Teuchos::ParameterList> solverParams =
> Teuchos::parameterList();
> solverParams->set ("Num Blocks", 100);
> solverParams->set ("Maximum Iterations", 1000);
> solverParams->set ("Convergence Tolerance", 10e-6);
> Teuchos::RCP< Belos::SolverManager<double,
> Tpetra::MultiVector<double, int, long long>,
> Tpetra::Operator<double, int, long long> > >
> testSolver = factory.create("GMRES", solverParams);
> Teuchos::RCP<Belos::LinearProblem<double,
> Tpetra::MultiVector<double, int, long long>,
> Tpetra::Operator<double, int, long long> > > problem
> = Teuchos::rcp(new Belos::LinearProblem<double,
> Tpetra::MultiVector<double, int, long long>,
> Tpetra::Operator<double, int, long long> >(testAOperator, testX,
> testB));
> problem->setLeftPrec(mueLuPreconditioner);
> success = problem->setProblem();
> if (success) testSolver->setProblem(problem);
> if (success) {
> Belos::ReturnType result = testSolver->solve();
> success = (result == Belos::Converged);
> }
> if (success) {
> for (long long int row=0; row<static_cast<long long
> int>(numEntries); ++row) {
> double thisValue = testX->getData(0)[row];
> long long int longValue = round(thisValue);
> if (longValue != row) {
> success = false;
> }
> }
> }
> return success;
> }
> int main(int argc, char **argv) {
> bool result = false;
> result = testIt(1000); //works
> std::cout << "Result for 1000: " << (result ? "True" :
> "False") << std::endl;
> result = testIt(2000); //works
> std::cout << "Result for 2000: " << (result ? "True" :
> "False") << std::endl;
> result = testIt(3000); //fails
> std::cout << "Result for 3000: " << (result ? "True" :
> "False") << std::endl;
>
> return 0;
> }
>
>
>
> --
> Corey A. Henderson
> PhD Candidate and NSF Graduate Fellow
> Dept. of Engineering Physics
> Univ. of Wisconsin - Madison
>
>
> _______________________________________________
> Trilinos-Users mailing list
> Trilinos-Users at trilinos.org
> https://trilinos.org/mailman/listinfo/trilinos-users
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://trilinos.org/pipermail/trilinos-users/attachments/20181219/74bce6ca/attachment.html>
More information about the Trilinos-Users
mailing list