41 MPI_Comm_rank(mpicomm_, &myrank);
42 MPI_Comm_size(mpicomm_, &commsize);
47 std::vector<IntersectionData> dummy(1);
48 intersections_.swap(dummy);
51 std::vector<Dune::FieldVector<ctype, dimworld> > patch0coords;
52 std::vector<unsigned int> patch0entities;
53 std::vector<Dune::GeometryType> patch0types;
54 std::vector<Dune::FieldVector<ctype,dimworld> > patch1coords;
55 std::vector<unsigned int> patch1entities;
56 std::vector<Dune::GeometryType> patch1types;
64 extractGrid(patch<0>(), patch0coords, patch0entities, patch0types);
65 extractGrid(patch<1>(), patch1coords, patch1entities, patch1types);
67 std::cout <<
">>>> rank " << myrank <<
" coords: "
68 << patch0coords.size() <<
" and " << patch1coords.size() << std::endl;
69 std::cout <<
">>>> rank " << myrank <<
" entities: "
70 << patch0entities.size() <<
" and " << patch1entities.size() << std::endl;
71 std::cout <<
">>>> rank " << myrank <<
" types: "
72 << patch0types.size() <<
" and " << patch1types.size() << std::endl;
75 const char prefix[] =
"GridGlue::Builder::build() : ";
77 sprintf(patch0surf,
"/tmp/vtk-patch0-test-%i", myrank);
79 sprintf(patch1surf,
"/tmp/vtk-patch1-test-%i", myrank);
99 patch0_is_.beginResize();
100 patch1_is_.beginResize();
105 const int mergingrank,
106 const std::vector<Dune::FieldVector<ctype,dimworld> >& remotePatch0coords,
107 const std::vector<unsigned int>& remotePatch0entities,
108 const std::vector<Dune::GeometryType>& remotePatch0types,
109 const std::vector<Dune::FieldVector<ctype,dimworld> >& remotePatch1coords,
110 const std::vector<unsigned int>& remotePatch1entities,
111 const std::vector<Dune::GeometryType>& remotePatch1types
114 if (remotePatch1entities.size() > 0 && patch0entities.size() > 0)
115 mergePatches(patch0coords, patch0entities, patch0types, myrank,
116 remotePatch1coords, remotePatch1entities, remotePatch1types, mergingrank);
117 if (mergingrank != myrank &&
118 remotePatch0entities.size() > 0 && patch1entities.size() > 0)
119 mergePatches(remotePatch0coords, remotePatch0entities, remotePatch0types, mergingrank,
120 patch1coords, patch1entities, patch1types, myrank);
123 patch0coords, patch0entities, patch0types,
124 patch1coords, patch1entities, patch1types
130 patch0_is_.endResize();
131 patch1_is_.endResize();
134 remoteIndices_.setIncludeSelf(
true);
136 remoteIndices_.setIndexSets(patch0_is_, patch1_is_, mpicomm_) ;
137 remoteIndices_.rebuild<
true>();
140#ifdef DEBUG_GRIDGLUE_PARALLELMERGE
141 for (
auto it = remoteIndices_.begin(); it != remoteIndices_.end(); it++)
143 std::cout << myrank <<
"\tri-list\t" << it->first << std::endl;
144 for (
auto xit = it->second.first->begin(); xit != it->second.first->end(); ++xit)
145 std::cout << myrank <<
"\tri-list 1 \t" << it->first <<
"\t" << *xit << std::endl;
146 for (
auto xit = it->second.second->begin(); xit != it->second.second->end(); ++xit)
147 std::cout << myrank <<
"\tri-list 2 \t" << it->first <<
"\t" << *xit << std::endl;
153 if (patch1entities.size() > 0 && patch0entities.size() > 0)
155 mergePatches(patch0coords, patch0entities, patch0types, myrank,
156 patch1coords, patch1entities, patch1types, myrank);
157#ifdef CALL_MERGER_TWICE
158 mergePatches(patch0coords, patch0entities, patch0types, myrank,
159 patch1coords, patch1entities, patch1types, myrank);
180 const std::vector<Dune::FieldVector<ctype,dimworld> >& patch0coords,
181 const std::vector<unsigned int>& patch0entities,
182 const std::vector<Dune::GeometryType>& patch0types,
183 const int patch0rank,
184 const std::vector<Dune::FieldVector<ctype,dimworld> >& patch1coords,
185 const std::vector<unsigned int>& patch1entities,
186 const std::vector<Dune::GeometryType>& patch1types,
187 const int patch1rank)
195 MPI_Comm_rank(mpicomm_, &myrank);
196 MPI_Comm_size(mpicomm_, &commsize);
200 const bool patch0local = (myrank == patch0rank);
201 const bool patch1local = (myrank == patch1rank);
204 const unsigned int offset = intersections_.size()-1;
207 <<
" GridGlue::mergePatches : rank " << patch0rank <<
" / " << patch1rank << std::endl;
210 merger_->build(patch0coords, patch0entities, patch0types,
211 patch1coords, patch1entities, patch1types);
214 intersections_.resize(merger_->nSimplices() + offset + 1);
215 for (
unsigned int i = 0; i < merger_->nSimplices(); ++i)
216 intersections_[offset + i] =
IntersectionData(*
this, i, offset, patch0local, patch1local);
218 index__sz = intersections_.size() - 1;
221 <<
" GridGlue::mergePatches : "
222 <<
"The number of remote intersections is " << intersections_.size()-1 << std::endl;
226 printVector(patch0entities,
"patch0entities",myrank);
229 printVector(patch1entities,
"patch1entities",myrank);
237 assert(Dune::RESIZE == patch0_is_.state());
238 assert(Dune::RESIZE == patch1_is_.state());
240 for (
unsigned int i = 0; i < merger_->nSimplices(); i++)
244 GlobalId gid(patch0rank, patch1rank, i);
245 if (it.template local<0>())
247 Dune::PartitionType ptype = patch<0>().element(it.template index<0>()).partitionType();
248 patch0_is_.add (gid, LocalIndex(offset+i, ptype) );
250 if (it.template local<1>())
252 Dune::PartitionType ptype = patch<1>().element(it.template index<1>()).partitionType();
253 patch1_is_.add (gid, LocalIndex(offset+i, ptype) );
266 std::vector<Dune::FieldVector<ctype, dimworld> > & coords,
267 std::vector<unsigned int> & entities,
268 std::vector<Dune::GeometryType>& geometryTypes)
const
270 std::vector<typename Extractor::Coords> tempcoords;
271 std::vector<typename Extractor::VertexVector> tempentities;
275 coords.reserve(tempcoords.size());
277 for (
unsigned int i = 0; i < tempcoords.size(); ++i)
280 coords.push_back(Dune::FieldVector<ctype, dimworld>());
281 for (
size_t j = 0; j <dimworld; ++j)
282 coords.back()[j] = tempcoords[i][j];
288 for (
unsigned int i = 0; i < tempentities.size(); ++i) {
289 for (
unsigned int j = 0; j < tempentities[i].size(); ++j)
290 entities.push_back(tempentities[i][j]);
302 Dune::InterfaceType iftype, Dune::CommunicationDirection dir)
const
305 typedef typename DataHandle::DataType DataType;
309 if (mpicomm_ != MPI_COMM_SELF)
315 Dune::dinfo <<
"GridGlue: parallel communication" << std::endl;
316 typedef Dune::EnumItem <Dune::PartitionType, Dune::InteriorEntity> InteriorFlags;
317 typedef Dune::EnumItem <Dune::PartitionType, Dune::OverlapEntity> OverlapFlags;
318 typedef Dune::EnumRange <Dune::PartitionType, Dune::InteriorEntity, Dune::GhostEntity> AllFlags;
319 Dune::Interface interface;
320 assert(remoteIndices_.isSynced());
323 case Dune::InteriorBorder_InteriorBorder_Interface :
324 interface.build (remoteIndices_, InteriorFlags(), InteriorFlags() );
326 case Dune::InteriorBorder_All_Interface :
327 if (dir == Dune::ForwardCommunication)
328 interface.build (remoteIndices_, InteriorFlags(), AllFlags() );
330 interface.build (remoteIndices_, AllFlags(), InteriorFlags() );
332 case Dune::Overlap_OverlapFront_Interface :
333 interface.build (remoteIndices_, OverlapFlags(), OverlapFlags() );
335 case Dune::Overlap_All_Interface :
336 if (dir == Dune::ForwardCommunication)
337 interface.build (remoteIndices_, OverlapFlags(), AllFlags() );
339 interface.build (remoteIndices_, AllFlags(), OverlapFlags() );
341 case Dune::All_All_Interface :
342 interface.build (remoteIndices_, AllFlags(), AllFlags() );
345 DUNE_THROW(Dune::NotImplemented,
"GridGlue::communicate for interface " << iftype <<
" not implemented");
353 commInfo.
data = &data;
356 Dune::BufferedCommunicator bComm ;
357 bComm.template build< CommInfo >(commInfo, commInfo, interface);
361 if (dir == Dune::ForwardCommunication)
372 Dune::dinfo <<
"GridGlue: sequential fallback communication" << std::endl;
375 int ssz = size() * 10;
376 int rsz = size() * 10;
379 auto sendbuffer = std::make_unique<DataType[]>(ssz);
380 auto receivebuffer = std::make_unique<DataType[]>(rsz);
389 if (dir == Dune::ForwardCommunication)
396 data.
gather(gatherbuffer, in.inside(), in);
406 data.
gather(gatherbuffer, in.outside(), in.flip());
412 for (
int i=0; i<ssz; i++)
413 receivebuffer[i] = sendbuffer[i];
422 if (dir == Dune::ForwardCommunication)
428 data.
scatter(scatterbuffer, in.outside(), in.flip(),
437 data.
scatter(scatterbuffer, in.inside(), in,
void mergePatches(const std::vector< Dune::FieldVector< ctype, dimworld > > &patch0coords, const std::vector< unsigned int > &patch0entities, const std::vector< Dune::GeometryType > &patch0types, const int patch0rank, const std::vector< Dune::FieldVector< ctype, dimworld > > &patch1coords, const std::vector< unsigned int > &patch1entities, const std::vector< Dune::GeometryType > &patch1types, const int patch1rank)
after building the merged grid the intersection can be updated through this method (for internal use)
Definition gridglue.cc:179