[Trilinos-Users] bug in Tpetra::CrsGraph()?

[Nico Schlömer]

Hi all,

I just noticed an oddity about the dynamic profile of a Tpetra::CrsGraph.
I'd have expected that, when using a dynamic profile (which is the
default), one can just insert as many indices as one wants in each row,
independent of the maxNumEntriesPerRow argument of the graph constructor.
As I can't tell beforehand what the maximum number of entries per row will
be, I generically set this value to 0. This worked out fine with
Epetra_CrsGraph.

In Tpetra::CrsGraph, inserting one set of indices is never a problem;
inserting a second set in the same row, though, blows away some of the
indices previously inserted. If maxNumEntriesPerRow is set to the correct
number of entries per row, things work out fine.

A minimal example code is attached to this mail. The output I get is
===================== *snip* =====================
The graph after the *first* insertion:
G[ 0 , 0 ]
G[ 0 , 1 ]
G[ 0 , 2 ]


The graph after the *second* insertion:
G[ 0 , -1572835640 ]
G[ 0 , 32515 ]
G[ 0 , -1572835640 ]
G[ 0 , 0 ]
G[ 0 , 3 ]
G[ 0 , 4 ]


The graph after fillComplete:
G[ 0 , -1572835640 ]
G[ 0 , 0 ]
G[ 0 , 3 ]
G[ 0 , 4 ]
G[ 0 , 32515 ]
===================== *snap* =====================

I worked around this now by setting maxNumEntriesPerRow=N, but that's not
very pretty.

Cheers,
Nico
-------------- next part --------------
#include <Teuchos_DefaultComm.hpp>
#include <Tpetra_CrsGraph.hpp>

typedef int ORD;

void showGraph( const Teuchos::RCP<const Tpetra::CrsGraph<ORD> > graph );

int
main (int argc, char** argv)
{
#ifdef HAVE_MPI
  MPI_Init( &argc, &argv );
#endif

  // Create a communicator for Tpetra objects
  const Teuchos::RCP<const Teuchos::Comm<ORD> > comm =
        Teuchos::DefaultComm<int>::getComm();

  // allocate map
  int NumGlobalElements = 10;
  Teuchos::RCP<Tpetra::Map<ORD> > map =
      Teuchos::rcp( new Tpetra::Map<ORD>( NumGlobalElements, 0, comm ) );

  // allocate graph
  int maxNumEntriesPerRow = 0; // all works well when this is set >=5.
  Teuchos::RCP<Tpetra::CrsGraph<ORD> > graph =
      Teuchos::rcp( new Tpetra::CrsGraph<ORD>( map, map, maxNumEntriesPerRow ) );

  // insert the first set of indices
  Teuchos::ArrayRCP<ORD> indices(3);
  indices[0] = 0;
  indices[1] = 1;
  indices[2] = 2;
  graph->insertLocalIndices( 0, indices() );
  
  // show graph
  std::cout << "The graph after the *first* insertion:" << std::endl;
  showGraph( graph );

  // insert the second set of indices
  indices[0] = 0;
  indices[1] = 3;
  indices[2] = 4;
  graph->insertLocalIndices( 0, indices() );

  // show graph
  std::cout << "\n\nThe graph after the *second* insertion:" << std::endl;
  showGraph( graph );
  
  graph->globalAssemble();
  graph->fillComplete();

  // show graph
  std::cout << "\n\nThe graph after fillComplete:" << std::endl;
  showGraph( graph );  
                           
#ifdef HAVE_MPI
  MPI_Finalize();
#endif

  return( EXIT_SUCCESS );
}

void
showGraph( const Teuchos::RCP<const Tpetra::CrsGraph<ORD> > graph )
{
  for ( int i=0; i<graph->getRowMap()->getNodeNumElements(); i++ )
  {
      Teuchos::ArrayRCP<const ORD> indices = graph->getLocalRowView( i ); 
      for ( int l=0; l<indices.size(); l++ )
          std::cout << "G[ " << i << " , " << indices[l] << " ]" << std::endl;
  }
  return;
}
-------------- next part --------------
CMAKE_MINIMUM_REQUIRED(VERSION 2.8)

FIND_PACKAGE( Trilinos REQUIRED )

SET( CMAKE_C_COMPILER ${Trilinos_C_COMPILER} )
SET( CMAKE_CXX_COMPILER ${Trilinos_CXX_COMPILER} )
SET( CMAKE_Fortran_COMPILER ${Trilinos_Fortran_COMPILER} )

PROJECT( extract CXX Fortran )

INCLUDE_DIRECTORIES ( ${Trilinos_INCLUDE_DIRS}
                      ${Trilinos_TPL_INCLUDE_DIRS} )
LINK_DIRECTORIES( ${Trilinos_LIBRARY_DIRS}
                  ${Trilinos_TPL_LIBRARY_DIRS} )

SET( MY_EXECUTABLE
     "crsgraph.exe" )
SET ( extract_SRCS
      main.cpp )
ADD_EXECUTABLE( ${MY_EXECUTABLE}
                ${extract_SRCS} )

SET_TARGET_PROPERTIES( ${MY_EXECUTABLE}
                       PROPERTIES LINKER_LANGUAGE Fortran )

TARGET_LINK_LIBRARIES( ${MY_EXECUTABLE}
                       ${Trilinos_LIBRARIES}
                       ${Trilinos_TPL_LIBRARIES} )