Hey guys,

the current MeshWorker::loop() uses the WorkStream::run() which copies in 
serial.
I tried to use the WorkStream::run() version which uses colored iterators 
and therefore assembles in parallel, in the sense of a tbb::parallel_for().

My code runs without using MPI as well as using an arbitrary amount of MPI 
procs with a single thread each. If I use 2 MPI procs with 2 threads each, 
it crashes due to a segmentation fault.
For the construction of the colored iterators I am using the 
make_graph_coloring() in such a way that two cells sharing a face have 
different colors.

The segmentation fault always appears while calling assemble() of the 
ResidualSimple in simple.h.

Does anybody have a clue why I get a segmentation fault after adding local 
data to the global dst vector inside assemble(), although I avoid write 
conflicts by colored iterators?

In the "error_message" I have attached the output my code writes for two 
MPI procs with two threads, respectively.

Best,
Julius

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to [email protected].
For more options, visit https://groups.google.com/d/optout.

Attachment: error_message
Description: Binary data

//dealii
#include <deal.II/base/graph_coloring.h>
#include <deal.II/base/multithread_info.h>
#include <deal.II/base/template_constraints.h>
#include <deal.II/base/work_stream.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/fe/fe_dgq.h>
#include <deal.II/fe/fe_values.h>
#include <deal.II/fe/mapping_q1.h>
#include <deal.II/grid/tria.h>
#include <deal.II/grid/grid_generator.h>
#include <deal.II/integrators/l2.h>
#include <deal.II/lac/vector.h>
#include <deal.II/meshworker/dof_info.h>
#include <deal.II/meshworker/integration_info.h>
#include <deal.II/meshworker/local_integrator.h>
#include <deal.II/meshworker/loop.h>
#include <deal.II/meshworker/simple.h>
#include <deal.II/integrators/laplace.h>

//mpi
#include <deal.II/distributed/tria.h>
#include <deal.II/lac/trilinos_vector.h>

//c++
#include <iostream>
#include <vector>
#include <cstdio>

template<int dim, int spacedim, typename ITERATOR, typename DOFINFO, typename INFOBOX, typename ASSEMBLER>
void colored_loop(const std::vector<std::vector<ITERATOR> > colored_iterators,
		  DOFINFO  &dof_info,
		  INFOBOX  &info,
		  const dealii::std_cxx11::function<void (DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker,
		  const dealii::std_cxx11::function<void (DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker,
		  const dealii::std_cxx11::function<void (DOFINFO &, DOFINFO &,
							  typename INFOBOX::CellInfo &,
							  typename INFOBOX::CellInfo &)> &face_worker,
		  ASSEMBLER  &assembler,
		  const dealii::MeshWorker::LoopControl &lctrl = dealii::MeshWorker::LoopControl())
{
  dealii::MeshWorker::DoFInfoBox<dim, DOFINFO> dinfo_box(dof_info);

  assembler.initialize_info(dinfo_box.cell, false);
  for (unsigned int i=0; i<dealii::GeometryInfo<dim>::faces_per_cell; ++i)
    {
      assembler.initialize_info(dinfo_box.interior[i], true);
      assembler.initialize_info(dinfo_box.exterior[i], true);
    }

  //  Loop over all cells                                                                                                              
  dealii::WorkStream::run(colored_iterators,
			  dealii::std_cxx11::bind(&dealii::MeshWorker::cell_action<INFOBOX, DOFINFO, dim, spacedim, ITERATOR>,
						  dealii::std_cxx11::_1, dealii::std_cxx11::_3, dealii::std_cxx11::_2,
						  cell_worker, boundary_worker, face_worker, lctrl),
			  dealii::std_cxx11::bind(&dealii::internal::assemble<dim,DOFINFO,ASSEMBLER>,
						  dealii::std_cxx11::_1, &assembler),
			  info, dinfo_box);
}

int main(int argc, char** argv)
{
  dealii::Utilities::MPI::MPI_InitFinalize mpi_initialization(argc,argv,2); 

  const unsigned int dim = 2;
  const unsigned int degree = 1;

  //template definitions
  typedef dealii::MeshWorker::DoFInfo<dim> DOFINFO ;
  typedef dealii::MeshWorker::IntegrationInfoBox<dim> INFOBOX ;
  typedef dealii::MeshWorker::Assembler::ResidualSimple<dealii::TrilinosWrappers::MPI::Vector> ASSEMBLER ;
  typedef typename dealii::DoFHandler<dim>::active_cell_iterator ITERATOR ;
  
  dealii::parallel::distributed::Triangulation<2> tr(MPI_COMM_WORLD);
  const dealii::FE_DGQ<dim> fe{degree};
  dealii::DoFHandler<dim> dof_handler{tr};
  const dealii::MappingQ1<dim> mapping;
  
  dealii::GridGenerator::hyper_cube (tr,0.,1.);
  tr.refine_global(4);
  dof_handler.distribute_dofs(fe);

  dealii::std_cxx11::function<std::vector<dealii::types::global_dof_index>
                              (const ITERATOR& )>
    conflicts
    = [] (const ITERATOR& cell)
    {
      std::vector< dealii::types::global_dof_index > conflict_indices;

      if( !(cell->is_artificial()) )
        {
          conflict_indices.resize(cell->get_fe().dofs_per_cell);
          cell->get_active_or_mg_dof_indices(conflict_indices);
                                                                                                    
	  std::vector< dealii::types::global_dof_index > cell_indices;
          for (unsigned int face_no=0; face_no < dealii::GeometryInfo<dim>::faces_per_cell; ++face_no)
            {
              if ( !(cell->at_boundary(face_no)) && !(cell->neighbor(face_no)->is_artificial()) )
                {
                  auto neighbor = cell->neighbor(face_no);
		  cell_indices.clear();
		  cell_indices.resize(cell->get_fe().dofs_per_cell);
		  neighbor->get_active_or_mg_dof_indices(cell_indices);
		  for( unsigned int i=0; i<cell_indices.size(); ++i)
		    conflict_indices.push_back(cell_indices[i]);
                }
            }
	}
      return conflict_indices;
    };

  std::vector<std::vector<ITERATOR> >
    colored_iterators = dealii::GraphColoring::make_graph_coloring ( dof_handler.begin_active(),
								     dof_handler.end(),
								     conflicts );

  std::cout << "n_colors= " << colored_iterators.size()  << std::endl;
  std::cout << "n_cells= " << tr.n_active_cells() << std::endl;
  std::cout << "n_dofs= " << dof_handler.n_dofs() << std::endl;

  const dealii::IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs();

  dealii::TrilinosWrappers::MPI::Vector dst(locally_owned_dofs,MPI_COMM_WORLD);
  dealii::MeshWorker::Assembler::ResidualSimple<dealii::TrilinosWrappers::MPI::Vector> assembler;
  dealii::AnyData dst_data;
  dst_data.add<dealii::TrilinosWrappers::MPI::Vector* >(&dst, "dst");
  assembler.initialize(dst_data);

  INFOBOX info_box;
  DOFINFO dof_info(dof_handler);
  
  colored_loop<dim,dim,ITERATOR,DOFINFO,INFOBOX,ASSEMBLER>
    (colored_iterators,
     dof_info, info_box,
     nullptr,nullptr,nullptr,
     assembler);
}

Reply via email to