Epetra Norm2 failure in parallel

[wd15]

The following code fails with trilinos 10.8.3 on bunter, sandbox and Jon's mac.

from fipy import *
L = 1.
nx = 6
m = Grid1D(nx=nx, dx=L / nx)
v = CellVariable(mesh=m, value=1.)
eqn = DiffusionTerm() - v

v.constrain(0.,  m.facesLeft)
v.constrain(1.,  m.facesRight)

res = 1.
sweep = 0
while res > 1e-8 and sweep < 100:
res = eqn.sweep(v)
sweep += 1
x = m.cellCenters[0]
answer = (numerix.exp(x) - numerix.exp(-x)) / (numerix.exp(L) - numerix.exp(-L))
print numerix.allclose(v, answer, rtol=2e-5)

in the following way

bunter)bunter[wd15]: mpirun -np 2 python tmp2.py
Traceback (most recent call last):
  File "tmp2.py", line 14, in <module>
Traceback (most recent call last):
  File "tmp2.py", line 14, in <module>
res = eqn.sweep(v)
  File "/users/wd15/Documents/python/fipy/bunter/fipy/terms/term.py", line 237, in sweep
res = eqn.sweep(v)
  File "/users/wd15/Documents/python/fipy/bunter/fipy/terms/term.py", line 237, in sweep
residual = solver._calcResidual(residualFn=residualFn)
  File "/users/wd15/Documents/python/fipy/bunter/fipy/solvers/trilinos/trilinosSolver.py", line 170, in _calcResidual
residual = solver._calcResidual(residualFn=residualFn)
  File "/users/wd15/Documents/python/fipy/bunter/fipy/solvers/trilinos/trilinosSolver.py", line 170, in _calcResidual
return comm.Norm2(self._calcResidualVector())
  File "/users/wd15/Documents/python/fipy/bunter/fipy/tools/commWrapper.py", line 97, in Norm2
return comm.Norm2(self._calcResidualVector())
  File "/users/wd15/Documents/python/fipy/bunter/fipy/tools/commWrapper.py", line 97, in Norm2
return vec.Norm2()
  File "/users/wd15/.virtualenvs/bunter/lib/python2.6/site-packages/PyTrilinos/Epetra.py", line 9594, in Norm2
return vec.Norm2()
  File "/users/wd15/.virtualenvs/bunter/lib/python2.6/site-packages/PyTrilinos/Epetra.py", line 9594, in Norm2
return _Epetra.NumPyVector_Norm2(self, *args)
RuntimeError: Norm2 returned error code -1
return _Epetra.NumPyVector_Norm2(self, *args)
RuntimeError: Norm2 returned error code -1

Imported from trac ticket #370, created by wd15 on 10-28-2011 at 11:32, last modified: 09-30-2013 at 21:44

[wd15]

As I suspected, the maps have changed. If the colMap is printed at line 156 in [bunter/fipy/solvers/trilinos/trilinosSolver.py#L156@4608 trilinosSolver.py]
it gives the following:

 globalMatrix.colMap 
Number of Global Elements  = 10
Number of Global Points    = 10
Maximum of all GIDs        = 5
Minimum of all GIDs        = 0
Index Base                 = 0
Constant Element Size      = 1

Number of Local Elements   = 5
Number of Local Points     = 5
Maximum of my GIDs         = 4
Minimum of my GIDs         = 0

MyPID           Local Index        Global Index  
0                 0                 0    
0                 1                 1    
0                 2                 2    
0                 3                 3    
0                 4                 4    

globalMatrix.colMap 
Number of Local Elements   = 5
Number of Local Points     = 5
Maximum of my GIDs         = 5
Minimum of my GIDs         = 1

MyPID           Local Index        Global Index  
1                 0                 1    
1                 1                 2    
1                 2                 3    
1                 3                 4    
1                 4                 5    

which is really bad news. There are some other test failures that might be related to this. Next job is to install an earlier version of trilinos and put together a simple example which shows the difference between the versions.

Trac comment by wd15 on 10-28-2011 at 13:01

[wd15]

I've also been trying to look at an Epetra C script.

#include "Epetra_SerialComm.h"
#include "Epetra_Map.h"
#include "Epetra_Vector.h"
#include "Epetra_Version.h"
#include "mpi.h"
#include "Epetra_MpiComm.h"
#include <vector>


int main(int argc, char *argv[])
{

  cout << Epetra_Version() << endl << endl;
  MPI_Init (&argc,&argv);
  Epetra_MpiComm Comm (MPI_COMM_WORLD);

  int MyPID = Comm.MyPID();

  cout << 'MyPID' << MyPID << endl;
  //int NumElements = 
  int N = 3;
  int rows[N];
  if (MyPID == 0)
{
rows[0] = 0;
rows[1] = 1;
rows[2] = 2;
}
  else
{
rows[0] = 1;
rows[1] = 2;
rows[2] = 3;
}
  // Construct a Map with NumElements and index base of 0
  //Epetra_Map Map(4, 3, rows, 0, Comm);
  Epetra_Map Map(6, 3, rows, 0, Comm);
  Epetra_Vector b(Map);

  b.Random();
  cout << "b" << b << endl;
  double bnorm;

  b.Norm2(&bnorm);

  cout << "2 norm of b = " << bnorm << endl;

  return 0;
}

To compile this I needed to do the following.

  $ export CPLUS_INCLUDE_PATH=~/.virtualenvs/bunter/include:/usr/include/mpi
  $ gcc tmp.cc -o tmp -lepetra -lmpi

Trac comment by wd15 on 10-28-2011 at 13:13

[wd15]

This issue has been fixed with r4816. Unfortunately, I have been unable to understand the changes in trilinos that caused Norm2 to start failing for overlapping vectors. It fails for versions of trilinos going back to 10.2.2.

Trac comment by wd15 on 10-31-2011 at 16:52

[fipymigrate]

In 84d55df:

#CommitTicketReference repository="fipy" revision="84d55df5ed448fe0eed90e96d2547aa7075043bc"
Fixed issue #274. Still 1 test failure in parallel.

git-svn-id: svn+ssh://code.matforge.org/fipy/branches/bunter@4816 d80e17d7-ff13-0410-a124-85740d801063

Trac comment by Daniel Wheeler daniel.wheeler@nist.gov on 06-12-2013 at 14:14

[guyer]

Marking milestone

Trac comment by guyer on 09-30-2013 at 21:44