[Trilinos-Users] Aztec seg fault.
Kurtis Nusbaum
klnusbaum at gmail.com
Fri Jan 15 16:31:00 MST 2010
I'm doing some research comparing the performance of the Belos and Aztec
solvers. Using stratimikos I set up a problem to solve various sizes of a
laplace_3d matrix (see attached source code if you whish). For the most part
the program works fine. But when ever I try to solve problems sizes greater
than 1000000 using aztec I get a segfault. The segfault seems to occur in a
function called AZ_manage_memory. Does anyone have an idea as to what is
going on? I've attached the settings i'm using for running Aztec via
Stratimikos and the gdb output from when I run the problem via gdb.
-Kurtis Nusbaum
-------------- next part --------------
An HTML attachment was scrubbed...
URL: https://software.sandia.gov/pipermail/trilinos-users/attachments/20100115/15e0bed2/attachment-0001.html
-------------- next part --------------
(gdb) run --linear-solver-params-file="AztecSettings.xml" --problem-size=1000000
Starting program: /net/home/f07/klnusbau/sandbox/BelosAztecTesting/belosaztecTester --linear-solver-params-file="AztecSettings.xml" --problem-size=1000000
[Thread debugging using libthread_db enabled]
Printing statistics of the Epetra linear system ...
Epetra_CrsMatrix epetra_A of dimension 1000000 x 1000000
||epetraA||inf = 12
||epetra_b||2 = 249.8
||epetra_x||2 = 0
Reading parameters from XML file "AztecSettings.xml" ...
Entering Thyra::AztecOOLinearOpWithSolveFactory::initializeOp_impl(...) ...
Entering Thyra::IfpackPreconditionerFactory::initializePrec(...) ...
Creating the initial Ifpack_Preconditioner object of type 'ILU' ...
=> Creation time = 7.7e-05 sec
Computing the factorization of the preconditioner ...
=> Factorization time = 8.78983 sec
Description of created preconditioner:
================================================================================
Ifpack_AdditiveSchwarz, overlap level = 0
Combine mode = Zero
Condition number estimate = 3.14192
Global number of rows = 1000000
Phase # calls Total Time (s) Total MFlops MFlops/s
----- ------- -------------- ------------ --------
Initialize() 1 10.4308 0 0
Compute() 1 8.27269 149.64 18.0884
ApplyInverse() 2 0.462207 172.184 372.526
================================================================================
Total time = 19.2225 sec
Leaving Thyra::IfpackPreconditionerFactory::initializePrec(...) ...
Leaving Thyra::AztecOOLinearOpWithSolveFactory::initializeOp_impl(...) ...
Solving block system using AztecOO ...
*******************************************************
***** Problem: Epetra::CrsMatrix
***** Preconditioned GMRES solution
***** Ifpack_AdditiveSchwarz, ov = 0, local solver =
***** `IFPACK ILU (fill=2, relax=0.000000, athr=0.000000, rthr=1.000000)'
***** Condition number estimate = 3.14192
***** No scaling
*******************************************************
Program received signal SIGSEGV, Segmentation fault.
0x00000000007bad4c in AZ_manage_memory (input_size=2408004816, action=0, type=-914901, name=0x7fffffffc450 "vblock in gmres0", status=0x7fffffffc4ec) at /home/f07/klnusbau/sandbox/Trilinos/packages/aztecoo/src/az_util.c:991
991 temp->name = (char *) &(dtmp[(aligned_str_mem+aligned_size)/
Missing separate debuginfos, use: debuginfo-install atlas-3.8.3-9.fc11.x86_64 blas-3.2.1-3.fc11.x86_64 glibc-2.10.2-1.x86_64 libgcc-4.4.1-2.fc11.x86_64 libgfortran-4.4.1-2.fc11.x86_64 libstdc++-4.4.1-2.fc11.x86_64
(gdb) bt
#0 0x00000000007bad4c in AZ_manage_memory (input_size=2408004816, action=0, type=-914901, name=0x7fffffffc450 "vblock in gmres0", status=0x7fffffffc4ec) at /home/f07/klnusbau/sandbox/Trilinos/packages/aztecoo/src/az_util.c:991
#1 0x00000000007df288 in AZ_pgmres (b=0x169f110, x=0x1e40430, weight=0x1411b588, options=0x1411b470, params=0x1411b540, proc_config=0x1411b3f0, status=0x1411b640, Amat=0x1411bb50, precond=0x1411cbd0, convergence_info=0x8d0a050)
at /home/f07/klnusbau/sandbox/Trilinos/packages/aztecoo/src/az_gmres.c:232
#2 0x00000000007ccc53 in AZ_oldsolve (x=0x1e40430, b=0x169f110, options=0x1411b470, params=0x1411b540, status=0x1411b640, proc_config=0x1411b3f0, Amat=0x1411bb50, precond=0x1411cbd0, scaling=0x1411d840)
at /home/f07/klnusbau/sandbox/Trilinos/packages/aztecoo/src/az_solve.c:540
#3 0x00000000007cc084 in AZ_iterate (x=0x1e40430, b=0x169f110, options=0x1411b470, params=0x1411b540, status=0x1411b640, proc_config=0x1411b3f0, Amat=0x1411bb50, precond=0x1411cbd0, scaling=0x0)
at /home/f07/klnusbau/sandbox/Trilinos/packages/aztecoo/src/az_solve.c:190
#4 0x00000000007a4a19 in AztecOO::Iterate (this=0x1411b2d0, MaxIters=400, Tolerance=1e-13) at /home/f07/klnusbau/sandbox/Trilinos/packages/aztecoo/src/AztecOO.cpp:902
#5 0x00000000005b0fba in Thyra::AztecOOLinearOpWithSolve::solve (this=0x2608070, M_trans=Thyra::NOTRANS, B=@0x25e2b60, X=0x25e2920, numBlocks=1, blockSolveCriteria=0x7fffffffd210, blockSolveStatus=0x7fffffffd1e0)
at /home/f07/klnusbau/sandbox/Trilinos/packages/stratimikos/adapters/aztecoo/src/Thyra_AztecOOLinearOpWithSolve.cpp:673
#6 0x00000000004a1bbb in Thyra::SingleScalarLinearOpWithSolveBase<double>::solve (this=0x26081b0, conj=Thyra::NONCONJ_ELE, B=@0x25e2b60, X=0x25e2920, numBlocks=1, blockSolveCriteria=0x7fffffffd210, blockSolveStatus=0x7fffffffd1e0)
at /home/f07/klnusbau/sandbox/Trilinos/packages/thyra/src/support/operator_solve/client_support/Thyra_SingleScalarLinearOpWithSolveBase.hpp:74
#7 0x000000000040f143 in Thyra::SolveStatus<Thyra::LinearOpWithSolveBase<double, double>::PromotedScalar> Thyra::solve<double, double>(Thyra::LinearOpWithSolveBase<double, double> const&, Thyra::EConj, Thyra::MultiVectorBase<double> const&, Thyra::MultiVectorBase<double>*, Thyra::SolveCriteria<Thyra::LinearOpWithSolveBase<double, double>::PromotedScalar> const*) ()
#8 0x000000000040c839 in Thyra::SolveStatus<double> Thyra::solve<double>(Thyra::LinearOpWithSolveBase<double, double> const&, Thyra::EOpTransp, Thyra::MultiVectorBase<double> const&, Thyra::MultiVectorBase<double>*, Thyra::SolveCriteria<double> const*) ()
#9 0x00000000004090cd in main ()
(gdb) list
986 (void) AZ_printf_err( " Asked for %ld bytes. Perhaps\n", size);
987 (void) AZ_printf_err( " a smaller problem should be run\n");
988 exit(-1);
989 }
990 temp = (struct mem_ptr *) &(dtmp[aligned_size/sizeof(double)]);
991 temp->name = (char *) &(dtmp[(aligned_str_mem+aligned_size)/
992 sizeof(double)]);
993 temp->address = dtmp;
994 for (i = 0 ; i < j; i++ ) (temp->name)[i] = name[i];
995
-------------- next part --------------
// @HEADER
// ***********************************************************************
//
// Stratimikos: Thyra-based strategies for linear solvers
// Copyright (2006) Sandia Corporation
//
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
//
// This library is free software; you can redistribute it and/or modify
// it under the terms of the GNU Lesser General Public License as
// published by the Free Software Foundation; either version 2.1 of the
// License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
// USA
// Questions? Contact Roscoe A. Bartlett (rabartl at sandia.gov)
//
// ***********************************************************************
// @HEADER
#include "Trilinos_Util_CrsMatrixGallery.h"
#include "Stratimikos_DefaultLinearSolverBuilder.hpp"
#include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
#include "Thyra_EpetraThyraWrappers.hpp"
#include "Thyra_EpetraLinearOp.hpp"
#include "Teuchos_VerboseObject.hpp"
#include "Teuchos_XMLParameterListHelpers.hpp"
#include "Teuchos_CommandLineProcessor.hpp"
#include "Teuchos_StandardCatchMacros.hpp"
#include "Epetra_SerialComm.h"
#include "Epetra_Vector.h"
#include "Epetra_Map.h"
#include "Epetra_LinearProblem.h"
namespace {
// Helper function to compute a single norm for a vector
double epetraNorm2( const Epetra_Vector &v )
{
double norm[1] = { -1.0 };
v.Norm2(&norm[0]);
return norm[0];
}
} // namespace
int main(int argc, char* argv[])
{
using Teuchos::rcp;
using Teuchos::RCP;
using Teuchos::CommandLineProcessor;
typedef Teuchos::ParameterList::PrintOptions PLPrintOptions;
bool success = true;
bool verbose = true;
Teuchos::RCP<Teuchos::FancyOStream>
out = Teuchos::VerboseObjectBase::getDefaultOStream();
try {
//
// A) Program setup code
//
//
// Read options from command-line
//
int problemSize = 1000;
double tol = 1e-5;
Stratimikos::DefaultLinearSolverBuilder linearSolverBuilder;
CommandLineProcessor clp(false); // Don't throw exceptions
// Set up command-line options for the linear solver that will be used!
linearSolverBuilder.setupCLP(&clp);
clp.setOption( "tol", &tol, "Tolerance to check against the scaled residual norm." );
clp.setOption( "problem-size", &problemSize, "The size of the problem");
clp.setDocString(
"Simple example for the use of the Stratimikos facade Stratimikos::DefaultLinearSolverBuilder.\n"
"\n"
"To solve a linear system from a Triutils gallery matrix"
" with options read from an XML file using --linear-solver-params-file=\"SomeFile.xml\"\n"
"Set the size of the problem using --problem-size\n"
);
CommandLineProcessor::EParseCommandLineReturn parse_return = clp.parse(argc,argv);
if( parse_return != CommandLineProcessor::PARSE_SUCCESSFUL ) return parse_return;
//
// B) Epetra-specific code that sets up the linear system to be solved
//
// While the below code reads in the Epetra objects from a file, you can
// setup the Epetra objects any way you would like. Note that this next
// set of code as nothing to do with Thyra at all, and it should not.
//
// *out << "\nReading linear system in Epetra format from the file \'"<<matrixFile<<"\' ...\n";
Epetra_SerialComm comm;
Trilinos_Util::CrsMatrixGallery matrixGallery("laplace_3d",comm);
matrixGallery.Set("problem_size", problemSize);
RCP<Epetra_CrsMatrix> epetra_A = RCP<Epetra_CrsMatrix>(matrixGallery.GetMatrix());
RCP<Epetra_Vector> epetra_x = RCP<Epetra_Vector>((*matrixGallery.GetLinearProblem()->GetLHS())(0));
RCP<Epetra_Vector> epetra_b = RCP<Epetra_Vector>((*matrixGallery.GetLinearProblem()->GetRHS())(0));
/*EpetraExt::readEpetraLinearSystem( matrixFile, comm, &epetra_A, NULL, &epetra_x, &epetra_b );
if(!epetra_b.get()) {
*out << "\nThe RHS b was not read in so generate a new random vector ...\n";*
epetra_b = rcp(new Epetra_Vector(epetra_A->OperatorRangeMap()));
epetra_b->Random();
}
if(!epetra_x.get()) {
*out << "\nThe LHS x was not read in so generate a new zero vector ...\n";
epetra_x = rcp(new Epetra_Vector(epetra_A->OperatorDomainMap()));
epetra_x->PutScalar(0.0); // Initial guess is critical!
}*/
*out << "\nPrinting statistics of the Epetra linear system ...\n";
*out
<< "\n Epetra_CrsMatrix epetra_A of dimension "
<< epetra_A->OperatorRangeMap().NumGlobalElements();
*out << " x " << epetra_A->OperatorDomainMap().NumGlobalElements()
<< "\n ||epetraA||inf = " << epetra_A->NormInf()
<< "\n ||epetra_b||2 = " << epetraNorm2(*epetra_b)
<< "\n ||epetra_x||2 = " << epetraNorm2(*epetra_x)
<< "\n";
//
// C) The "Glue" code that takes Epetra objects and wraps them as Thyra
// objects
//
// This next set of code wraps the Epetra objects that define the linear
// system to be solved as Thyra objects so that they can be passed to the
// linear solver.
//
RCP<const Thyra::LinearOpBase<double> >
A = Thyra::epetraLinearOp( epetra_A );
RCP<Thyra::VectorBase<double> >
x = Thyra::create_Vector( epetra_x, A->domain() );
RCP<const Thyra::VectorBase<double> >
b = Thyra::create_Vector( epetra_b, A->range() );
// Note that above Thyra is only interacted with in the most trival of
// ways. For most users, Thyra will only be seen as a thin wrapper that
// they need know little about in order to wrap their objects in order to
// pass them to Thyra-enabled solvers.
//
// D) Thyra-specific code for solving the linear system
//
// Note that this code has no mention of any concrete implementation and
// therefore can be used in any use case.
//
// Reading in the solver parameters from the parameters file and/or from
// the command line. This was setup by the command-line options
// set by the setupCLP(...) function above.
linearSolverBuilder.readParameters(out.get());
// Create a linear solver factory given information read from the
// parameter list.
RCP<Thyra::LinearOpWithSolveFactoryBase<double> >
lowsFactory = linearSolverBuilder.createLinearSolveStrategy("");
// Setup output stream and the verbosity level
lowsFactory->setOStream(out);
lowsFactory->setVerbLevel(Teuchos::VERB_MEDIUM);
// Create a linear solver based on the forward operator A
RCP<Thyra::LinearOpWithSolveBase<double> >
lows = Thyra::linearOpWithSolve(*lowsFactory, A);
// Solve the linear system (note: the initial guess in 'x' is critical)
Thyra::SolveStatus<double>
status = Thyra::solve(*lows, Thyra::NOTRANS, *b, &*x);
*out << "\nSolve status:\n" << status;
// Write the linear solver parameters after they were read
linearSolverBuilder.writeParamsFile(*lowsFactory);
//
// E) Post process the solution and check the error
//
// Note that the below code is based only on the Epetra objects themselves
// and does not in any way depend or interact with any Thyra-based
// objects. The point is that most users of Thyra can largely gloss over
// the fact that Thyra is really being used for anything.
//
// Wipe out the Thyra wrapper for x to guarantee that the solution will be
// written back to epetra_x! At the time of this writting this is not
// really needed but the behavior may change at some point so this is a
// good idea.
x = Teuchos::null;
*out
<< "\nSolution ||epetra_x||2 = " << epetraNorm2(*epetra_x) << "\n";
*out << "\nTesting the solution error ||b-A*x||/||b|| computed through the Epetra objects ...\n";
// r = b - A*x
Epetra_Vector epetra_r(*epetra_b);
{
Epetra_Vector epetra_A_x(epetra_A->OperatorRangeMap());
epetra_A->Apply(*epetra_x,epetra_A_x);
epetra_r.Update(-1.0,epetra_A_x,1.0);
}
const double
nrm_r = epetraNorm2(epetra_r),
nrm_b = epetraNorm2(*epetra_b),
rel_err = ( nrm_r / nrm_b );
const bool
passed = (rel_err <= tol);
*out
<< "||b-A*x||/||b|| = " << nrm_r << "/" << nrm_b << " = " << rel_err
<< " < tol = " << tol << " ? " << ( passed ? "passed" : "failed" ) << "\n";
if(!passed) success = false;
epetra_A.release();
epetra_x.release();
epetra_b.release();
}
TEUCHOS_STANDARD_CATCH_STATEMENTS(verbose,std::cerr,success)
if(success) *out << "\nCongratulations! All of the tests checked out!\n";
else *out << "\nOh no! At least one of the tests failed!\n";
return ( success ? 0 : 1 );
}
-------------- next part --------------
A non-text attachment was scrubbed...
Name: AztecSettings.xml
Type: text/xml
Size: 1287 bytes
Desc: not available
Url : https://software.sandia.gov/pipermail/trilinos-users/attachments/20100115/15e0bed2/attachment-0001.xml
More information about the Trilinos-Users
mailing list