Actual source code: Sifter.hh
petsc-3.4.5 2014-06-29
1: #ifndef included_ALE_Sifter_hh
2: #define included_ALE_Sifter_hh
4: /*
5: #include <boost/multi_index_container.hpp>
6: #include <boost/multi_index/member.hpp>
7: #include <boost/multi_index/ordered_index.hpp>
8: #include <boost/multi_index/composite_key.hpp>
9: */
10: #include <iostream>
12: // ALE extensions
14: #ifndef included_ALE_hh
15: #include <sieve/ALE.hh>
16: #endif
18: extern PetscErrorCode PetscObjectDestroy_PetscObject(PetscObject*);
20: namespace ALE {
22: namespace SifterDef {
23: // Defines the traits of a sequence representing a subset of a multi_index container Index_.
24: // A sequence defines output (input in std terminology) iterators for traversing an Index_ object.
25: // Upon dereferencing values are extracted from each result record using a ValueExtractor_ object.
26: template <typename Index_, typename ValueExtractor_>
27: struct IndexSequenceTraits {
28: typedef Index_ index_type;
29: class iterator_base {
30: public:
31: // Standard iterator typedefs
32: typedef ValueExtractor_ extractor_type;
33: typedef std::input_iterator_tag iterator_category;
34: typedef typename extractor_type::result_type value_type;
35: typedef int difference_type;
36: typedef value_type* pointer;
37: typedef value_type& reference;
39: // Underlying iterator type
40: typedef typename index_type::iterator itor_type;
41: protected:
42: // Underlying iterator
43: itor_type _itor;
44: // Member extractor
45: extractor_type _ex;
46: public:
47: iterator_base(itor_type itor) {
48: this->_itor = itor_type(itor);
49: };
50: virtual ~iterator_base() {};
51: virtual bool operator==(const iterator_base& iter) const {return this->_itor == iter._itor;};
52: virtual bool operator!=(const iterator_base& iter) const {return this->_itor != iter._itor;};
53: // FIX: operator*() should return a const reference, but it won't compile that way, because _ex() returns const value_type
54: virtual const value_type operator*() const {return _ex(*(this->_itor));};
55: };// class iterator_base
56: class iterator : public iterator_base {
57: public:
58: // Standard iterator typedefs
59: typedef typename iterator_base::iterator_category iterator_category;
60: typedef typename iterator_base::value_type value_type;
61: typedef typename iterator_base::extractor_type extractor_type;
62: typedef typename iterator_base::difference_type difference_type;
63: typedef typename iterator_base::pointer pointer;
64: typedef typename iterator_base::reference reference;
65: // Underlying iterator type
66: typedef typename iterator_base::itor_type itor_type;
67: public:
68: iterator(const itor_type& itor) : iterator_base(itor) {};
69: virtual ~iterator() {};
70: //
71: virtual iterator operator++() {++this->_itor; return *this;};
72: virtual iterator operator++(int n) {iterator tmp(this->_itor); ++this->_itor; return tmp;};
73: };// class iterator
74: }; // struct IndexSequenceTraits
76: template <typename Index_, typename ValueExtractor_>
77: struct ReversibleIndexSequenceTraits {
78: typedef IndexSequenceTraits<Index_, ValueExtractor_> base_traits;
79: typedef typename base_traits::iterator_base iterator_base;
80: typedef typename base_traits::iterator iterator;
81: typedef typename base_traits::index_type index_type;
83: // reverse_iterator is the reverse of iterator
84: class reverse_iterator : public iterator_base {
85: public:
86: // Standard iterator typedefs
87: typedef typename iterator_base::iterator_category iterator_category;
88: typedef typename iterator_base::value_type value_type;
89: typedef typename iterator_base::extractor_type extractor_type;
90: typedef typename iterator_base::difference_type difference_type;
91: typedef typename iterator_base::pointer pointer;
92: typedef typename iterator_base::reference reference;
93: // Underlying iterator type
94: typedef typename iterator_base::itor_type itor_type;
95: public:
96: reverse_iterator(const itor_type& itor) : iterator_base(itor) {};
97: virtual ~reverse_iterator() {};
98: //
99: virtual reverse_iterator operator++() {--this->_itor; return *this;};
100: virtual reverse_iterator operator++(int n) {reverse_iterator tmp(this->_itor); --this->_itor; return tmp;};
101: };
102: }; // class ReversibleIndexSequenceTraits
105: //
106: // Rec & RecContainer definitions.
107: // Rec is intended to denote a graph point record.
108: //
109: template <typename Point_>
110: struct Rec {
111: typedef Point_ point_type;
112: template<typename OtherPoint_>
113: struct rebind {
114: typedef Rec<OtherPoint_> type;
115: };
116: point_type point;
117: int degree;
118: // Basic interface
119: Rec() : degree(0){};
120: Rec(const Rec& r) : point(r.point), degree(r.degree) {}
121: //Rec(const point_type& p) : point(p), degree(0) {};
122: Rec(const point_type& p, const int d) : point(p), degree(d) {};
123: // Printing
124: friend std::ostream& operator<<(std::ostream& os, const Rec& p) {
125: os << "<" << p.point << ", "<< p.degree << ">";
126: return os;
127: };
129: struct degreeAdjuster {
130: degreeAdjuster(int newDegree) : _newDegree(newDegree) {};
131: void operator()(Rec& r) { r.degree = this->_newDegree; }
132: private:
133: int _newDegree;
134: };// degreeAdjuster()
136: };// class Rec
138: template <typename Point_, typename Rec_>
139: struct RecContainerTraits {
140: typedef Rec_ rec_type;
141: // Index tags
142: struct pointTag{};
143: // Rec set definition
144: typedef ::boost::multi_index::multi_index_container<
145: rec_type,
146: ::boost::multi_index::indexed_by<
147: ::boost::multi_index::ordered_unique<
148: ::boost::multi_index::tag<pointTag>, BOOST_MULTI_INDEX_MEMBER(rec_type, typename rec_type::point_type, point)
149: >
150: >,
151: ALE_ALLOCATOR<rec_type>
152: > set_type;
153: //
154: // Return types
155: //
157: class PointSequence {
158: public:
159: typedef IndexSequenceTraits<typename ::boost::multi_index::index<set_type, pointTag>::type,
160: BOOST_MULTI_INDEX_MEMBER(rec_type, typename rec_type::point_type,point)>
161: traits;
162: protected:
163: const typename traits::index_type& _index;
164: public:
166: // Need to extend the inherited iterator to be able to extract the degree
167: class iterator : public traits::iterator {
168: public:
169: iterator(const typename traits::iterator::itor_type& itor) : traits::iterator(itor) {};
170: virtual const int& degree() const {return this->_itor->degree;};
171: };
173: PointSequence(const PointSequence& seq) : _index(seq._index) {};
174: PointSequence(typename traits::index_type& index) : _index(index) {};
175: virtual ~PointSequence(){};
177: virtual bool empty(){return this->_index.empty();};
179: virtual typename traits::index_type::size_type size() {return this->_index.size();};
181: virtual iterator begin() {
182: // Retrieve the beginning iterator of the index
183: return iterator(this->_index.begin());
184: };
185: virtual iterator end() {
186: // Retrieve the ending iterator of the index
187: // Since the elements in this index are ordered by degree, this amounts to the end() of the index.
188: return iterator(this->_index.end());
189: };
190: virtual bool contains(const typename rec_type::point_type& p) {
191: // Check whether a given point is in the index
192: return (this->_index.find(p) != this->_index.end());
193: }
194: }; // class PointSequence
195: };// struct RecContainerTraits
198: template <typename Point_, typename Rec_>
199: struct RecContainer {
200: typedef RecContainerTraits<Point_, Rec_> traits;
201: typedef typename traits::set_type set_type;
202: template <typename OtherPoint_, typename OtherRec_>
203: struct rebind {
204: typedef RecContainer<OtherPoint_, OtherRec_> type;
205: };
206: set_type set;
207: //
208: void removePoint(const typename traits::rec_type::point_type& p) {
209: /*typename ::boost::multi_index::index<set_type, typename traits::pointTag>::type& index =
210: ::boost::multi_index::get<typename traits::pointTag>(this->set);
211: typename ::boost::multi_index::index<set_type, typename traits::pointTag>::type::iterator i = index.find(p);
212: if (i != index.end()) { // Point exists
213: i = index.erase(i);
214: }*/
215: this->set.erase(p);
216: };
217: //
218: void adjustDegree(const typename traits::rec_type::point_type& p, int delta) {
219: typename ::boost::multi_index::index<set_type, typename traits::pointTag>::type& index =
220: ::boost::multi_index::get<typename traits::pointTag>(this->set);
221: typename ::boost::multi_index::index<set_type, typename traits::pointTag>::type::iterator i = index.find(p);
222: if (i == index.end()) { // No such point exists
223: if(delta < 0) { // Cannot decrease degree of a non-existent point
224: ostringstream err;
225: err << "ERROR: adjustDegree: Non-existent point " << p;
226: std::cout << err << std::endl;
227: throw(Exception(err.str().c_str()));
228: }
229: else { // We CAN INCREASE the degree of a non-existent point: simply insert a new element with degree == delta
230: std::pair<typename ::boost::multi_index::index<set_type, typename traits::pointTag>::type::iterator, bool> ii;
231: typename traits::rec_type r(p,delta);
232: ii = index.insert(r);
233: if(ii.second == false) {
234: ostringstream err;
235: err << "ERROR: adjustDegree: Failed to insert a rec " << r;
236: std::cout << err << std::endl;
237: throw(Exception(err.str().c_str()));
238: }
239: }
240: }
241: else { // Point exists, so we try to modify its degree
242: // If the adjustment is zero, there is nothing to do, otherwise ...
243: if(delta != 0) {
244: int newDegree = i->degree + delta;
245: if(newDegree < 0) {
246: ostringstream ss;
247: ss << "adjustDegree: Adjustment of " << *i << " by " << delta << " would result in negative degree: " << newDegree;
248: throw Exception(ss.str().c_str());
249: }
250: index.modify(i, typename traits::rec_type::degreeAdjuster(newDegree));
251: }
252: }
253: }; // adjustDegree()
254: }; // struct RecContainer
256: //
257: // Arrow & ArrowContainer definitions
258: //
259: template<typename Source_, typename Target_, typename Color_>
260: struct Arrow { //: public ALE::def::Arrow<Source_, Target_, Color_> {
261: typedef Arrow arrow_type;
262: typedef Source_ source_type;
263: typedef Target_ target_type;
264: typedef Color_ color_type;
265: source_type source;
266: target_type target;
267: color_type color;
268: Arrow(const source_type& s, const target_type& t, const color_type& c) : source(s), target(t), color(c) {};
269: // Flipping
270: template <typename OtherSource_, typename OtherTarget_, typename OtherColor_>
271: struct rebind {
272: typedef Arrow<OtherSource_, OtherTarget_, OtherColor_> type;
273: };
274: struct flip {
275: typedef Arrow<target_type, source_type, color_type> type;
276: type arrow(const arrow_type& a) { return type(a.target, a.source, a.color);};
277: };
279: // Printing
280: friend std::ostream& operator<<(std::ostream& os, const Arrow& a) {
281: os << a.source << " --(" << a.color << ")--> " << a.target;
282: return os;
283: }
285: // Arrow modifiers
286: struct sourceChanger {
287: sourceChanger(const source_type& newSource) : _newSource(newSource) {};
288: void operator()(arrow_type& a) {a.source = this->_newSource;}
289: private:
290: source_type _newSource;
291: };
293: struct targetChanger {
294: targetChanger(const target_type& newTarget) : _newTarget(newTarget) {};
295: void operator()(arrow_type& a) { a.target = this->_newTarget;}
296: private:
297: const target_type _newTarget;
298: };
299: };// struct Arrow
302: template<typename Source_, typename Target_, typename Color_, typename SupportCompare_>
303: struct ArrowContainerTraits {
304: public:
305: //
306: // Encapsulated types
307: //
308: typedef Arrow<Source_,Target_,Color_> arrow_type;
309: typedef typename arrow_type::source_type source_type;
310: typedef typename arrow_type::target_type target_type;
311: typedef typename arrow_type::color_type color_type;
312: typedef SupportCompare_ support_compare_type;
313: // Index tags
314: struct sourceColorTag{};
315: struct targetColorTag{};
316: struct sourceTargetTag{};
318: // Sequence traits and sequence types
319: template <typename Index_, typename Key_, typename SubKey_, typename ValueExtractor_>
320: class ArrowSequence {
321: // ArrowSequence implements ReversibleIndexSequencTraits with Index_ and ValueExtractor_ types.
322: // A Key_ object and an optional SubKey_ object are used to extract the index subset.
323: public:
324: typedef ReversibleIndexSequenceTraits<Index_, ValueExtractor_> traits;
325: //typedef source_type source_type;
326: //typedef target_type target_type;
327: //typedef arrow_type arrow_type;
328: //
329: typedef Key_ key_type;
330: typedef SubKey_ subkey_type;
331: protected:
332: typename traits::index_type& _index;
333: key_type key;
334: subkey_type subkey;
335: bool useSubkey;
336: public:
337: // Need to extend the inherited iterators to be able to extract arrow color
338: class iterator : public traits::iterator {
339: public:
340: iterator(const typename traits::iterator::itor_type& itor) : traits::iterator(itor) {};
341: virtual const source_type& source() const {return this->_itor->source;};
342: virtual const color_type& color() const {return this->_itor->color;};
343: virtual const target_type& target() const {return this->_itor->target;};
344: virtual const arrow_type& arrow() const {return *(this->_itor);};
345: };
346: class reverse_iterator : public traits::reverse_iterator {
347: public:
348: reverse_iterator(const typename traits::reverse_iterator::itor_type& itor) : traits::reverse_iterator(itor) {};
349: virtual const source_type& source() const {return this->_itor->source;};
350: virtual const color_type& color() const {return this->_itor->color;};
351: virtual const target_type& target() const {return this->_itor->target;};
352: virtual const arrow_type& arrow() const {return *(this->_itor);};
353: };
354: public:
355: //
356: // Basic ArrowSequence interface
357: //
358: ArrowSequence(const ArrowSequence& seq) : _index(seq._index), key(seq.key), subkey(seq.subkey), useSubkey(seq.useSubkey) {};
359: ArrowSequence(typename traits::index_type& index, const key_type& k) :
360: _index(index), key(k), subkey(subkey_type()), useSubkey(0) {};
361: ArrowSequence(typename traits::index_type& index, const key_type& k, const subkey_type& kk) :
362: _index(index), key(k), subkey(kk), useSubkey(1){};
363: virtual ~ArrowSequence() {};
365: void setKey(const key_type& key) {this->key = key;};
366: void setSubkey(const subkey_type& subkey) {this->subkey = subkey;};
367: void setUseSubkey(const bool& useSubkey) {this->useSubkey = useSubkey;};
369: virtual bool empty() {return this->_index.empty();};
371: virtual typename traits::index_type::size_type size() {
372: if (this->useSubkey) {
373: return this->_index.count(::boost::make_tuple(this->key,this->subkey));
374: } else {
375: return this->_index.count(::boost::make_tuple(this->key));
376: }
377: };
379: virtual iterator begin() {
380: if (this->useSubkey) {
381: return iterator(this->_index.lower_bound(::boost::make_tuple(this->key,this->subkey)));
382: } else {
383: return iterator(this->_index.lower_bound(::boost::make_tuple(this->key)));
384: }
385: };
387: virtual iterator end() {
388: if (this->useSubkey) {
389: return iterator(this->_index.upper_bound(::boost::make_tuple(this->key,this->subkey)));
390: } else {
391: return iterator(this->_index.upper_bound(::boost::make_tuple(this->key)));
392: }
393: };
395: virtual reverse_iterator rbegin() {
396: if (this->useSubkey) {
397: return reverse_iterator(--this->_index.upper_bound(::boost::make_tuple(this->key,this->subkey)));
398: } else {
399: return reverse_iterator(--this->_index.upper_bound(::boost::make_tuple(this->key)));
400: }
401: };
403: virtual reverse_iterator rend() {
404: if (this->useSubkey) {
405: return reverse_iterator(--this->_index.lower_bound(::boost::make_tuple(this->key,this->subkey)));
406: } else {
407: return reverse_iterator(--this->_index.lower_bound(::boost::make_tuple(this->key)));
408: }
409: };
411: template<typename ostream_type>
412: void view(ostream_type& os, const bool& useColor = false, const char* label = NULL){
413: if(label != NULL) {
414: os << "Viewing " << label << " sequence:" << std::endl;
415: }
416: os << "[";
417: for(iterator i = this->begin(); i != this->end(); i++) {
418: os << " (" << *i;
419: if(useColor) {
420: os << "," << i.color();
421: }
422: os << ")";
423: }
424: os << " ]" << std::endl;
425: }
426: };// class ArrowSequence
427: };// class ArrowContainerTraits
430: // The specialized ArrowContainer types distinguish the cases of unique and multiple colors of arrows on
431: // for each (source,target) pair (i.e., a single arrow, or multiple arrows between each pair of points).
432: typedef enum {multiColor, uniColor} ColorMultiplicity;
434: template<typename Source_, typename Target_, typename Color_, ColorMultiplicity colorMultiplicity, typename SupportCompare_>
435: struct ArrowContainer {};
437: template<typename Source_, typename Target_, typename Color_, typename SupportCompare_>
438: struct ArrowContainer<Source_, Target_, Color_, multiColor, SupportCompare_> {
439: // Define container's encapsulated types
440: typedef ArrowContainerTraits<Source_, Target_, Color_, SupportCompare_> traits;
441: // need to def arrow_type locally, since BOOST_MULTI_INDEX_MEMBER barfs when first template parameter starts with 'typename'
442: typedef typename traits::arrow_type arrow_type;
443: // Container set type
444: typedef ::boost::multi_index::multi_index_container<
445: typename traits::arrow_type,
446: ::boost::multi_index::indexed_by<
447: ::boost::multi_index::ordered_non_unique<
448: ::boost::multi_index::tag<typename traits::sourceTargetTag>,
449: ::boost::multi_index::composite_key<
450: typename traits::arrow_type,
451: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::source_type, source),
452: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::target_type, target),
453: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::color_type, color)
454: >
455: >,
456: ::boost::multi_index::ordered_non_unique<
457: ::boost::multi_index::tag<typename traits::sourceColorTag>,
458: ::boost::multi_index::composite_key<
459: typename traits::arrow_type,
460: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::source_type, source),
461: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::color_type, color),
462: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::target_type, target)
463: >
464: >,
465: ::boost::multi_index::ordered_non_unique<
466: ::boost::multi_index::tag<typename traits::targetColorTag>,
467: ::boost::multi_index::composite_key<
468: typename traits::arrow_type,
469: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::target_type, target),
470: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::color_type, color),
471: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::source_type, source)
472: >
473: >
474: >,
475: ALE_ALLOCATOR<typename traits::arrow_type>
476: > set_type;
477: // multi-index set of multicolor arrows
478: set_type set;
479: }; // class ArrowContainer<multiColor>
481: template<typename Source_, typename Target_, typename Color_, typename SupportCompare_>
482: struct ArrowContainer<Source_, Target_, Color_, uniColor, SupportCompare_> {
483: // Define container's encapsulated types
484: typedef ArrowContainerTraits<Source_, Target_, Color_, SupportCompare_> traits;
485: // need to def arrow_type locally, since BOOST_MULTI_INDEX_MEMBER barfs when first template parameter starts with 'typename'
486: typedef typename traits::arrow_type arrow_type;
488: // multi-index set type -- arrow set
489: typedef ::boost::multi_index::multi_index_container<
490: typename traits::arrow_type,
491: ::boost::multi_index::indexed_by<
492: ::boost::multi_index::ordered_unique<
493: ::boost::multi_index::tag<typename traits::sourceTargetTag>,
494: ::boost::multi_index::composite_key<
495: typename traits::arrow_type,
496: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::source_type, source),
497: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::target_type, target),
498: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::color_type, color)
499: >
500: >,
501: ::boost::multi_index::ordered_non_unique<
502: ::boost::multi_index::tag<typename traits::sourceColorTag>,
503: ::boost::multi_index::composite_key<
504: typename traits::arrow_type,
505: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::source_type, source),
506: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::color_type, color),
507: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::target_type, target)
508: >,
509: SupportCompare_
510: >,
511: ::boost::multi_index::ordered_non_unique<
512: ::boost::multi_index::tag<typename traits::targetColorTag>,
513: ::boost::multi_index::composite_key<
514: typename traits::arrow_type,
515: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::target_type, target),
516: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::color_type, color),
517: BOOST_MULTI_INDEX_MEMBER(arrow_type, typename traits::source_type, source)
518: >
519: >
520: >,
521: ALE_ALLOCATOR<typename traits::arrow_type>
522: > set_type;
523: // multi-index set of unicolor arrow records
524: set_type set;
525: }; // class ArrowContainer<uniColor>
526: }; // namespace SifterDef
528: //
529: // ASifter (short for Abstract Sifter, structurally a bipartite graph with colored arrows) implements a sequential interface
530: // similar to that of Sieve, except the source and target points may have different types and iterated operations (e.g., nCone,
531: // closure) are not available.
532: //
533: template<typename Source_, typename Target_, typename Color_, SifterDef::ColorMultiplicity colorMultiplicity, typename SupportCompare_ = ::boost::multi_index::composite_key_compare<std::less<Source_>, std::less<Color_>, std::less<Target_> >, typename SourceCtnr_ = SifterDef::RecContainer<Source_, SifterDef::Rec<Source_> >, typename TargetCtnr_ = SifterDef::RecContainer<Target_, SifterDef::Rec<Target_> > >
534: class ASifter { // class ASifter
535: public:
536: typedef struct {
537: typedef ASifter<Source_, Target_, Color_, colorMultiplicity, SupportCompare_, SourceCtnr_, TargetCtnr_> graph_type;
538: // Encapsulated container types
539: typedef SifterDef::ArrowContainer<Source_, Target_, Color_, colorMultiplicity, SupportCompare_> arrow_container_type;
540: typedef SourceCtnr_ cap_container_type;
541: typedef TargetCtnr_ base_container_type;
542: // Types associated with records held in containers
543: typedef typename arrow_container_type::traits::arrow_type arrow_type;
544: typedef typename arrow_container_type::traits::source_type source_type;
545: typedef typename cap_container_type::traits::rec_type sourceRec_type;
546: typedef typename arrow_container_type::traits::target_type target_type;
547: typedef typename base_container_type::traits::rec_type targetRec_type;
548: typedef typename arrow_container_type::traits::color_type color_type;
549: typedef typename arrow_container_type::traits::support_compare_type support_compare_type;
550: // Convenient tag names
551: typedef typename arrow_container_type::traits::sourceColorTag supportInd;
552: typedef typename arrow_container_type::traits::targetColorTag coneInd;
553: typedef typename arrow_container_type::traits::sourceTargetTag arrowInd;
554: typedef typename base_container_type::traits::pointTag baseInd;
555: typedef typename cap_container_type::traits::pointTag capInd;
556: //
557: // Return types
558: //
559: typedef typename
560: arrow_container_type::traits::template ArrowSequence<typename ::boost::multi_index::index<typename arrow_container_type::set_type,arrowInd>::type, source_type, target_type, BOOST_MULTI_INDEX_MEMBER(arrow_type, color_type, color)>
561: arrowSequence;
563: // FIX: This is a temp fix to include addArrow into the interface; should probably be pushed up to ArrowSequence
564: struct coneSequence : public arrow_container_type::traits::template ArrowSequence<typename ::boost::multi_index::index<typename arrow_container_type::set_type,coneInd>::type, target_type, color_type, BOOST_MULTI_INDEX_MEMBER(arrow_type, source_type, source)> {
565: protected:
566: graph_type& _graph;
567: public:
568: typedef typename
569: arrow_container_type::traits::template ArrowSequence<typename ::boost::multi_index::index<typename arrow_container_type::set_type,coneInd>::type, target_type, color_type, BOOST_MULTI_INDEX_MEMBER(arrow_type, source_type, source)> base_type;
570: // Encapsulated types
571: typedef typename base_type::traits traits;
572: typedef typename base_type::iterator iterator;
573: typedef typename base_type::reverse_iterator reverse_iterator;
574: // Basic interface
575: coneSequence(const coneSequence& seq) : base_type(seq), _graph(seq._graph) {};
576: coneSequence(graph_type& graph, typename traits::index_type& index, const typename base_type::key_type& k) : base_type(index, k), _graph(graph){};
577: coneSequence(graph_type& graph, typename traits::index_type& index, const typename base_type::key_type& k, const typename base_type::subkey_type& kk) : base_type(index, k, kk), _graph(graph) {};
578: virtual ~coneSequence() {};
580: // Extended interface
581: void addArrow(const arrow_type& a) {
582: // if(a.target != this->key) {
583: // throw ALE::Exception("Arrow target mismatch in a coneSequence");
584: // }
585: this->_graph.addArrow(a);
586: };
587: void addArrow(const source_type& s, const color_type& c){
588: this->_graph.addArrow(arrow_type(s,this->key,c));
589: };
591: virtual bool contains(const source_type& s) {
592: // Check whether a given point is in the index
593: typename ::boost::multi_index::index<typename ASifter::traits::arrow_container_type::set_type,typename ASifter::traits::arrowInd>::type& index = ::boost::multi_index::get<typename ASifter::traits::arrowInd>(this->_graph._arrows.set);
594: return (index.find(::boost::make_tuple(s,this->key)) != index.end());
595: };
596: };// struct coneSequence
598: // FIX: This is a temp fix to include addArrow into the interface; should probably be pushed up to ArrowSequence
599: struct supportSequence : public arrow_container_type::traits::template ArrowSequence<typename ::boost::multi_index::index<typename arrow_container_type::set_type,supportInd>::type, source_type, color_type, BOOST_MULTI_INDEX_MEMBER(arrow_type, target_type, target)> {
600: protected:
601: graph_type& _graph;
602: public:
603: typedef typename
604: arrow_container_type::traits::template ArrowSequence<typename ::boost::multi_index::index<typename arrow_container_type::set_type,supportInd>::type, source_type, color_type, BOOST_MULTI_INDEX_MEMBER(arrow_type, target_type, target)> base_type;
605: // Encapsulated types
606: typedef typename base_type::traits traits;
607: typedef typename base_type::iterator iterator;
608: typedef typename base_type::reverse_iterator reverse_iterator;
609: // Basic interface
610: supportSequence(const supportSequence& seq) : base_type(seq), _graph(seq._graph) {};
611: supportSequence(graph_type& graph, typename traits::index_type& index, const typename base_type::key_type& k) : base_type(index, k), _graph(graph){};
612: supportSequence(graph_type& graph, typename traits::index_type& index, const typename base_type::key_type& k, const typename base_type::subkey_type& kk) : base_type(index, k, kk), _graph(graph) {};
613: virtual ~supportSequence() {};
615: // FIX: WARNING: (or a HACK?): we flip the arrow on addition here.
616: // Fancy interface
617: void addArrow(const typename arrow_type::flip::type& af) {
618: this->_graph.addArrow(af.target, af.source, af.color);
619: };
620: void addArrow(const target_type& t, const color_type& c){
621: this->_graph.addArrow(arrow_type(this->key,t,c));
622: };
623: };// struct supportSequence
626: typedef typename base_container_type::traits::PointSequence baseSequence;
627: typedef typename cap_container_type::traits::PointSequence capSequence;
628: typedef std::set<source_type> coneSet;
629: typedef ALE::array<source_type> coneArray;
630: typedef std::set<target_type> supportSet;
631: typedef ALE::array<target_type> supportArray;
632: } traits;
634: template <typename OtherSource_, typename OtherTarget_, typename OtherColor_, SifterDef::ColorMultiplicity otherColorMultiplicity,
635: typename OtherSupportCompare_ = ::boost::multi_index::composite_key_compare<std::less<OtherSource_>, std::less<OtherColor_>, std::less<OtherTarget_> >,
636: typename OtherSourceCtnr_ = SifterDef::RecContainer<OtherSource_, SifterDef::Rec<OtherSource_> >,
637: typename OtherTargetCtnr_ = SifterDef::RecContainer<OtherTarget_, SifterDef::Rec<OtherTarget_> > >
638: struct rebind {
639: typedef ASifter<OtherSource_, OtherTarget_, OtherColor_, otherColorMultiplicity, OtherSupportCompare_, OtherSourceCtnr_, OtherTargetCtnr_> type;
640: };
642: public:
643: // Debug level
644: int _debug;
645: //protected:
646: typename traits::arrow_container_type _arrows;
647: typename traits::base_container_type _base;
648: typename traits::cap_container_type _cap;
649: protected:
650: MPI_Comm _comm;
651: int _commRank;
652: int _commSize;
653: PetscObject _petscObj;
654: void __init(MPI_Comm comm) {
655: static PetscClassId sifterType = -1;
656: //const char *id_name = ALE::getClassName<T>();
657: const char *id_name = "Sifter";
658: PetscErrorCode ierr;
660: if (sifterType < 0) {
661: PetscClassIdRegister(id_name,&sifterType);CHKERROR(ierr, "Error in MPI_Comm_rank");
662: }
663: this->_comm = comm;
664: MPI_Comm_rank(this->_comm, &this->_commRank);CHKERROR(ierr, "Error in MPI_Comm_rank");
665: MPI_Comm_size(this->_comm, &this->_commSize);CHKERROR(ierr, "Error in MPI_Comm_rank");
666: #ifdef USE_PETSC_OBJ
667: PetscObjectCreateGeneric(this->_comm, sifterType, id_name, &this->_petscObj);CHKERROR(ierr, "Error in PetscObjectCreate");
668: #endif
669: //ALE::restoreClassName<T>(id_name);
670: };
671: // We store these sequence objects to avoid creating them each query
672: Obj<typename traits::coneSequence> _coneSeq;
673: Obj<typename traits::supportSequence> _supportSeq;
674: public:
675: //
676: // Basic interface
677: //
678: ASifter(MPI_Comm comm = PETSC_COMM_SELF, const int& debug = 0) : _debug(debug), _petscObj(NULL) {
679: __init(comm);
680: this->_coneSeq = new typename traits::coneSequence(*this, ::boost::multi_index::get<typename traits::coneInd>(this->_arrows.set), typename traits::target_type());
681: this->_supportSeq = new typename traits::supportSequence(*this, ::boost::multi_index::get<typename traits::supportInd>(this->_arrows.set), typename traits::source_type());
682: }
683: virtual ~ASifter() {
684: #ifdef USE_PETSC_OBJ
686: PetscObjectDestroy(&this->_petscObj);CHKERROR(ierr, "Failed in PetscObjectDestroy");
687: #endif
688: };
689: //
690: // Query methods
691: //
692: int debug() const {return this->_debug;};
693: void setDebug(const int debug) {this->_debug = debug;};
694: MPI_Comm comm() const {return this->_comm;};
695: int commSize() const {return this->_commSize;};
696: int commRank() const {return this->_commRank;}
697: #ifdef USE_PETSC_OBJ
698: PetscObject petscObj() const {return this->_petscObj;};
699: #endif
701: // Added to allow optimized versions
702: void assemble() {};
703: void assemblePoints() {};
704: // FIX: need const_cap, const_base returning const capSequence etc, but those need to have const_iterators, const_begin etc.
705: Obj<typename traits::capSequence> cap() {
706: return typename traits::capSequence(::boost::multi_index::get<typename traits::capInd>(this->_cap.set));
707: };
708: Obj<typename traits::baseSequence> base() {
709: return typename traits::baseSequence(::boost::multi_index::get<typename traits::baseInd>(this->_base.set));
710: };
711: typename traits::capSequence::iterator capBegin() {
712: return this->cap()->begin();
713: };
714: typename traits::capSequence::iterator capEnd() {
715: return this->cap()->end();
716: };
717: typename traits::baseSequence::iterator baseBegin() {
718: return this->base()->begin();
719: };
720: typename traits::baseSequence::iterator baseEnd() {
721: return this->base()->end();
722: };
723: int getBaseSize() {return this->base()->size();};
724: void setBaseSize(int size) {};
725: int getCapSize() {return this->cap()->size();};
726: void setCapSize(int size) {};
727: bool capContains(const typename traits::source_type& p) {
728: typename traits::capSequence cap(::boost::multi_index::get<typename traits::capInd>(this->_cap.set));
730: //for(typename traits::capSequence::iterator c_iter = cap.begin(); c_iter != cap.end(); ++c_iter) {
731: //}
732: return cap.contains(p);
733: };
734: bool baseContains(const typename traits::target_type& p) {
735: typename traits::baseSequence base(::boost::multi_index::get<typename traits::baseInd>(this->_base.set));
737: //for(typename traits::capSequence::iterator c_iter = cap.begin(); c_iter != cap.end(); ++c_iter) {
738: //}
739: return base.contains(p);
740: };
741: // FIX: should probably have cone and const_cone etc, since arrows can be modified through an iterator (modifyColor).
742: Obj<typename traits::arrowSequence>
743: arrows(const typename traits::source_type& s, const typename traits::target_type& t) {
744: return typename traits::arrowSequence(::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set), s, t);
745: };
746: Obj<typename traits::arrowSequence>
747: arrows(const typename traits::source_type& s) {
748: return typename traits::arrowSequence(::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set), s);
749: };
750: #ifdef SLOW
751: Obj<typename traits::coneSequence>
752: cone(const typename traits::target_type& p) {
753: return typename traits::coneSequence(*this, ::boost::multi_index::get<typename traits::coneInd>(this->_arrows.set), p);
754: };
755: #else
756: const Obj<typename traits::coneSequence>&
757: cone(const typename traits::target_type& p) {
758: this->_coneSeq->setKey(p);
759: this->_coneSeq->setUseSubkey(false);
760: return this->_coneSeq;
761: };
762: const typename traits::coneSequence::iterator
763: coneBegin(const typename traits::target_type& p) {
764: return this->cone(p)->begin();
765: };
766: const typename traits::coneSequence::iterator
767: coneEnd(const typename traits::target_type& p) {
768: return this->cone(p)->end();
769: };
770: void setConeSize(const typename traits::target_type& p, int size) {};
771: #endif
772: template<class InputSequence>
773: Obj<typename traits::coneSet>
774: cone(const Obj<InputSequence>& points) {
775: return this->cone(points, typename traits::color_type(), false);
776: }
777: #ifdef SLOW
778: Obj<typename traits::coneSequence>
779: cone(const typename traits::target_type& p, const typename traits::color_type& color) {
780: return typename traits::coneSequence(*this, ::boost::multi_index::get<typename traits::coneInd>(this->_arrows.set), p, color);
781: }
782: #else
783: const Obj<typename traits::coneSequence>&
784: cone(const typename traits::target_type& p, const typename traits::color_type& color) {
785: this->_coneSeq->setKey(p);
786: this->_coneSeq->setSubkey(color);
787: this->_coneSeq->setUseSubkey(true);
788: return this->_coneSeq;
789: };
790: #endif
791: template<class InputSequence>
792: Obj<typename traits::coneSet>
793: cone(const Obj<InputSequence>& points, const typename traits::color_type& color, bool useColor = true) {
794: Obj<typename traits::coneSet> cone = typename traits::coneSet();
795: for(typename InputSequence::iterator p_itor = points->begin(); p_itor != points->end(); ++p_itor) {
796: Obj<typename traits::coneSequence> pCone;
797: if (useColor) {
798: pCone = this->cone(*p_itor, color);
799: } else {
800: pCone = this->cone(*p_itor);
801: }
802: cone->insert(pCone->begin(), pCone->end());
803: }
804: return cone;
805: }
806: template<typename PointCheck>
807: bool coneContains(const typename traits::target_type& p, const PointCheck& checker) {
808: typename traits::coneSequence cone(*this, ::boost::multi_index::get<typename traits::coneInd>(this->_arrows.set), p);
810: for(typename traits::coneSequence::iterator c_iter = cone.begin(); c_iter != cone.end(); ++c_iter) {
811: if (checker(*c_iter, p)) return true;
812: }
813: return false;
814: }
815: template<typename PointProcess>
816: void coneApply(const typename traits::target_type& p, PointProcess& processor) {
817: typename traits::coneSequence cone(*this, ::boost::multi_index::get<typename traits::coneInd>(this->_arrows.set), p);
819: for(typename traits::coneSequence::iterator c_iter = cone.begin(); c_iter != cone.end(); ++c_iter) {
820: processor(*c_iter, p);
821: }
822: }
823: int getConeSize(const typename traits::target_type& p) {
824: return this->cone(p)->size();
825: }
826: #ifdef SLOW
827: Obj<typename traits::supportSequence>
828: support(const typename traits::source_type& p) {
829: return typename traits::supportSequence(*this, ::boost::multi_index::get<typename traits::supportInd>(this->_arrows.set), p);
830: }
831: #else
832: const Obj<typename traits::supportSequence>&
833: support(const typename traits::source_type& p) {
834: this->_supportSeq->setKey(p);
835: this->_supportSeq->setUseSubkey(false);
836: return this->_supportSeq;
837: };
838: const typename traits::supportSequence::iterator
839: supportBegin(const typename traits::source_type& p) {
840: return this->support(p)->begin();
841: };
842: const typename traits::supportSequence::iterator
843: supportEnd(const typename traits::source_type& p) {
844: return this->support(p)->end();
845: };
846: const typename traits::supportSequence::iterator
847: supportBegin(const typename traits::source_type& p, const typename traits::color_type& color) {
848: return this->support(p, color)->begin();
849: };
850: const typename traits::supportSequence::iterator
851: supportEnd(const typename traits::source_type& p, const typename traits::color_type& color) {
852: return this->support(p, color)->end();
853: };
854: void setSupportSize(const typename traits::source_type& p, int size) {};
855: #endif
856: #ifdef SLOW
857: Obj<typename traits::supportSequence>
858: support(const typename traits::source_type& p, const typename traits::color_type& color) {
859: return typename traits::supportSequence(*this, ::boost::multi_index::get<typename traits::supportInd>(this->_arrows.set), p, color);
860: };
861: #else
862: const Obj<typename traits::supportSequence>&
863: support(const typename traits::source_type& p, const typename traits::color_type& color) {
864: this->_supportSeq->setKey(p);
865: this->_supportSeq->setSubkey(color);
866: this->_supportSeq->setUseSubkey(true);
867: return this->_supportSeq;
868: };
869: #endif
870: template<class InputSequence>
871: Obj<typename traits::supportSet>
872: support(const Obj<InputSequence>& sources) {
873: return this->support(sources, typename traits::color_type(), false);
874: }
875: template<class InputSequence>
876: Obj<typename traits::supportSet>
877: support(const Obj<InputSequence>& points, const typename traits::color_type& color, bool useColor = true){
878: Obj<typename traits::supportSet> supp = typename traits::supportSet();
879: for(typename InputSequence::iterator p_itor = points->begin(); p_itor != points->end(); ++p_itor) {
880: Obj<typename traits::supportSequence> pSupport;
881: if (useColor) {
882: pSupport = this->support(*p_itor, color);
883: } else {
884: pSupport = this->support(*p_itor);
885: }
886: supp->insert(pSupport->begin(), pSupport->end());
887: }
888: return supp;
889: }
890: template<typename PointCheck>
891: bool supportContains(const typename traits::source_type& p, const PointCheck& checker) {
892: typename traits::supportSequence support(*this, ::boost::multi_index::get<typename traits::supportInd>(this->_arrows.set), p);
894: for(typename traits::supportSequence::iterator s_iter = support.begin(); s_iter != support.end(); ++s_iter) {
895: if (checker(*s_iter, p)) return true;
896: }
897: return false;
898: }
899: template<typename PointProcess>
900: void supportApply(const typename traits::source_type& p, PointProcess& processor) {
901: typename traits::supportSequence support(*this, ::boost::multi_index::get<typename traits::supportInd>(this->_arrows.set), p);
903: for(typename traits::supportSequence::iterator s_iter = support.begin(); s_iter != support.end(); ++s_iter) {
904: processor(*s_iter, p);
905: }
906: }
907: int getSupportSize(const typename traits::source_type& p) {
908: return this->support(p)->size();
909: }
910: int getSupportSize(const typename traits::source_type& p, const typename traits::color_type& color) {
911: return this->support(p, color)->size();
912: }
914: template<typename ostream_type>
915: void view(ostream_type& os, const char* label = NULL, bool rawData = false){
916: int rank = this->commRank();
918: if(label != NULL) {
919: os << "["<<rank<<"]Viewing Sifter '" << label << "':" << std::endl;
920: }
921: else {
922: os << "["<<rank<<"]Viewing a Sifter:" << std::endl;
923: }
924: if(!rawData) {
925: os << "cap --> base:" << std::endl;
926: Obj<typename traits::capSequence> cap = this->cap();
927: for(typename traits::capSequence::iterator capi = cap->begin(); capi != cap->end(); capi++) {
928: const Obj<typename traits::supportSequence>& supp = this->support(*capi);
930: for(typename traits::supportSequence::iterator suppi = supp->begin(); suppi != supp->end(); suppi++) {
931: os << *capi << "--(" << suppi.color() << ")-->" << *suppi << std::endl;
932: }
933: }
934: os << "base <-- cap:" << std::endl;
935: Obj<typename traits::baseSequence> base = this->base();
936: for(typename traits::baseSequence::iterator basei = base->begin(); basei != base->end(); basei++) {
937: const Obj<typename traits::coneSequence>& cone = this->cone(*basei);
939: for(typename traits::coneSequence::iterator conei = cone->begin(); conei != cone->end(); conei++) {
940: os << *basei << "<--(" << conei.color() << ")--" << *conei << std::endl;
941: }
942: }
943: os << "cap --> outdegrees:" << std::endl;
944: for(typename traits::capSequence::iterator capi = cap->begin(); capi != cap->end(); capi++) {
945: os << *capi << "-->" << capi.degree() << std::endl;
946: }
947: os << "base <-- indegrees:" << std::endl;
948: for(typename traits::baseSequence::iterator basei = base->begin(); basei != base->end(); basei++) {
949: os << *basei << "<--" << basei.degree() << std::endl;
950: }
951: }
952: else {
953: os << "'raw' arrow set:" << std::endl;
954: for(typename traits::arrow_container_type::set_type::iterator ai = _arrows.set.begin(); ai != _arrows.set.end(); ai++)
955: {
956: typename traits::arrow_type arr = *ai;
957: os << arr << std::endl;
958: }
959: os << "'raw' base set:" << std::endl;
960: for(typename traits::base_container_type::set_type::iterator bi = _base.set.begin(); bi != _base.set.end(); bi++)
961: {
962: typename traits::base_container_type::traits::rec_type bp = *bi;
963: os << bp << std::endl;
964: }
965: os << "'raw' cap set:" << std::endl;
966: for(typename traits::cap_container_type::set_type::iterator ci = _cap.set.begin(); ci != _cap.set.end(); ci++)
967: {
968: typename traits::cap_container_type::traits::rec_type cp = *ci;
969: os << cp << std::endl;
970: }
971: }
972: }
973: // A parallel viewer
976: PetscErrorCode view(const char* label = NULL, bool raw = false){
978: ostringstream txt;
980: if(this->_debug) {
981: std::cout << "viewing a Sifter, comm = " << this->comm() << ", PETSC_COMM_SELF = " << PETSC_COMM_SELF << ", commRank = " << this->commRank() << std::endl;
982: }
983: if(label != NULL) {
984: PetscPrintf(this->comm(), "viewing Sifter: '%s'\n", label);
985: } else {
986: PetscPrintf(this->comm(), "viewing a Sifter: \n");
987: }
988: if(!raw) {
989: ostringstream txt;
990: if(this->commRank() == 0) {
991: txt << "cap --> base:\n";
992: }
993: typename traits::capSequence cap = this->cap();
994: typename traits::baseSequence base = this->base();
995: if(cap.empty()) {
996: txt << "[" << this->commRank() << "]: empty" << std::endl;
997: }
998: for(typename traits::capSequence::iterator capi = cap.begin(); capi != cap.end(); capi++) {
999: const Obj<typename traits::supportSequence>& supp = this->support(*capi);
1000: for(typename traits::supportSequence::iterator suppi = supp->begin(); suppi != supp->end(); suppi++) {
1001: txt << "[" << this->commRank() << "]: " << *capi << "--(" << suppi.color() << ")-->" << *suppi << std::endl;
1002: }
1003: }
1004: //
1005: PetscSynchronizedPrintf(this->comm(), txt.str().c_str());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1006: PetscSynchronizedFlush(this->comm());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1007: #if 0
1008: //
1009: ostringstream txt1;
1010: if(this->commRank() == 0) {
1011: //txt1 << "cap <point,degree>:\n";
1012: txt1 << "cap:\n";
1013: }
1014: txt1 << "[" << this->commRank() << "]: [";
1015: for(typename traits::capSequence::iterator capi = cap.begin(); capi != cap.end(); capi++) {
1016: //txt1 << " <" << *capi << "," << capi.degree() << ">";
1017: txt1 << " " << *capi;
1018: }
1019: txt1 << " ]" << std::endl;
1020: //
1021: PetscSynchronizedPrintf(this->comm(), txt1.str().c_str());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1022: PetscSynchronizedFlush(this->comm());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1023: //
1024: ostringstream txt2;
1025: if(this->commRank() == 0) {
1026: //txt2 << "base <point,degree>:\n";
1027: txt2 << "base:\n";
1028: }
1029: txt2 << "[" << this->commRank() << "]: [";
1030: for(typename traits::baseSequence::iterator basei = base.begin(); basei != base.end(); basei++) {
1031: txt2 << " " << *basei;
1032: //txt2 << " <" << *basei << "," << basei.degree() << ">";
1033: }
1034: txt2 << " ]" << std::endl;
1035: //
1036: PetscSynchronizedPrintf(this->comm(), txt2.str().c_str());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1037: PetscSynchronizedFlush(this->comm());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1038: #endif
1039: }
1040: else { // if(raw)
1041: ostringstream txt;
1042: if(this->commRank() == 0) {
1043: txt << "'raw' arrow set:" << std::endl;
1044: }
1045: for(typename traits::arrow_container_type::set_type::iterator ai = _arrows.set.begin(); ai != _arrows.set.end(); ai++)
1046: {
1047: typename traits::arrow_type arr = *ai;
1048: txt << "[" << this->commRank() << "]: " << arr << std::endl;
1049: }
1050: PetscSynchronizedPrintf(this->comm(), txt.str().c_str());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1051: PetscSynchronizedFlush(this->comm());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1052: //
1053: ostringstream txt1;
1054: if(this->commRank() == 0) {
1055: txt1 << "'raw' base set:" << std::endl;
1056: }
1057: for(typename traits::base_container_type::set_type::iterator bi = _base.set.begin(); bi != _base.set.end(); bi++)
1058: {
1059: typename traits::base_container_type::traits::rec_type bp = *bi;
1060: txt1 << "[" << this->commRank() << "]: " << bp << std::endl;
1061: }
1062: PetscSynchronizedPrintf(this->comm(), txt1.str().c_str());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1063: PetscSynchronizedFlush(this->comm());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1064: //
1065: ostringstream txt2;
1066: if(this->commRank() == 0) {
1067: txt2 << "'raw' cap set:" << std::endl;
1068: }
1069: for(typename traits::cap_container_type::set_type::iterator ci = _cap.set.begin(); ci != _cap.set.end(); ci++)
1070: {
1071: typename traits::cap_container_type::traits::rec_type cp = *ci;
1072: txt2 << "[" << this->commRank() << "]: " << cp << std::endl;
1073: }
1074: PetscSynchronizedPrintf(this->comm(), txt2.str().c_str());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1075: PetscSynchronizedFlush(this->comm());CHKERROR(ierr, "Error in PetscSynchronizedFlush");
1076: }// if(raw)
1078: return(0);
1079: };
1080: public:
1081: void copy(ASifter *newSifter) {
1082: const typename traits::baseSequence::iterator sBegin = this->baseBegin();
1083: const typename traits::baseSequence::iterator sEnd = this->baseEnd();
1085: for(typename traits::baseSequence::iterator r_iter = sBegin; r_iter != sEnd; ++r_iter) {
1086: const typename traits::coneSequence::iterator pBegin = this->coneBegin(*r_iter);
1087: const typename traits::coneSequence::iterator pEnd = this->coneEnd(*r_iter);
1089: for(typename traits::coneSequence::iterator p_iter = pBegin; p_iter != pEnd; ++p_iter) {
1090: const Obj<typename traits::supportSequence>& support = this->support(*p_iter);
1091: const typename traits::supportSequence::iterator supBegin = support->begin();
1092: const typename traits::supportSequence::iterator supEnd = support->end();
1094: for(typename traits::supportSequence::iterator s_iter = supBegin; s_iter != supEnd; ++s_iter) {
1095: newSifter->addArrow(*p_iter, *s_iter, s_iter.color());
1096: }
1097: }
1098: }
1099: };
1100: //
1101: // Lattice queries
1102: //
1103: template<class targetInputSequence>
1104: Obj<typename traits::coneSequence> meet(const Obj<targetInputSequence>& targets);
1105: // unimplemented
1106: template<class targetInputSequence>
1107: Obj<typename traits::coneSequence> meet(const Obj<targetInputSequence>& targets, const typename traits::color_type& color);
1108: // unimplemented
1109: template<class sourceInputSequence>
1110: Obj<typename traits::coneSequence> join(const Obj<sourceInputSequence>& sources);
1111: // unimplemented
1112: template<class sourceInputSequence>
1113: Obj<typename traits::coneSequence> join(const Obj<sourceInputSequence>& sources, const typename traits::color_type& color);
1114: public:
1115: //
1116: // Structural manipulation
1117: //
1118: void clear() {
1119: this->_arrows.set.clear(); this->_base.set.clear(); this->_cap.set.clear();
1120: };
1121: void addBasePoint(const typename traits::target_type t) {
1122: /* // Increase degree by 0, which won't affect an existing point and will insert a new point, if necessery
1123: this->_base.adjustDegree(t,0); */
1124: this->_base.set.insert(typename traits::targetRec_type(t,0));
1125: };
1126: void addBasePoint(const typename traits::targetRec_type b) {
1127: this->_base.set.insert(b);
1128: };
1129: void removeBasePoint(const typename traits::target_type t) {
1130: if (this->_debug) {std::cout << " Removing " << t << " from the base" << std::endl;}
1131: // Clear the cone and remove the point from _base
1132: this->clearCone(t);
1133: this->_base.removePoint(t);
1134: };
1135: void addCapPoint(const typename traits::source_type s) {
1136: /* // Increase degree by 0, which won't affect an existing point and will insert a new point, if necessery
1137: this->_cap.adjustDegree(s,0); */
1138: this->_cap.set.insert(typename traits::sourceRec_type(s,0));
1139: };
1140: void addCapPoint(const typename traits::sourceRec_type c) {
1141: this->_cap.set.insert(c);
1142: };
1143: void removeCapPoint(const typename traits::source_type s) {
1144: if (this->_debug) {std::cout << " Removing " << s << " from the cap" << std::endl;}
1145: // Clear the support and remove the point from _cap
1146: this->clearSupport(s);
1147: this->_cap.removePoint(s);
1148: };
1149: virtual void addArrow(const typename traits::source_type& p, const typename traits::target_type& q) {
1150: this->addArrow(p, q, typename traits::color_type());
1151: };
1152: virtual void addArrow(const typename traits::source_type& p, const typename traits::target_type& q, const typename traits::color_type& color) {
1153: this->addArrow(typename traits::arrow_type(p, q, color));
1154: //std::cout << "Added " << arrow_type(p, q, color);
1155: };
1156: virtual bool checkArrow(const typename traits::arrow_type& a) {
1157: if (this->_cap.set.find(a.source) == this->_cap.set.end()) return false;
1158: if (this->_base.set.find(a.target) == this->_base.set.end()) return false;
1159: return true;
1160: };
1161: virtual void addArrow(const typename traits::arrow_type& a, bool noNewPoints = false) {
1162: if (noNewPoints && !this->checkArrow(a)) return;
1163: this->_arrows.set.insert(a);
1164: this->addBasePoint(a.target);
1165: this->addCapPoint(a.source);
1166: };
1167: virtual void removeArrow(const typename traits::source_type& p, const typename traits::target_type& q) {
1168: this->removeArrow(typename traits::arrow_type(p, q, typename traits::color_type()));
1169: };
1170: virtual void removeArrow(const typename traits::arrow_type& a) {
1171: // First, produce an arrow sequence for the given source, target combination.
1172: typename traits::arrowSequence::traits::index_type& arrowIndex =
1173: ::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set);
1174: typename traits::arrowSequence::traits::index_type::iterator i,ii,j;
1175: i = arrowIndex.lower_bound(::boost::make_tuple(a.source,a.target));
1176: ii = arrowIndex.upper_bound(::boost::make_tuple(a.source, a.target));
1177: if (this->_debug) {
1178: std::cout << "removeArrow: attempting to remove arrow:" << a << std::endl;
1179: std::cout << "removeArrow: candidate arrows are:" << std::endl;
1180: }
1181: for(j = i; j != ii; j++) {
1182: if (this->_debug) {
1183: std::cout << " " << *j;
1184: }
1185: // Find the arrow of right color and remove it
1186: if(j->color == a.color) {
1187: if (this->_debug) {
1188: std::cout << std::endl << "removeArrow: found:" << *j << std::endl;
1189: }
1190: /* this->_base.adjustDegree(a.target, -1); this->_cap.adjustDegree(a.source,-1); */
1191: arrowIndex.erase(j);
1192: break;
1193: }
1194: }
1195: };
1197: void addCone(const typename traits::source_type& source, const typename traits::target_type& target){
1198: this->addArrow(source, target);
1199: }
1200: template<class sourceInputSequence>
1201: void addCone(const Obj<sourceInputSequence>& sources, const typename traits::target_type& target) {
1202: this->addCone(sources, target, typename traits::color_type());
1203: }
1204: void addCone(const typename traits::source_type& source, const typename traits::target_type& target, const typename traits::color_type& color) {
1205: this->addArrow(source, target, color);
1206: }
1207: template<class sourceInputSequence>
1208: void
1209: addCone(const Obj<sourceInputSequence>& sources, const typename traits::target_type& target, const typename traits::color_type& color){
1210: if (this->_debug > 1) {std::cout << "Adding a cone " << std::endl;}
1211: for(typename sourceInputSequence::iterator iter = sources->begin(); iter != sources->end(); ++iter) {
1212: if (this->_debug > 1) {std::cout << "Adding arrow from " << *iter << " to " << target << "(" << color << ")" << std::endl;}
1213: this->addArrow(*iter, target, color);
1214: }
1215: }
1216: void clearCone(const typename traits::target_type& t) {
1217: clearCone(t, typename traits::color_type(), false);
1218: };
1220: void clearCone(const typename traits::target_type& t, const typename traits::color_type& color, bool useColor = true) {
1221: // Use the cone sequence types to clear the cone
1222: typename traits::coneSequence::traits::index_type& coneIndex =
1223: ::boost::multi_index::get<typename traits::coneInd>(this->_arrows.set);
1224: typename traits::coneSequence::traits::index_type::iterator i, ii, j;
1225: if (this->_debug > 20) {
1226: std::cout << "clearCone: removing cone over " << t;
1227: if(useColor) {
1228: std::cout << " with color" << color << std::endl;
1229: const Obj<typename traits::coneSequence>& cone = this->cone(t,color);
1230: std::cout << "[";
1231: for(typename traits::coneSequence::iterator ci = cone->begin(); ci != cone->end(); ci++) {
1232: std::cout << " " << ci.arrow();
1233: }
1234: std::cout << "]" << std::endl;
1235: }
1236: else {
1237: std::cout << std::endl;
1238: const Obj<typename traits::coneSequence>& cone = this->cone(t);
1239: std::cout << "[";
1240: for(typename traits::coneSequence::iterator ci = cone->begin(); ci != cone->end(); ci++) {
1241: std::cout << " " << ci.arrow();
1242: }
1243: std::cout << "]" << std::endl;
1244: }
1245: }
1246: if (useColor) {
1247: i = coneIndex.lower_bound(::boost::make_tuple(t,color));
1248: ii = coneIndex.upper_bound(::boost::make_tuple(t,color));
1249: } else {
1250: i = coneIndex.lower_bound(::boost::make_tuple(t));
1251: ii = coneIndex.upper_bound(::boost::make_tuple(t));
1252: }
1253: for(j = i; j != ii; j++){
1254: // Adjust the degrees before removing the arrow; use a different iterator, since we'll need i,ii to do the arrow erasing.
1255: if (this->_debug) {
1256: std::cout << "clearCone: adjusting degrees for endpoints of arrow: " << *j << std::endl;
1257: }
1258: /* this->_cap.adjustDegree(j->source, -1);
1259: this->_base.adjustDegree(j->target, -1); */
1260: }
1261: coneIndex.erase(i,ii);
1262: };// clearCone()
1264: template<class InputSequence>
1265: void
1266: restrictBase(const Obj<InputSequence>& points) {
1267: typename traits::baseSequence base = this->base();
1268: typename std::set<typename traits::target_type> remove;
1270: for(typename traits::baseSequence::iterator bi = base.begin(); bi != base.end(); bi++) {
1271: // Check whether *bi is in points, if it is NOT, remove it
1272: // if (!points->contains(*bi)) {
1273: if (points->find(*bi) == points->end()) {
1274: // this->removeBasePoint(*bi);
1275: remove.insert(*bi);
1276: }
1277: }
1278: //FIX
1279: for(typename std::set<typename traits::target_type>::iterator r_iter = remove.begin(); r_iter != remove.end(); ++r_iter) {
1280: this->removeBasePoint(*r_iter);
1281: }
1282: }
1284: template<class InputSequence>
1285: void
1286: excludeBase(const Obj<InputSequence>& points) {
1287: for(typename InputSequence::iterator pi = points->begin(); pi != points->end(); pi++) {
1288: this->removeBasePoint(*pi);
1289: }
1290: }
1292: template<class InputSequence>
1293: void
1294: restrictCap(const Obj<InputSequence>& points) {
1295: typename traits::capSequence cap = this->cap();
1296: for(typename traits::capSequence::iterator ci = cap.begin(); ci != cap.end(); ci++) {
1297: // Check whether *ci is in points, if it is NOT, remove it
1298: if(points->find(*ci) == points->end()) {
1299: this->removeCapPoint(*ci);
1300: }
1301: }
1302: }
1304: template<class InputSequence>
1305: void
1306: excludeCap(const Obj<InputSequence>& points) {
1307: for(typename InputSequence::iterator pi = points->begin(); pi != points->end(); pi++) {
1308: this->removeCapPoint(*pi);
1309: }
1310: }
1312: void clearSupport(const typename traits::source_type& s) {
1313: clearSupport(s, typename traits::color_type(), false);
1314: };
1315: void clearSupport(const typename traits::source_type& s, const typename traits::color_type& color, bool useColor = true) {
1316: // Use the cone sequence types to clear the cone
1317: typename
1318: traits::supportSequence::traits::index_type& suppIndex = ::boost::multi_index::get<typename traits::supportInd>(this->_arrows.set);
1319: typename traits::supportSequence::traits::index_type::iterator i, ii, j;
1320: if (useColor) {
1321: i = suppIndex.lower_bound(::boost::make_tuple(s,color));
1322: ii = suppIndex.upper_bound(::boost::make_tuple(s,color));
1323: } else {
1324: i = suppIndex.lower_bound(::boost::make_tuple(s));
1325: ii = suppIndex.upper_bound(::boost::make_tuple(s));
1326: }
1327: for(j = i; j != ii; j++){
1328: // Adjust the degrees before removing the arrow
1329: /* this->_cap.adjustDegree(j->source, -1);
1330: this->_base.adjustDegree(j->target, -1); */
1331: }
1332: suppIndex.erase(i,ii);
1333: }
1334: void setCone(const typename traits::source_type& source, const typename traits::target_type& target){
1335: this->clearCone(target, typename traits::color_type(), false); this->addCone(source, target);
1336: }
1337: template<class sourceInputSequence>
1338: void setCone(const Obj<sourceInputSequence>& sources, const typename traits::target_type& target) {
1339: this->clearCone(target, typename traits::color_type(), false); this->addCone(sources, target, typename traits::color_type());
1340: }
1341: void setCone(const typename traits::source_type& source, const typename traits::target_type& target, const typename traits::color_type& color) {
1342: this->clearCone(target, color, true); this->addCone(source, target, color);
1343: }
1344: template<class sourceInputSequence>
1345: void setCone(const Obj<sourceInputSequence>& sources, const typename traits::target_type& target, const typename traits::color_type& color){
1346: this->clearCone(target, color, true); this->addCone(sources, target, color);
1347: }
1348: template<class targetInputSequence>
1349: void addSupport(const typename traits::source_type& source, const Obj<targetInputSequence >& targets) {
1350: this->addSupport(source, targets, typename traits::color_type());
1351: }
1352: template<class targetInputSequence>
1353: void addSupport(const typename traits::source_type& source, const Obj<targetInputSequence>& targets, const typename traits::color_type& color) {
1354: const typename targetInputSequence::iterator end = targets->end();
1356: for(typename targetInputSequence::iterator iter = targets->begin(); iter != end; ++iter) {
1357: this->addArrow(source, *iter, color);
1358: }
1359: }
1360: template<typename Sifter_>
1361: void add(const Obj<Sifter_>& cbg, bool noNewPoints = false) {
1362: typename ::boost::multi_index::index<typename Sifter_::traits::arrow_container_type::set_type, typename Sifter_::traits::arrowInd>::type& aInd = ::boost::multi_index::get<typename Sifter_::traits::arrowInd>(cbg->_arrows.set);
1364: for(typename ::boost::multi_index::index<typename Sifter_::traits::arrow_container_type::set_type, typename Sifter_::traits::arrowInd>::type::iterator a_iter = aInd.begin(); a_iter != aInd.end(); ++a_iter) {
1365: this->addArrow(*a_iter, noNewPoints);
1366: }
1367: if (!noNewPoints) {
1368: typename ::boost::multi_index::index<typename Sifter_::traits::base_container_type::set_type, typename Sifter_::traits::baseInd>::type& bInd = ::boost::multi_index::get<typename Sifter_::traits::baseInd>(this->_base.set);
1370: for(typename ::boost::multi_index::index<typename Sifter_::traits::base_container_type::set_type, typename Sifter_::traits::baseInd>::type::iterator b_iter = bInd.begin(); b_iter != bInd.end(); ++b_iter) {
1371: this->addBasePoint(*b_iter);
1372: }
1373: typename ::boost::multi_index::index<typename Sifter_::traits::cap_container_type::set_type, typename Sifter_::traits::capInd>::type& cInd = ::boost::multi_index::get<typename Sifter_::traits::capInd>(this->_cap.set);
1375: for(typename ::boost::multi_index::index<typename Sifter_::traits::cap_container_type::set_type, typename Sifter_::traits::capInd>::type::iterator c_iter = cInd.begin(); c_iter != cInd.end(); ++c_iter) {
1376: this->addCapPoint(*c_iter);
1377: }
1378: }
1379: }
1380: }; // class ASifter
1382: // A UniSifter aka Sifter
1383: template <typename Source_, typename Target_, typename Color_,
1384: typename SupportCompare_ = ::boost::multi_index::composite_key_compare<std::less<Source_>, std::less<Color_>, std::less<Target_> >,
1385: typename SourceCtnr_ = SifterDef:: RecContainer<Source_, SifterDef::Rec<Source_> >, typename TargetCtnr_= SifterDef::RecContainer<Target_, SifterDef::Rec<Target_> > >
1386: class Sifter : public ASifter<Source_, Target_, Color_, SifterDef::uniColor, SupportCompare_, SourceCtnr_, TargetCtnr_> {
1387: public:
1388: typedef typename ASifter<Source_, Target_, Color_, SifterDef::uniColor, SupportCompare_, SourceCtnr_, TargetCtnr_>::traits traits;
1389: template <typename OtherSource_, typename OtherTarget_, typename OtherColor_,
1390: typename OtherSupportCompare_ = ::boost::multi_index::composite_key_compare<std::less<OtherSource_>, std::less<OtherColor_>, std::less<OtherTarget_> >,
1391: typename OtherSourceCtnr_ = SifterDef::RecContainer<OtherSource_, SifterDef::Rec<OtherSource_> >,
1392: typename OtherTargetCtnr_ = SifterDef::RecContainer<OtherTarget_, SifterDef::Rec<OtherTarget_> > >
1393: struct rebind {
1394: typedef Sifter<OtherSource_, OtherTarget_, OtherColor_, OtherSupportCompare_, OtherSourceCtnr_, OtherTargetCtnr_> type;
1395: };
1396: // Re-export some typedefs expected by CoSifter
1397: typedef typename traits::source_type source_type;
1398: typedef typename traits::target_type target_type;
1399: typedef typename traits::color_type color_type;
1400: typedef typename traits::arrow_type Arrow_;
1401: typedef typename traits::coneSequence coneSequence;
1402: typedef typename traits::supportSequence supportSequence;
1403: typedef typename traits::baseSequence baseSequence;
1404: typedef typename traits::capSequence capSequence;
1405: // Basic interface
1406: Sifter(MPI_Comm comm = PETSC_COMM_SELF, const int& debug = 0) :
1407: ASifter<Source_, Target_, Color_, SifterDef::uniColor, SupportCompare_, SourceCtnr_, TargetCtnr_>(comm, debug) {};
1409: const typename traits::color_type&
1410: getColor(const typename traits::source_type& s, const typename traits::target_type& t, bool fail = true) {
1411: typedef typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type,typename traits::arrowInd>::type index_type;
1413: const index_type& _index = ::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set);
1414: #if 0
1415: ::boost::tuple<typename traits::source_type, typename traits::target_type> key = ::boost::make_tuple(s, t);
1416: typename index_type::iterator begin = _index.lower_bound(key);
1417: if(begin != _index.upper_bound(key)) {
1418: return begin->color;
1419: }
1420: #else
1421: const typename index_type::iterator begin = _index.find(::boost::make_tuple(s, t));
1422: if (begin != _index.end()) {
1423: return begin->color;
1424: }
1425: #endif
1426: // typename traits::arrowSequence arr(::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set), s, t);
1427: // if(arr.begin() != arr.end()) {
1428: // return arr.begin().color();
1429: // }
1430: if (fail) {
1431: ostringstream o;
1432: o << "Arrow " << s << " --> " << t << " not present";
1433: throw ALE::Exception(o.str().c_str());
1434: } else {
1435: static typename traits::color_type c;
1436: return c;
1437: }
1438: };
1440: template<typename ColorChanger>
1441: void modifyColor(const typename traits::source_type& s, const typename traits::target_type& t, const ColorChanger& changeColor) {
1442: typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type, typename traits::arrowInd>::type& index =
1443: ::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set);
1444: typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type, typename traits::arrowInd>::type::iterator i =
1445: index.find(::boost::make_tuple(s,t));
1446: if (i != index.end()) {
1447: index.modify(i, changeColor);
1448: } else {
1449: typename traits::arrow_type a(s, t, typename traits::color_type());
1450: changeColor(a);
1451: this->addArrow(a);
1452: }
1453: }
1455: struct ColorSetter {
1456: ColorSetter(const typename traits::color_type& color) : _color(color) {};
1457: void operator()(typename traits::arrow_type& p) const {
1458: p.color = _color;
1459: }
1460: private:
1461: const typename traits::color_type& _color;
1462: };
1464: void setColor(const typename traits::source_type& s, const typename traits::target_type& t, const typename traits::color_type& color) {
1465: ColorSetter colorSetter(color);
1466: typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type, typename traits::arrowInd>::type& index =
1467: ::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set);
1468: typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type, typename traits::arrowInd>::type::iterator i =
1469: index.find(::boost::make_tuple(s,t));
1470: if (i != index.end()) {
1471: index.modify(i, colorSetter);
1472: } else {
1473: typename traits::arrow_type a(s, t, color);
1474: this->addArrow(a);
1475: }
1476: };
1477: template<typename Labeling, typename AnotherSifter>
1478: void relabel(Labeling& relabeling, AnotherSifter& newLabel) {
1479: typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type, typename traits::arrowInd>::type& aInd = ::boost::multi_index::get<typename traits::arrowInd>(this->_arrows.set);
1481: for(typename ::boost::multi_index::index<typename traits::arrow_container_type::set_type, typename traits::arrowInd>::type::iterator a_iter = aInd.begin(); a_iter != aInd.end(); ++a_iter) {
1482: const typename traits::source_type newSource = relabeling.restrictPoint(a_iter->source)[0];
1483: const typename traits::target_type newTarget = relabeling.restrictPoint(a_iter->target)[0];
1485: newLabel.addArrow(newSource, newTarget);
1486: }
1487: }
1488: };// class Sifter
1490: class SifterSerializer {
1491: public:
1492: template<typename Sifter>
1493: static void writeSifter(std::ofstream& fs, Sifter& sifter) {
1494: typename Sifter::traits::arrow_container_type::set_type::size_type numArrows;
1496: if (sifter.commRank() == 0) {
1497: // Write local
1498: fs << sifter._arrows.set.size() << std::endl;
1499: for(typename Sifter::traits::arrow_container_type::set_type::iterator ai = sifter._arrows.set.begin(); ai != sifter._arrows.set.end(); ai++) {
1500: fs << ai->source << " " << ai->target << " " << ai->color << std::endl;
1501: }
1502: // Receive and write remote
1503: for(int p = 1; p < sifter.commSize(); ++p) {
1504: PetscInt size;
1505: PetscInt *arrows;
1506: MPI_Status status;
1509: MPI_Recv(&size, 1, MPIU_INT, p, 1, sifter.comm(), &status);CHKERRXX(ierr);
1510: numArrows = size;
1511: fs << numArrows << std::endl;
1512: PetscMalloc(size*3 * sizeof(PetscInt), &arrows);CHKERRXX(ierr);
1513: MPI_Recv(arrows, size*3, MPIU_INT, p, 1, sifter.comm(), &status);CHKERRXX(ierr);
1514: for(PetscInt a = 0; a < size; ++a) {
1515: typename Sifter::traits::arrow_type::source_type source = arrows[a*3+0];
1516: typename Sifter::traits::arrow_type::target_type target = arrows[a*3+1];
1517: typename Sifter::traits::arrow_type::color_type color = arrows[a*3+2];
1519: fs << source << " " << target << " " << color << std::endl;
1520: }
1521: PetscFree(arrows);CHKERRXX(ierr);
1522: }
1523: } else {
1524: // Send remote
1525: PetscInt size = sifter._arrows.set.size();
1526: PetscInt a = 0;
1527: PetscInt *arrows;
1530: MPI_Send(&size, 1, MPIU_INT, 0, 1, sifter.comm());CHKERRXX(ierr);
1531: // There is no nice way to make a generic MPI type here. Really sucky
1532: PetscMalloc(size*3 * sizeof(PetscInt), &arrows);CHKERRXX(ierr);
1533: for(typename Sifter::traits::arrow_container_type::set_type::iterator ai = sifter._arrows.set.begin(); ai != sifter._arrows.set.end(); ai++, ++a) {
1534: arrows[a*3+0] = ai->source;
1535: arrows[a*3+1] = ai->target;
1536: arrows[a*3+2] = ai->color;
1537: }
1538: MPI_Send(arrows, size*3, MPIU_INT, 0, 1, sifter.comm());CHKERRXX(ierr);
1539: PetscFree(arrows);CHKERRXX(ierr);
1540: }
1541: }
1542: template<typename Sifter>
1543: static void loadSifter(std::ifstream& fs, Sifter& sifter) {
1544: typedef typename Sifter::traits::arrow_container_type::set_type::size_type size_type;
1545: if (sifter.commRank() == 0) {
1546: // Load local
1547: size_type numArrows;
1549: fs >> numArrows;
1550: for(size_type a = 0; a < numArrows; ++a) {
1551: typename Sifter::traits::arrow_type::source_type source;
1552: typename Sifter::traits::arrow_type::target_type target;
1553: typename Sifter::traits::arrow_type::color_type color;
1555: fs >> source;
1556: fs >> target;
1557: fs >> color;
1558: sifter.addArrow(typename Sifter::traits::arrow_type(source, target, color));
1559: }
1560: // Load and send remote
1561: for(int p = 1; p < sifter.commSize(); ++p) {
1562: PetscInt size;
1563: PetscInt *arrows;
1566: fs >> numArrows;
1567: size = numArrows;
1568: MPI_Send(&size, 1, MPIU_INT, p, 1, sifter.comm());CHKERRXX(ierr);
1569: PetscMalloc(size*3 * sizeof(PetscInt), &arrows);CHKERRXX(ierr);
1570: for(PetscInt a = 0; a < size; ++a) {
1571: typename Sifter::traits::arrow_type::source_type source;
1572: typename Sifter::traits::arrow_type::target_type target;
1573: typename Sifter::traits::arrow_type::color_type color;
1575: fs >> source;
1576: fs >> target;
1577: fs >> color;
1578: arrows[a*3+0] = source;
1579: arrows[a*3+1] = target;
1580: arrows[a*3+2] = color;
1581: }
1582: MPI_Send(arrows, size*3, MPIU_INT, p, 1, sifter.comm());CHKERRXX(ierr);
1583: PetscFree(arrows);CHKERRXX(ierr);
1584: }
1585: } else {
1586: // Load remote
1587: PetscInt size;
1588: PetscInt *arrows;
1589: MPI_Status status;
1592: MPI_Recv(&size, 1, MPIU_INT, 0, 1, sifter.comm(), &status);CHKERRXX(ierr);
1593: PetscMalloc(size*3 * sizeof(PetscInt), &arrows);CHKERRXX(ierr);
1594: MPI_Recv(arrows, size*3, MPIU_INT, 0, 1, sifter.comm(), &status);CHKERRXX(ierr);
1595: for(PetscInt a = 0; a < size; ++a) {
1596: sifter.addArrow(typename Sifter::traits::arrow_type(arrows[a*3+0], arrows[a*3+1], arrows[a*3+2]));
1597: }
1598: PetscFree(arrows);CHKERRXX(ierr);
1599: }
1600: }
1601: };
1602: } // namespace ALE
1604: #endif