[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