Dear libmesh developers

VTK output crashes when run in parallel (mpirun -n 2 ...). The following 
assertion fails:

Assertion `it!=_global_to_local_map.end()' failed.
[...]/libmesh/include/numerics/petsc_vector.h, line 1073, [...]

The problem is reproducible with the following minimal working example:

//START MINIMAL WORKING EXAMPLE CODE
#include "libmesh.h"
#include "mesh.h"
#include "mesh_generation.h"
#include "equation_systems.h"
#include "explicit_system.h"
#include "vtk_io.h"

using namespace libMesh;

int main(int argc, char** argv)
{
  LibMeshInit init(argc, argv);
  {
    Mesh mesh;
    MeshTools::Generation::build_line(mesh, 4, 0, 1, EDGE2);
    EquationSystems equation_systems(mesh);
    ExplicitSystem& system = equation_systems.add_system<ExplicitSystem>("sys");
    system.add_variable("u");
    equation_systems.init();
    VTKIO(mesh).write_equation_systems("out.pvtu", equation_systems); // causes 
assertion failure
  }
  return 0;
}
//END MINIMAL WORKING EXAMPLE CODE


The stack trace looks as follows:


0: libMesh::print_trace(std::ostream&)
1: libMesh::PetscVector<double>::map_global_to_local_index(unsigned int) const
2: libMesh::PetscVector<double>::operator()(unsigned int) const
3: libMesh::System::current_solution(unsigned int) const
4: libMesh::VTKIO::solution_to_vtk(libMesh::EquationSystems const&, 
vtkUnstructuredGrid*&)
5: libMesh::VTKIO::write_equation_systems(std::string const&, 
libMesh::EquationSystems const&)
6: main
7: __libc_start_main


The assertion does NOT fail when running in serial mode (without mpirun) or 
when using ExodusII_IO. Since I'm bound to VTK and parallel processing, I will 
need to fix this. The problem seems to be that current_solution.operator() is 
called for ALL dofs in line 333 of vtk_io.C, instead of only the local ones. 
Therefore, the global-to-local index map doesn't contain the sought indices, 
when processor 0 attempts to write data owned by processor 1.

Solution output appears to be more sophisticated in ExodusII_IO, where 
MeshOutput<MT>::_build_variable_names_and_solution_vector is called, which does 
a EquationSystems::allgather().

What's the expert's ansatz to fixing VTK output? Should I have a go at 
duplicating the ExodusII output procedure, or does anybody see a simple(r) way 
of fixing the problem? Perhaps a localize_to_one() call on the solution vector? 
Would that compromise parallel efficiency?
------------------------------------------------------------------------------
The ultimate all-in-one performance toolkit: Intel(R) Parallel Studio XE:
Pinpoint memory and threading errors before they happen.
Find and fix more than 250 security defects in the development cycle.
Locate bottlenecks in serial and parallel code that limit performance.
http://p.sf.net/sfu/intel-dev2devfeb
_______________________________________________
Libmesh-devel mailing list
Libmesh-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/libmesh-devel

Reply via email to