Branch data Line data Source code
1 : : #include "moab/ScdInterface.hpp"
2 : : #include "moab/Core.hpp"
3 : : #include "SequenceManager.hpp"
4 : : #include "EntitySequence.hpp"
5 : : #include "StructuredElementSeq.hpp"
6 : : #include "VertexSequence.hpp"
7 : : #include "ScdVertexData.hpp"
8 : : #include "MBTagConventions.hpp"
9 : : #ifdef MOAB_HAVE_MPI
10 : : #include "moab/ParallelComm.hpp"
11 : : #include "moab/TupleList.hpp"
12 : : #include "moab/gs.hpp"
13 : : #endif
14 : : #include "assert.h"
15 : : #include <iostream>
16 : : #include <functional>
17 : :
18 : : #define ERRORR( rval, str ) \
19 : : { \
20 : : if( MB_SUCCESS != rval ) \
21 : : { \
22 : : std::cerr << str; \
23 : : return rval; \
24 : : } \
25 : : }
26 : :
27 : : namespace moab
28 : : {
29 : :
30 : : const char* ScdParData::PartitionMethodNames[] = { "alljorkori", "alljkbal", "sqij", "sqjk",
31 : : "sqijk", "trivial", "rcbzoltan", "nopart" };
32 : :
33 : 29 : ScdInterface::ScdInterface( Interface* imp, bool boxes )
34 : : : mbImpl( imp ), searchedBoxes( false ), boxPeriodicTag( 0 ), boxDimsTag( 0 ), globalBoxDimsTag( 0 ),
35 : 29 : partMethodTag( 0 ), boxSetTag( 0 )
36 : : {
37 [ - + ][ # # ]: 29 : if( boxes ) find_boxes( scdBoxes );
38 : 29 : }
39 : :
40 : : // Destructor
41 : 58 : ScdInterface::~ScdInterface()
42 : : {
43 : 29 : std::vector< ScdBox* > tmp_boxes;
44 : 29 : tmp_boxes.swap( scdBoxes );
45 : :
46 [ + + ]: 75 : for( std::vector< ScdBox* >::iterator rit = tmp_boxes.begin(); rit != tmp_boxes.end(); ++rit )
47 [ + - ]: 46 : delete *rit;
48 : :
49 [ + + ]: 29 : if( box_set_tag( false ) ) mbImpl->tag_delete( box_set_tag() );
50 : 29 : }
51 : :
52 : 14761 : Interface* ScdInterface::impl() const
53 : : {
54 : 14761 : return mbImpl;
55 : : }
56 : :
57 : 2 : ErrorCode ScdInterface::find_boxes( std::vector< ScdBox* >& scd_boxes )
58 : : {
59 [ + - ]: 2 : Range tmp_boxes;
60 [ + - ]: 2 : ErrorCode rval = find_boxes( tmp_boxes );
61 [ - + ]: 2 : if( MB_SUCCESS != rval ) return rval;
62 : :
63 [ + - ][ + - ]: 4 : for( Range::iterator rit = tmp_boxes.begin(); rit != tmp_boxes.end(); ++rit )
[ + - ][ + - ]
[ + + ]
64 : : {
65 [ + - ][ + - ]: 2 : ScdBox* tmp_box = get_scd_box( *rit );
66 [ + - ]: 2 : if( tmp_box )
67 [ + - ]: 2 : scd_boxes.push_back( tmp_box );
68 : : else
69 : 0 : rval = MB_FAILURE;
70 : : }
71 : :
72 : 2 : return rval;
73 : : }
74 : :
75 : 19 : ErrorCode ScdInterface::find_boxes( Range& scd_boxes )
76 : : {
77 : 19 : ErrorCode rval = MB_SUCCESS;
78 [ + - ]: 19 : box_dims_tag();
79 [ + - ]: 19 : Range boxes;
80 [ + + ]: 19 : if( !searchedBoxes )
81 : : {
82 [ + - ]: 18 : rval = mbImpl->get_entities_by_type_and_tag( 0, MBENTITYSET, &boxDimsTag, NULL, 1, boxes, Interface::UNION );
83 : 18 : searchedBoxes = true;
84 [ + - ][ + + ]: 18 : if( !boxes.empty() )
85 : : {
86 [ + - ][ + - ]: 1 : scdBoxes.resize( boxes.size() );
87 [ + - ][ + - ]: 1 : rval = mbImpl->tag_get_data( boxSetTag, boxes, &scdBoxes[0] );
88 : 1 : ScdBox* dum = NULL;
89 : : // std::remove_if(scdBoxes.begin(), scdBoxes.end(),
90 : : // std::bind2nd(std::equal_to<ScdBox*>(), dum) ) ;
91 : : std::remove_if( scdBoxes.begin(), scdBoxes.end(),
92 [ + - ][ + - ]: 18 : std::bind( std::equal_to< ScdBox* >(), std::placeholders::_1, dum ) );
93 : : // https://stackoverflow.com/questions/32739018/a-replacement-for-stdbind2nd
94 : : }
95 : : }
96 : :
97 [ + - ][ + - ]: 21 : for( std::vector< ScdBox* >::iterator vit = scdBoxes.begin(); vit != scdBoxes.end(); ++vit )
[ + + ]
98 [ + - ][ + - ]: 2 : scd_boxes.insert( ( *vit )->box_set() );
[ + - ]
99 : :
100 : 19 : return rval;
101 : : }
102 : :
103 : 2 : ScdBox* ScdInterface::get_scd_box( EntityHandle eh )
104 : : {
105 : 2 : ScdBox* scd_box = NULL;
106 [ + - ][ - + ]: 2 : if( !box_set_tag( false ) ) return scd_box;
107 : :
108 [ + - ][ + - ]: 2 : mbImpl->tag_get_data( box_set_tag(), &eh, 1, &scd_box );
109 : 2 : return scd_box;
110 : : }
111 : :
112 : 13 : ErrorCode ScdInterface::construct_box( HomCoord low, HomCoord high, const double* const coords, unsigned int num_coords,
113 : : ScdBox*& new_box, int* const lperiodic, ScdParData* par_data, bool assign_gids,
114 : : int tag_shared_ents )
115 : : {
116 : : // create a rectangular structured mesh block
117 : : ErrorCode rval;
118 : :
119 : 13 : int tmp_lper[3] = { 0, 0, 0 };
120 [ + + ][ + - ]: 13 : if( lperiodic ) std::copy( lperiodic, lperiodic + 3, tmp_lper );
121 : :
122 : : #ifndef MOAB_HAVE_MPI
123 : : if( -1 != tag_shared_ents ) ERRORR( MB_FAILURE, "Parallel capability requested but MOAB not compiled parallel." );
124 : : if( -1 == tag_shared_ents && !assign_gids ) assign_gids = true; // need to assign gids in order to tag shared verts
125 : : #else
126 [ - + ][ # # ]: 13 : if( par_data && low == high && ScdParData::NOPART != par_data->partMethod )
[ # # ][ # # ]
[ - + ]
127 : : {
128 : : // requesting creation of parallel mesh, so need to compute partition
129 [ # # ]: 0 : if( !par_data->pComm )
130 : : {
131 : : // this is a really boneheaded way to have to create a PC
132 [ # # ]: 0 : par_data->pComm = ParallelComm::get_pcomm( mbImpl, 0 );
133 [ # # ][ # # ]: 0 : if( NULL == par_data->pComm ) par_data->pComm = new ParallelComm( mbImpl, MPI_COMM_WORLD );
[ # # ]
134 : : }
135 : : int ldims[6];
136 [ # # ][ # # ]: 0 : rval = compute_partition( par_data->pComm->size(), par_data->pComm->rank(), *par_data, ldims, tmp_lper,
137 [ # # ][ # # ]: 0 : par_data->pDims );ERRORR( rval, "Error returned from compute_partition." );
[ # # ]
138 [ # # ]: 0 : low.set( ldims[0], ldims[1], ldims[2] );
139 [ # # ]: 0 : high.set( ldims[3], ldims[4], ldims[5] );
140 [ # # ][ # # ]: 0 : if( par_data->pComm->get_debug_verbosity() > 0 )
141 : : {
142 [ # # ][ # # ]: 0 : std::cout << "Proc " << par_data->pComm->rank() << ": " << *par_data;
[ # # ][ # # ]
[ # # ]
143 [ # # ][ # # ]: 0 : std::cout << "Proc " << par_data->pComm->rank() << " local dims: " << low << "-" << high << std::endl;
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
144 : : }
145 : : }
146 : : #endif
147 : :
148 [ + - ][ + - ]: 13 : HomCoord tmp_size = high - low + HomCoord( 1, 1, 1, 0 );
[ + - ]
149 [ + - ][ + - ]: 52 : if( ( tmp_size[1] && num_coords && (int)num_coords < tmp_size[0] ) ||
[ - + ][ # # ]
[ # # ][ + - ]
[ - + ][ + - ]
[ - + # #
# # ]
150 [ + - ][ - + ]: 39 : ( tmp_size[2] && num_coords && (int)num_coords < tmp_size[0] * tmp_size[1] ) )
[ # # ][ # # ]
[ # # ][ - + ]
[ - + ][ + - ]
[ # # # #
# # ]
151 : 0 : return MB_FAILURE;
152 : :
153 [ + - ][ - + ]: 13 : rval = create_scd_sequence( low, high, MBVERTEX, 0, new_box );ERRORR( rval, "Trouble creating scd vertex sequence." );
[ # # ]
154 : :
155 : : // set the vertex coordinates
156 : : double *xc, *yc, *zc;
157 [ + - ][ - + ]: 13 : rval = new_box->get_coordinate_arrays( xc, yc, zc );ERRORR( rval, "Couldn't get vertex coordinate arrays." );
[ # # ]
158 : :
159 [ - + ][ # # ]: 13 : if( coords && num_coords )
160 : : {
161 : 0 : unsigned int i = 0;
162 [ # # ][ # # ]: 0 : for( int kl = low[2]; kl <= high[2]; kl++ )
[ # # ]
163 : : {
164 [ # # ][ # # ]: 0 : for( int jl = low[1]; jl <= high[1]; jl++ )
[ # # ]
165 : : {
166 [ # # ][ # # ]: 0 : for( int il = low[0]; il <= high[0]; il++ )
[ # # ]
167 : : {
168 : 0 : xc[i] = coords[3 * i];
169 [ # # ][ # # ]: 0 : if( new_box->box_size()[1] ) yc[i] = coords[3 * i + 1];
[ # # ]
170 [ # # ][ # # ]: 0 : if( new_box->box_size()[2] ) zc[i] = coords[3 * i + 2];
[ # # ]
171 : 0 : i++;
172 : : }
173 : : }
174 : 0 : }
175 : : }
176 : : else
177 : : {
178 : 13 : unsigned int i = 0;
179 [ + - ][ + - ]: 75 : for( int kl = low[2]; kl <= high[2]; kl++ )
[ + + ]
180 : : {
181 [ + - ][ + - ]: 772 : for( int jl = low[1]; jl <= high[1]; jl++ )
[ + + ]
182 : : {
183 [ + - ][ + - ]: 15412 : for( int il = low[0]; il <= high[0]; il++ )
[ + + ]
184 : : {
185 : 14702 : xc[i] = (double)il;
186 [ + - ][ + - ]: 14702 : if( new_box->box_size()[1] )
[ + - ]
187 : 14702 : yc[i] = (double)jl;
188 : : else
189 : 0 : yc[i] = 0.0;
190 [ + - ][ + - ]: 14702 : if( new_box->box_size()[2] )
[ + - ]
191 : 14702 : zc[i] = (double)kl;
192 : : else
193 : 0 : zc[i] = 0.0;
194 : 14702 : i++;
195 : : }
196 : : }
197 : : }
198 : : }
199 : :
200 : : // create element sequence
201 [ - + ]: 13 : Core* mbcore = dynamic_cast< Core* >( mbImpl );
202 [ + - ]: 13 : SequenceManager* seq_mgr = mbcore->sequence_manager();
203 : :
204 : : EntitySequence* tmp_seq;
205 : : EntityHandle start_ent;
206 : :
207 : : // construct the sequence
208 : 13 : EntityType this_tp = MBHEX;
209 [ + - ][ + + ]: 13 : if( 1 >= tmp_size[2] ) this_tp = MBQUAD;
210 [ + - ][ + + ]: 13 : if( 1 >= tmp_size[2] && 1 >= tmp_size[1] ) this_tp = MBEDGE;
[ + - ][ - + ]
[ + + ][ + - ]
[ - + ]
[ # # # # ]
211 [ + - ][ - + ]: 13 : rval = seq_mgr->create_scd_sequence( low, high, this_tp, 0, start_ent, tmp_seq, tmp_lper );ERRORR( rval, "Trouble creating scd element sequence." );
[ # # ]
212 : :
213 [ + - ]: 13 : new_box->elem_seq( tmp_seq );
214 [ + - ]: 13 : new_box->start_element( start_ent );
215 : :
216 : : // add vertex seq to element seq, forward orientation, unity transform
217 : : rval = new_box->add_vbox( new_box,
218 : : // p1: imin,jmin
219 : : low, low,
220 : : // p2: imax,jmin
221 [ + - ][ + - ]: 13 : low + HomCoord( 1, 0, 0 ), low + HomCoord( 1, 0, 0 ),
[ + - ][ + - ]
222 : : // p3: imin,jmax
223 [ + - ][ + - ]: 26 : low + HomCoord( 0, 1, 0 ), low + HomCoord( 0, 1, 0 ) );ERRORR( rval, "Error constructing structured element sequence." );
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ - + ]
[ # # ]
224 : :
225 : : // add the new hexes to the scd box set; vertices were added in call to create_scd_sequence
226 [ + - ][ + - ]: 13 : Range tmp_range( new_box->start_element(), new_box->start_element() + new_box->num_elements() - 1 );
[ + - ][ + - ]
227 [ + - ][ + - ]: 13 : rval = mbImpl->add_entities( new_box->box_set(), tmp_range );ERRORR( rval, "Couldn't add new hexes to box set." );
[ - + ][ # # ]
228 : :
229 [ - + ][ # # ]: 13 : if( par_data ) new_box->par_data( *par_data );
230 : :
231 [ + + ]: 13 : if( assign_gids )
232 : : {
233 [ + - ][ - + ]: 2 : rval = assign_global_ids( new_box );ERRORR( rval, "Trouble assigning global ids" );
[ # # ]
234 : : }
235 : :
236 : : #ifdef MOAB_HAVE_MPI
237 [ - + ][ # # ]: 13 : if( par_data && -1 != tag_shared_ents ) { rval = tag_shared_vertices( par_data->pComm, new_box ); }
[ # # ]
238 : : #endif
239 : :
240 : 13 : return MB_SUCCESS;
241 : : }
242 : :
243 : 2 : ErrorCode ScdInterface::assign_global_ids( ScdBox* box )
244 : : {
245 : : // Get a ptr to global id memory
246 : : void* data;
247 : 2 : int count = 0;
248 [ + - ]: 2 : Tag gid_tag = mbImpl->globalId_tag();
249 [ + - ][ + - ]: 2 : Range tmp_range( box->start_vertex(), box->start_vertex() + box->num_vertices() );
[ + - ][ + - ]
250 [ + - ][ + - ]: 2 : ErrorCode rval = mbImpl->tag_iterate( gid_tag, tmp_range.begin(), tmp_range.end(), count, data );ERRORR( rval, "Failed to get tag iterator." );
[ + - ][ - + ]
[ # # ]
251 [ + - ][ - + ]: 2 : assert( count == box->num_vertices() );
252 : 2 : int* gid_data = (int*)data;
253 [ + - ][ + - ]: 2 : int di = box->par_data().gDims[3] - box->par_data().gDims[0] + 1;
254 [ + - ][ + - ]: 2 : int dj = box->par_data().gDims[4] - box->par_data().gDims[1] + 1;
255 : :
256 [ + - ][ + - ]: 4 : for( int kl = box->box_dims()[2]; kl <= box->box_dims()[5]; kl++ )
[ + + ]
257 : : {
258 [ + - ][ + - ]: 132 : for( int jl = box->box_dims()[1]; jl <= box->box_dims()[4]; jl++ )
[ + + ]
259 : : {
260 [ + - ][ + - ]: 8580 : for( int il = box->box_dims()[0]; il <= box->box_dims()[3]; il++ )
[ + + ]
261 : : {
262 : : int itmp =
263 [ + - ][ + - ]: 16900 : ( !box->locally_periodic()[0] && box->par_data().gPeriodic[0] && il == box->par_data().gDims[3]
[ - + ][ # # ]
[ # # ]
264 [ # # ]: 0 : ? box->par_data().gDims[0]
265 [ + - ]: 16900 : : il );
266 [ + - ]: 8450 : *gid_data = ( -1 != kl ? kl * di * dj : 0 ) + jl * di + itmp + 1;
267 : 8450 : gid_data++;
268 : : }
269 : : }
270 : : }
271 : :
272 : 2 : return MB_SUCCESS;
273 : : }
274 : :
275 : 46 : ErrorCode ScdInterface::create_scd_sequence( const HomCoord& low, const HomCoord& high, EntityType tp, int starting_id,
276 : : ScdBox*& new_box, int* is_periodic )
277 : : {
278 [ + - ][ + - ]: 46 : HomCoord tmp_size = high - low + HomCoord( 1, 1, 1, 0 );
[ + - ]
279 [ + + ][ + - ]: 97 : if( ( tp == MBHEX && 1 >= tmp_size[2] ) || ( tp == MBQUAD && 1 >= tmp_size[1] ) ||
[ + - ][ + + ]
[ + - ][ + - ]
[ + + ][ + + ]
[ + + ][ - +
# # # # ]
280 [ + - ][ - + ]: 51 : ( tp == MBEDGE && 1 >= tmp_size[0] ) )
[ + + ][ # # ]
281 : 0 : return MB_TYPE_OUT_OF_RANGE;
282 : :
283 [ - + ]: 46 : Core* mbcore = dynamic_cast< Core* >( mbImpl );
284 [ - + ]: 46 : assert( mbcore != NULL );
285 [ + - ]: 46 : SequenceManager* seq_mgr = mbcore->sequence_manager();
286 : :
287 : : EntitySequence* tmp_seq;
288 : : EntityHandle start_ent, scd_set;
289 : :
290 : : // construct the sequence
291 [ + - ]: 46 : ErrorCode rval = seq_mgr->create_scd_sequence( low, high, tp, starting_id, start_ent, tmp_seq, is_periodic );
292 [ - + ]: 46 : if( MB_SUCCESS != rval ) return rval;
293 : :
294 : : // create the set for this rectangle
295 [ + - ]: 46 : rval = create_box_set( low, high, scd_set );
296 [ - + ]: 46 : if( MB_SUCCESS != rval ) return rval;
297 : :
298 : : // make the ScdBox
299 [ + - ][ + - ]: 46 : new_box = new ScdBox( this, scd_set, tmp_seq );
300 [ - + ]: 46 : if( !new_box ) return MB_FAILURE;
301 : :
302 : : // set the start vertex/element
303 [ + - ]: 46 : Range new_range;
304 [ + + ][ + - ]: 46 : if( MBVERTEX == tp ) { new_range.insert( start_ent, start_ent + new_box->num_vertices() - 1 ); }
[ + - ]
305 : : else
306 : : {
307 [ + - ][ + - ]: 15 : new_range.insert( start_ent, start_ent + new_box->num_elements() - 1 );
308 : : }
309 : :
310 : : // put the entities in the box set
311 [ + - ]: 46 : rval = mbImpl->add_entities( scd_set, new_range );
312 [ - + ]: 46 : if( MB_SUCCESS != rval ) return rval;
313 : :
314 : : // tag the set with the box
315 [ + - ][ + - ]: 46 : rval = mbImpl->tag_set_data( box_set_tag(), &scd_set, 1, &new_box );
316 [ - + ]: 46 : if( MB_SUCCESS != rval ) return rval;
317 : :
318 : 46 : return MB_SUCCESS;
319 : : }
320 : :
321 : 46 : ErrorCode ScdInterface::create_box_set( const HomCoord& low, const HomCoord& high, EntityHandle& scd_set,
322 : : int* is_periodic )
323 : : {
324 : : // create the set and put the entities in it
325 [ + - ]: 46 : ErrorCode rval = mbImpl->create_meshset( MESHSET_SET, scd_set );
326 [ - + ]: 46 : if( MB_SUCCESS != rval ) return rval;
327 : :
328 : : // tag the set with parameter extents
329 : : int boxdims[6];
330 [ + + ]: 184 : for( int i = 0; i < 3; i++ )
331 [ + - ]: 138 : boxdims[i] = low[i];
332 [ + + ]: 184 : for( int i = 0; i < 3; i++ )
333 [ + - ]: 138 : boxdims[3 + i] = high[i];
334 [ + - ][ + - ]: 46 : rval = mbImpl->tag_set_data( box_dims_tag(), &scd_set, 1, boxdims );
335 [ - + ]: 46 : if( MB_SUCCESS != rval ) return rval;
336 : :
337 [ - + ]: 46 : if( is_periodic )
338 : : {
339 [ # # ][ # # ]: 0 : rval = mbImpl->tag_set_data( box_periodic_tag(), &scd_set, 1, is_periodic );
340 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
341 : : }
342 : :
343 : 46 : return rval;
344 : : }
345 : :
346 : 0 : Tag ScdInterface::box_periodic_tag( bool create_if_missing )
347 : : {
348 : : // Reset boxPeriodicTag in case it has been deleted (e.g. by Core::clean_up_failed_read)
349 [ # # ]: 0 : if( boxPeriodicTag )
350 : : {
351 [ # # ]: 0 : std::string tag_name;
352 [ # # ][ # # ]: 0 : if( MB_TAG_NOT_FOUND == mbImpl->tag_get_name( boxPeriodicTag, tag_name ) ) boxPeriodicTag = NULL;
353 : : }
354 : :
355 [ # # ][ # # ]: 0 : if( boxPeriodicTag || !create_if_missing ) return boxPeriodicTag;
356 : :
357 : : ErrorCode rval =
358 : 0 : mbImpl->tag_get_handle( "BOX_PERIODIC", 3, MB_TYPE_INTEGER, boxPeriodicTag, MB_TAG_SPARSE | MB_TAG_CREAT );
359 [ # # ]: 0 : if( MB_SUCCESS != rval ) return 0;
360 : 0 : return boxPeriodicTag;
361 : : }
362 : :
363 : 65 : Tag ScdInterface::box_dims_tag( bool create_if_missing )
364 : : {
365 : : // Reset boxDimsTag in case it has been deleted (e.g. by clean_up_failed_read)
366 [ + + ]: 65 : if( boxDimsTag )
367 : : {
368 [ + - ]: 38 : std::string tag_name;
369 [ + - ][ - + ]: 38 : if( MB_TAG_NOT_FOUND == mbImpl->tag_get_name( boxDimsTag, tag_name ) ) boxDimsTag = NULL;
370 : : }
371 : :
372 [ + + ][ - + ]: 65 : if( boxDimsTag || !create_if_missing ) return boxDimsTag;
373 : :
374 : 27 : ErrorCode rval = mbImpl->tag_get_handle( "BOX_DIMS", 6, MB_TYPE_INTEGER, boxDimsTag, MB_TAG_SPARSE | MB_TAG_CREAT );
375 [ - + ]: 27 : if( MB_SUCCESS != rval ) return 0;
376 : 65 : return boxDimsTag;
377 : : }
378 : :
379 : 0 : Tag ScdInterface::global_box_dims_tag( bool create_if_missing )
380 : : {
381 : : // Reset globalBoxDimsTag in case it has been deleted (e.g. by Core::clean_up_failed_read)
382 [ # # ]: 0 : if( globalBoxDimsTag )
383 : : {
384 [ # # ]: 0 : std::string tag_name;
385 [ # # ][ # # ]: 0 : if( MB_TAG_NOT_FOUND == mbImpl->tag_get_name( globalBoxDimsTag, tag_name ) ) globalBoxDimsTag = NULL;
386 : : }
387 : :
388 [ # # ][ # # ]: 0 : if( globalBoxDimsTag || !create_if_missing ) return globalBoxDimsTag;
389 : :
390 : : ErrorCode rval =
391 : 0 : mbImpl->tag_get_handle( "GLOBAL_BOX_DIMS", 6, MB_TYPE_INTEGER, globalBoxDimsTag, MB_TAG_SPARSE | MB_TAG_CREAT );
392 [ # # ]: 0 : if( MB_SUCCESS != rval ) return 0;
393 : 0 : return globalBoxDimsTag;
394 : : }
395 : :
396 : 0 : Tag ScdInterface::part_method_tag( bool create_if_missing )
397 : : {
398 : : // Reset partMethodTag in case it has been deleted (e.g. by Core::clean_up_failed_read)
399 [ # # ]: 0 : if( partMethodTag )
400 : : {
401 [ # # ]: 0 : std::string tag_name;
402 [ # # ][ # # ]: 0 : if( MB_TAG_NOT_FOUND == mbImpl->tag_get_name( partMethodTag, tag_name ) ) partMethodTag = NULL;
403 : : }
404 : :
405 [ # # ][ # # ]: 0 : if( partMethodTag || !create_if_missing ) return partMethodTag;
406 : :
407 : : ErrorCode rval =
408 : 0 : mbImpl->tag_get_handle( "PARTITION_METHOD", 1, MB_TYPE_INTEGER, partMethodTag, MB_TAG_SPARSE | MB_TAG_CREAT );
409 [ # # ]: 0 : if( MB_SUCCESS != rval ) return 0;
410 : 0 : return partMethodTag;
411 : : }
412 : :
413 : 135 : Tag ScdInterface::box_set_tag( bool create_if_missing )
414 : : {
415 : : // Reset boxSetTag in case it has been deleted (e.g. by Core::clean_up_failed_read)
416 [ + + ]: 135 : if( boxSetTag )
417 : : {
418 [ + - ]: 106 : std::string tag_name;
419 [ + - ][ - + ]: 106 : if( MB_TAG_NOT_FOUND == mbImpl->tag_get_name( boxSetTag, tag_name ) ) boxSetTag = NULL;
420 : : }
421 : :
422 [ + + ][ + + ]: 135 : if( boxSetTag || !create_if_missing ) return boxSetTag;
423 : :
424 : : ErrorCode rval = mbImpl->tag_get_handle( "__BOX_SET", sizeof( ScdBox* ), MB_TYPE_OPAQUE, boxSetTag,
425 : 10 : MB_TAG_SPARSE | MB_TAG_CREAT );
426 [ - + ]: 10 : if( MB_SUCCESS != rval ) return 0;
427 : 135 : return boxSetTag;
428 : : }
429 : :
430 : : //! Remove the box from the list on ScdInterface
431 : 46 : ErrorCode ScdInterface::remove_box( ScdBox* box )
432 : : {
433 [ + - ]: 46 : std::vector< ScdBox* >::iterator vit = std::find( scdBoxes.begin(), scdBoxes.end(), box );
434 [ + - ][ - + ]: 46 : if( vit != scdBoxes.end() )
435 : : {
436 [ # # ]: 0 : scdBoxes.erase( vit );
437 : 0 : return MB_SUCCESS;
438 : : }
439 : : else
440 : 46 : return MB_FAILURE;
441 : : }
442 : :
443 : : //! Add the box to the list on ScdInterface
444 : 46 : ErrorCode ScdInterface::add_box( ScdBox* box )
445 : : {
446 : 46 : scdBoxes.push_back( box );
447 : 46 : return MB_SUCCESS;
448 : : }
449 : :
450 : 0 : ErrorCode ScdInterface::get_boxes( std::vector< ScdBox* >& boxes )
451 : : {
452 : 0 : std::copy( scdBoxes.begin(), scdBoxes.end(), std::back_inserter( boxes ) );
453 : 0 : return MB_SUCCESS;
454 : : }
455 : :
456 : 46 : ScdBox::ScdBox( ScdInterface* impl, EntityHandle bset, EntitySequence* seq1, EntitySequence* seq2 )
457 : 46 : : scImpl( impl ), boxSet( bset ), vertDat( NULL ), elemSeq( NULL ), startVertex( 0 ), startElem( 0 )
458 : : {
459 [ + + ]: 322 : for( int i = 0; i < 6; i++ )
460 : 276 : boxDims[i] = 0;
461 [ + + ]: 184 : for( int i = 0; i < 3; i++ )
462 : 138 : locallyPeriodic[i] = false;
463 [ - + ]: 46 : VertexSequence* vseq = dynamic_cast< VertexSequence* >( seq1 );
464 [ + + ][ - + ]: 46 : if( vseq ) vertDat = dynamic_cast< ScdVertexData* >( vseq->data() );
465 [ + + ]: 46 : if( vertDat )
466 : : {
467 : : // retrieve the parametric space
468 [ + + ]: 124 : for( int i = 0; i < 3; i++ )
469 : : {
470 [ + - ][ + - ]: 93 : boxDims[i] = vertDat->min_params()[i];
471 [ + - ][ + - ]: 93 : boxDims[3 + i] = vertDat->max_params()[i];
472 : : }
473 : 31 : startVertex = vertDat->start_handle();
474 : : }
475 [ + - ]: 15 : else if( impl->boxDimsTag )
476 : : {
477 : : // look for parametric space info on set
478 : 15 : ErrorCode rval = impl->mbImpl->tag_get_data( impl->boxDimsTag, &bset, 1, boxDims );
479 [ + - ]: 15 : if( MB_SUCCESS == rval )
480 : : {
481 [ + - ]: 15 : Range verts;
482 [ + - ]: 15 : impl->mbImpl->get_entities_by_dimension( bset, 0, verts );
483 [ + - ][ - + ]: 15 : if( !verts.empty() ) startVertex = *verts.begin();
[ # # ][ # # ]
484 : : }
485 : : }
486 : :
487 [ + - ]: 46 : elemSeq = dynamic_cast< StructuredElementSeq* >( seq2 );
488 [ + - ][ - + ]: 46 : if( !elemSeq ) elemSeq = dynamic_cast< StructuredElementSeq* >( seq1 );
489 : :
490 [ + + ]: 46 : if( elemSeq )
491 : : {
492 [ - + ]: 15 : if( vertDat )
493 : : {
494 : : // check the parametric space to make sure it's consistent
495 [ # # ][ # # ]: 0 : assert( elemSeq->sdata()->min_params() == HomCoord( boxDims, 3 ) &&
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # #
# # # # #
# ]
496 [ # # ]: 0 : ( elemSeq->sdata()->max_params() + HomCoord( 1, 1, 1 ) ) == HomCoord( boxDims, 3 ) );
497 : : }
498 : : else
499 : : {
500 : : // get the parametric space from the element sequence
501 [ + + ]: 60 : for( int i = 0; i < 3; i++ )
502 : : {
503 [ + - ][ + - ]: 45 : boxDims[i] = elemSeq->sdata()->min_params()[i];
[ + - ]
504 [ + - ][ + - ]: 45 : boxDims[3 + i] = elemSeq->sdata()->max_params()[i];
[ + - ]
505 : : }
506 : : }
507 : :
508 : 15 : startElem = elemSeq->start_handle();
509 : : }
510 : : else
511 : : {
512 [ + - ]: 31 : Range elems;
513 : : impl->mbImpl->get_entities_by_dimension(
514 [ + + ][ + + ]: 31 : bset, ( boxDims[2] == boxDims[5] ? ( boxDims[1] == boxDims[4] ? 1 : 2 ) : 3 ), elems );
[ + - ]
515 [ + - ][ - + ]: 31 : if( !elems.empty() ) startElem = *elems.begin();
[ # # ][ # # ]
516 : : // call the following w/o looking at return value, since it doesn't really need to be there
517 [ - + ][ # # ]: 31 : if( impl->boxPeriodicTag ) impl->mbImpl->tag_get_data( impl->boxPeriodicTag, &bset, 1, locallyPeriodic );
518 : : }
519 : :
520 [ + + ][ - + ]: 46 : assert( vertDat || elemSeq || boxDims[0] != boxDims[3] || boxDims[1] != boxDims[4] || boxDims[2] != boxDims[5] );
[ # # ][ # # ]
[ # # ]
521 : :
522 [ + - ][ + - ]: 46 : boxSize = HomCoord( boxDims + 3, 3 ) - HomCoord( boxDims, 3 ) + HomCoord( 1, 1, 1 );
[ + - ][ + - ]
[ + - ]
523 [ + - ][ + - ]: 46 : boxSizeIJ = ( boxSize[1] ? boxSize[1] : 1 ) * boxSize[0];
[ + - ][ + - ]
[ + - ][ # # ]
524 [ + - ]: 46 : boxSizeIM1 = boxSize[0] - ( locallyPeriodic[0] ? 0 : 1 );
525 [ + - ][ + - ]: 46 : boxSizeIJM1 = ( boxSize[1] ? ( boxSize[1] - ( locallyPeriodic[1] ? 0 : 1 ) ) : 1 ) * boxSizeIM1;
[ + - ][ + - ]
[ # # ]
526 : :
527 : 46 : scImpl->add_box( this );
528 : 46 : }
529 : :
530 : 46 : ScdBox::~ScdBox()
531 : : {
532 : : // Reset the tag on the set
533 [ + - ]: 46 : if( boxSet )
534 : : {
535 : : // It is possible that the box set entity has been deleted (e.g. by
536 : : // Core::clean_up_failed_read)
537 [ - + ]: 46 : Core* mbcore = dynamic_cast< Core* >( scImpl->mbImpl );
538 [ - + ]: 46 : assert( mbcore != NULL );
539 [ + - ]: 46 : if( mbcore->is_valid( boxSet ) )
540 : : {
541 : 46 : ScdBox* tmp_ptr = NULL;
542 : 46 : scImpl->mbImpl->tag_set_data( scImpl->box_set_tag(), &boxSet, 1, &tmp_ptr );
543 : : }
544 : : else
545 : 46 : boxSet = 0;
546 : : }
547 : :
548 : 46 : scImpl->remove_box( this );
549 : 46 : }
550 : :
551 : 6520 : EntityHandle ScdBox::get_vertex_from_seq( int i, int j, int k ) const
552 : : {
553 [ - + ]: 6520 : assert( elemSeq );
554 : 6520 : return elemSeq->get_vertex( i, j, k );
555 : : }
556 : :
557 : 13229 : int ScdBox::box_dimension() const
558 : : {
559 [ + + ]: 13229 : return ( startElem ? scImpl->mbImpl->dimension_from_handle( startElem ) : -1 );
560 : : }
561 : :
562 : 39 : ErrorCode ScdBox::add_vbox( ScdBox* vbox, HomCoord from1, HomCoord to1, HomCoord from2, HomCoord to2, HomCoord from3,
563 : : HomCoord to3, bool bb_input, const HomCoord& bb_min, const HomCoord& bb_max )
564 : : {
565 [ - + ]: 39 : if( !vbox->vertDat ) return MB_FAILURE;
566 : 39 : ScdVertexData* dum_data = dynamic_cast< ScdVertexData* >( vbox->vertDat );
567 : : ErrorCode rval =
568 : 39 : elemSeq->sdata()->add_vsequence( dum_data, from1, to1, from2, to2, from3, to3, bb_input, bb_min, bb_max );
569 : 39 : return rval;
570 : : }
571 : :
572 : 15 : bool ScdBox::boundary_complete() const
573 : : {
574 : 15 : return elemSeq->boundary_complete();
575 : : }
576 : :
577 : 15 : ErrorCode ScdBox::get_coordinate_arrays( double*& xc, double*& yc, double*& zc )
578 : : {
579 [ - + ]: 15 : if( !vertDat ) return MB_FAILURE;
580 : :
581 : 15 : xc = reinterpret_cast< double* >( vertDat->get_sequence_data( 0 ) );
582 : 15 : yc = reinterpret_cast< double* >( vertDat->get_sequence_data( 1 ) );
583 : 15 : zc = reinterpret_cast< double* >( vertDat->get_sequence_data( 2 ) );
584 : 15 : return MB_SUCCESS;
585 : : }
586 : :
587 : 0 : ErrorCode ScdBox::get_coordinate_arrays( const double*& xc, const double*& yc, const double*& zc ) const
588 : : {
589 [ # # ]: 0 : if( !vertDat ) return MB_FAILURE;
590 : 0 : xc = reinterpret_cast< const double* >( vertDat->get_sequence_data( 0 ) );
591 : 0 : yc = reinterpret_cast< const double* >( vertDat->get_sequence_data( 1 ) );
592 : 0 : zc = reinterpret_cast< const double* >( vertDat->get_sequence_data( 2 ) );
593 : 0 : return MB_SUCCESS;
594 : : }
595 : :
596 : 0 : ErrorCode ScdBox::vert_dat( ScdVertexData* vert_dt )
597 : : {
598 : 0 : vertDat = vert_dt;
599 : 0 : return MB_SUCCESS;
600 : : }
601 : :
602 : 13 : ErrorCode ScdBox::elem_seq( EntitySequence* elem_sq )
603 : : {
604 [ - + ]: 13 : elemSeq = dynamic_cast< StructuredElementSeq* >( elem_sq );
605 [ + - ]: 13 : if( elemSeq ) elemSeq->is_periodic( locallyPeriodic );
606 : :
607 [ + + ][ + - ]: 13 : if( locallyPeriodic[0] ) boxSizeIM1 = boxSize[0] - ( locallyPeriodic[0] ? 0 : 1 );
608 [ + + ][ + + ]: 13 : if( locallyPeriodic[0] || locallyPeriodic[1] )
609 [ + - ][ + - ]: 6 : boxSizeIJM1 = ( boxSize[1] ? ( boxSize[1] - ( locallyPeriodic[1] ? 0 : 1 ) ) : 1 ) * boxSizeIM1;
[ + - ][ + - ]
[ # # ]
610 : :
611 [ + - ]: 13 : return ( elemSeq ? MB_SUCCESS : MB_FAILURE );
612 : : }
613 : :
614 : 13229 : ErrorCode ScdBox::get_params( EntityHandle ent, HomCoord& ijkd ) const
615 : : {
616 : : // check first whether this is an intermediate entity, so we know what to do
617 : 13229 : int dimension = box_dimension();
618 : 13229 : int this_dim = scImpl->impl()->dimension_from_handle( ent );
619 : :
620 [ + + ][ + - ]: 13229 : if( ( 0 == this_dim && !vertDat ) || ( this_dim && this_dim == dimension ) )
[ + + ][ + - ]
621 : : {
622 [ - + ]: 3506 : assert( elemSeq );
623 [ + - ][ + - ]: 3506 : return elemSeq->get_params( ent, ijkd[0], ijkd[1], ijkd[2] );
[ + - ][ + - ]
624 : : }
625 : :
626 [ + - ][ + - ]: 9723 : else if( !this_dim && vertDat )
627 [ + - ][ + - ]: 9723 : return vertDat->get_params( ent, ijkd[0], ijkd[1], ijkd[2] );
[ + - ][ + - ]
628 : :
629 : : else
630 : 13229 : return MB_NOT_IMPLEMENTED;
631 : : }
632 : :
633 : : //! Get the entity of specified dimension adjacent to parametric element
634 : : /**
635 : : * \param dim Dimension of adjacent entity being requested
636 : : * \param i Parametric coordinates of cell being evaluated
637 : : * \param j Parametric coordinates of cell being evaluated
638 : : * \param k Parametric coordinates of cell being evaluated
639 : : * \param dir Direction (0, 1, or 2), for getting adjacent edges (2d, 3d) or faces (3d)
640 : : * \param ent EntityHandle of adjacent entity
641 : : * \param create_if_missing If true, creates the entity if it doesn't already exist
642 : : */
643 : 914 : ErrorCode ScdBox::get_adj_edge_or_face( int dim, int i, int j, int k, int dir, EntityHandle& ent,
644 : : bool create_if_missing ) const
645 : : {
646 : : // describe connectivity of sub-element in static array
647 : : // subconnect[dim-1][dir][numv][ijk] where dimensions are:
648 : : // [dim-1]: dim=1 or 2, so this is 0 or 1
649 : : // [dir]: one of 0..2, for ijk directions in a hex
650 : : // [numv]: number of vertices describing sub entity = 2*dim <= 4
651 : : // [ijk]: 3 values for i, j, k
652 : : int subconnect[2][3][4][3] = { { { { 0, 0, 0 }, { 1, 0, 0 }, { -1, -1, -1 }, { -1, -1, -1 } }, // i edge
653 : : { { 0, 0, 0 }, { 0, 1, 0 }, { -1, -1, -1 }, { -1, -1, -1 } }, // j edge
654 : : { { 0, 0, 0 }, { 0, 0, 1 }, { -1, -1, -1 }, { -1, -1, -1 } } }, // k edge
655 : :
656 : : { { { 0, 0, 0 }, { 0, 1, 0 }, { 0, 1, 1 }, { 0, 0, 1 } }, // i face
657 : : { { 0, 0, 0 }, { 1, 0, 0 }, { 1, 0, 1 }, { 0, 0, 1 } }, // j face
658 : 914 : { { 0, 0, 0 }, { 1, 0, 0 }, { 1, 1, 0 }, { 0, 1, 0 } } } }; // k face
659 : :
660 : : // check proper input dimensions and lower bound
661 [ + - ][ + - ]: 914 : if( dim < 1 || dim > 2 || i < boxDims[0] || j < boxDims[1] || k < boxDims[2] ) return MB_FAILURE;
[ + - ][ + - ]
[ - + ]
662 : :
663 : : // now check upper bound; parameters must be <= upper corner, since edges/faces
664 : : // follow element parameterization, not vertex parameterization
665 [ + - ][ - + ]: 914 : else if( ( boxDims[3] != boxDims[0] && i > ( locallyPeriodic[0] ? boxDims[3] + 1 : boxDims[3] ) ) ||
[ + - ][ + - ]
666 [ - + ][ + - ]: 914 : ( boxDims[4] != boxDims[1] && j > ( locallyPeriodic[1] ? boxDims[4] + 1 : boxDims[4] ) ) ||
[ + + ]
667 [ - + ]: 872 : ( boxDims[5] != boxDims[2] && k > boxDims[5] ) )
668 : 0 : return MB_FAILURE;
669 : :
670 : : // get the vertices making up this entity
671 : : EntityHandle verts[4];
672 [ + + ]: 4414 : for( int ind = 0; ind < 2 * dim; ind++ )
673 : : {
674 : 3500 : int i1 = i + subconnect[dim - 1][dir][ind][0];
675 : 3500 : int j1 = j + subconnect[dim - 1][dir][ind][1];
676 : : // if periodic in i and i1 is boxDims[3]+1, wrap around
677 [ - + ][ # # ]: 3500 : if( locallyPeriodic[0] && i1 == boxDims[3] + 1 ) i1 = boxDims[0];
678 : : // if periodic in j and j1 is boxDims[4]+1, wrap around
679 [ - + ][ # # ]: 3500 : if( locallyPeriodic[1] && j1 == boxDims[4] + 1 ) j1 = boxDims[1];
680 [ + - ]: 3500 : verts[ind] = get_vertex( i1, j1, k + subconnect[dim - 1][dir][ind][2] );
681 [ - + ]: 3500 : if( !verts[ind] ) return MB_FAILURE;
682 : : }
683 : :
684 [ + - ]: 914 : Range ents;
685 [ + - ][ + - ]: 914 : ErrorCode rval = scImpl->impl()->get_adjacencies( verts, 2 * dim, dim, false, ents );
686 [ - + ]: 914 : if( MB_SUCCESS != rval ) return rval;
687 : :
688 [ + - ][ - + ]: 914 : if( ents.size() > 1 )
689 : 0 : return MB_FAILURE;
690 : :
691 [ + - ][ + + ]: 914 : else if( ents.size() == 1 )
692 : : {
693 [ + - ][ + - ]: 411 : ent = *ents.begin();
694 : : }
695 [ + - ]: 503 : else if( create_if_missing )
696 [ + - ][ + + ]: 503 : rval = scImpl->impl()->create_element( ( 1 == dim ? MBEDGE : MBQUAD ), verts, 2 * dim, ent );
[ + - ]
697 : :
698 : 914 : return rval;
699 : : }
700 : :
701 : : #ifndef MOAB_HAVE_MPI
702 : : ErrorCode ScdInterface::tag_shared_vertices( ParallelComm*, ScdBox* )
703 : : {
704 : : return MB_FAILURE;
705 : : #else
706 : 0 : ErrorCode ScdInterface::tag_shared_vertices( ParallelComm* pcomm, ScdBox* box )
707 : : {
708 [ # # ]: 0 : EntityHandle seth = box->box_set();
709 : :
710 : : // check the # ents in the box against the num in the set, to make sure it's only 1 box;
711 : : // reuse tmp_range
712 [ # # ]: 0 : Range tmp_range;
713 [ # # ][ # # ]: 0 : ErrorCode rval = mbImpl->get_entities_by_dimension( seth, box->box_dimension(), tmp_range );
714 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
715 [ # # ][ # # ]: 0 : if( box->num_elements() != (int)tmp_range.size() ) return MB_FAILURE;
[ # # ]
716 : :
717 [ # # ]: 0 : const int* gdims = box->par_data().gDims;
718 [ # # ][ # # ]: 0 : if( ( gdims[0] == gdims[3] && gdims[1] == gdims[4] && gdims[2] == gdims[5] ) || -1 == box->par_data().partMethod )
[ # # ][ # # ]
[ # # ][ # # ]
719 : 0 : return MB_FAILURE;
720 : :
721 : : // ok, we have a partitioned box; get the vertices shared with other processors
722 [ # # ][ # # ]: 0 : std::vector< int > procs, offsets, shared_indices;
[ # # ]
723 [ # # ]: 0 : rval = get_shared_vertices( pcomm, box, procs, offsets, shared_indices );
724 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
725 : :
726 : : // post receives for start handles once we know how many to look for
727 [ # # ][ # # ]: 0 : std::vector< MPI_Request > recv_reqs( procs.size(), MPI_REQUEST_NULL ), send_reqs( procs.size(), MPI_REQUEST_NULL );
728 [ # # ][ # # ]: 0 : std::vector< EntityHandle > rhandles( 4 * procs.size() ), shandles( 4 );
729 [ # # ]: 0 : for( unsigned int i = 0; i < procs.size(); i++ )
730 : : {
731 [ # # ][ # # ]: 0 : int success = MPI_Irecv( (void*)&rhandles[4 * i], 4 * sizeof( EntityHandle ), MPI_UNSIGNED_CHAR, procs[i], 1,
732 [ # # ][ # # ]: 0 : pcomm->proc_config().proc_comm(), &recv_reqs[i] );
[ # # ][ # # ]
733 [ # # ]: 0 : if( success != MPI_SUCCESS ) return MB_FAILURE;
734 : : }
735 : :
736 : : // send our own start handles
737 [ # # ][ # # ]: 0 : shandles[0] = box->start_vertex();
738 [ # # ]: 0 : shandles[1] = 0;
739 [ # # ][ # # ]: 0 : if( box->box_dimension() == 1 )
740 : : {
741 [ # # ][ # # ]: 0 : shandles[1] = box->start_element();
742 [ # # ]: 0 : shandles[2] = 0;
743 [ # # ]: 0 : shandles[3] = 0;
744 : : }
745 [ # # ][ # # ]: 0 : else if( box->box_dimension() == 2 )
746 : : {
747 [ # # ][ # # ]: 0 : shandles[2] = box->start_element();
748 [ # # ]: 0 : shandles[3] = 0;
749 : : }
750 : : else
751 : : {
752 [ # # ]: 0 : shandles[2] = 0;
753 [ # # ][ # # ]: 0 : shandles[3] = box->start_element();
754 : : }
755 [ # # ]: 0 : for( unsigned int i = 0; i < procs.size(); i++ )
756 : : {
757 [ # # ][ # # ]: 0 : int success = MPI_Isend( (void*)&shandles[0], 4 * sizeof( EntityHandle ), MPI_UNSIGNED_CHAR, procs[i], 1,
758 [ # # ][ # # ]: 0 : pcomm->proc_config().proc_comm(), &send_reqs[i] );
[ # # ][ # # ]
759 [ # # ]: 0 : if( success != MPI_SUCCESS ) return MB_FAILURE;
760 : : }
761 : :
762 : : // receive start handles and save info to a tuple list
763 : 0 : int incoming = procs.size();
764 : : int p, j, k;
765 : : MPI_Status status;
766 [ # # ]: 0 : TupleList shared_data;
767 [ # # ]: 0 : shared_data.initialize( 1, 0, 2, 0, shared_indices.size() / 2 );
768 [ # # ]: 0 : shared_data.enableWriteAccess();
769 : :
770 : 0 : j = 0;
771 : 0 : k = 0;
772 [ # # ]: 0 : while( incoming )
773 : : {
774 [ # # ][ # # ]: 0 : int success = MPI_Waitany( procs.size(), &recv_reqs[0], &p, &status );
775 [ # # ]: 0 : if( MPI_SUCCESS != success ) return MB_FAILURE;
776 [ # # ][ # # ]: 0 : unsigned int num_indices = ( offsets[p + 1] - offsets[p] ) / 2;
777 [ # # ][ # # ]: 0 : int *lh = &shared_indices[offsets[p]], *rh = lh + num_indices;
778 [ # # ]: 0 : for( unsigned int i = 0; i < num_indices; i++ )
779 : : {
780 [ # # ]: 0 : shared_data.vi_wr[j++] = procs[p];
781 [ # # ]: 0 : shared_data.vul_wr[k++] = shandles[0] + lh[i];
782 [ # # ]: 0 : shared_data.vul_wr[k++] = rhandles[4 * p] + rh[i];
783 [ # # ]: 0 : shared_data.inc_n();
784 : : }
785 : 0 : incoming--;
786 : : }
787 : :
788 : : // still need to wait for the send requests
789 [ # # ]: 0 : std::vector< MPI_Status > mult_status( procs.size() );
790 [ # # ][ # # ]: 0 : int success = MPI_Waitall( procs.size(), &send_reqs[0], &mult_status[0] );
[ # # ]
791 [ # # ][ # # ]: 0 : if( MPI_SUCCESS != success ) { MB_SET_ERR( MB_FAILURE, "Failed in waitall in ScdInterface::tag_shared_vertices" ); }
[ # # ][ # # ]
[ # # ][ # # ]
792 : : // sort by local handle
793 [ # # ]: 0 : TupleList::buffer sort_buffer;
794 [ # # ]: 0 : sort_buffer.buffer_init( shared_indices.size() / 2 );
795 [ # # ]: 0 : shared_data.sort( 1, &sort_buffer );
796 [ # # ]: 0 : sort_buffer.reset();
797 : :
798 : : // process into sharing data
799 [ # # ]: 0 : std::map< std::vector< int >, std::vector< EntityHandle > > proc_nvecs;
800 [ # # ]: 0 : Range dum;
801 [ # # ]: 0 : rval = pcomm->tag_shared_verts( shared_data, proc_nvecs, dum, 0 );
802 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
803 : :
804 : : // create interface sets
805 [ # # ]: 0 : rval = pcomm->create_interface_sets( proc_nvecs );
806 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
807 : :
808 : : // add the box to the PComm's partitionSets
809 [ # # ][ # # ]: 0 : pcomm->partition_sets().insert( box->box_set() );
[ # # ]
810 : :
811 : : // make sure buffers are allocated for communicating procs
812 [ # # ][ # # ]: 0 : for( std::vector< int >::iterator pit = procs.begin(); pit != procs.end(); ++pit )
[ # # ]
813 [ # # ][ # # ]: 0 : pcomm->get_buffers( *pit );
814 : :
815 [ # # ][ # # ]: 0 : if( pcomm->get_debug_verbosity() > 1 ) pcomm->list_entities( NULL, 1 );
[ # # ]
816 : :
817 : : #ifndef NDEBUG
818 [ # # ]: 0 : rval = pcomm->check_all_shared_handles();
819 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
820 : : #endif
821 : :
822 : 0 : return MB_SUCCESS;
823 : :
824 : : #endif
825 : : }
826 : :
827 : 299 : ErrorCode ScdInterface::get_neighbor_alljkbal( int np, int pfrom, const int* const gdims, const int* const gperiodic,
828 : : const int* const dijk, int& pto, int* rdims, int* facedims,
829 : : int* across_bdy )
830 : : {
831 [ + + ]: 299 : if( dijk[0] != 0 )
832 : : {
833 : 180 : pto = -1;
834 : 180 : return MB_SUCCESS;
835 : : }
836 : :
837 : 119 : pto = -1;
838 : 119 : across_bdy[0] = across_bdy[1] = across_bdy[2] = 0;
839 : :
840 : : int ldims[6], pijk[3], lperiodic[3];
841 [ + - ]: 119 : ErrorCode rval = compute_partition_alljkbal( np, pfrom, gdims, gperiodic, ldims, lperiodic, pijk );
842 [ - + ]: 119 : if( MB_SUCCESS != rval ) return rval;
843 [ - + ]: 119 : assert( pijk[1] * pijk[2] == np );
844 : 119 : pto = -1;
845 : 119 : bool bot_j = pfrom< pijk[2], top_j = pfrom > np - pijk[2];
846 [ + + ][ + + ]: 119 : if( ( 1 == pijk[2] && dijk[2] ) || // 1d in j means no neighbors with dk != 0
[ + + ]
847 [ + + ][ + + ]: 107 : ( !( pfrom % pijk[2] ) && -1 == dijk[2] ) || // at -k bdy
848 [ + + ][ + + ]: 95 : ( pfrom % pijk[2] == pijk[2] - 1 && 1 == dijk[2] ) || // at +k bdy
849 [ + + ][ - + ]: 83 : ( pfrom < pijk[2] && -1 == dijk[1] && !gperiodic[1] ) || // down and not periodic
[ - + ]
850 [ # # ][ # # ]: 0 : ( pfrom >= np - pijk[2] && 1 == dijk[1] && !gperiodic[1] ) ) // up and not periodic
851 : 41 : return MB_SUCCESS;
852 : :
853 : 78 : pto = pfrom;
854 [ + - ]: 78 : std::copy( ldims, ldims + 6, rdims );
855 [ + - ]: 78 : std::copy( ldims, ldims + 6, facedims );
856 : :
857 [ + + ]: 78 : if( 0 != dijk[1] )
858 : : {
859 : 62 : pto = ( pto + dijk[1] * pijk[2] + np ) % np;
860 [ + - ][ - + ]: 62 : assert( pto >= 0 && pto < np );
861 : 62 : int dj = ( gdims[4] - gdims[1] ) / pijk[1], extra = ( gdims[4] - gdims[1] ) % pijk[1];
862 [ + + ]: 62 : if( -1 == dijk[1] )
863 : : {
864 : 31 : facedims[4] = facedims[1];
865 [ - + ]: 31 : if( bot_j )
866 : : {
867 : : // going across periodic lower bdy in j
868 : 0 : rdims[4] = gdims[4];
869 : 0 : across_bdy[1] = -1;
870 : : }
871 : : else
872 : : {
873 : 31 : rdims[4] = ldims[1];
874 : : }
875 : 31 : rdims[1] = rdims[4] - dj;
876 [ - + ]: 31 : if( pto < extra ) rdims[1]--;
877 : : }
878 : : else
879 : : {
880 [ - + ]: 31 : if( pfrom > np - pijk[2] ) facedims[4] = gdims[1];
881 : 31 : facedims[1] = facedims[4];
882 [ - + ]: 31 : if( top_j )
883 : : {
884 : : // going across periodic upper bdy in j
885 : 0 : rdims[1] = gdims[1];
886 : 0 : across_bdy[1] = 1;
887 : : }
888 : : else
889 : : {
890 : 31 : rdims[1] = ldims[4];
891 : : }
892 : 31 : rdims[4] = rdims[1] + dj;
893 [ - + ]: 62 : if( pto < extra ) rdims[4]++;
894 : : }
895 : : }
896 [ + + ]: 78 : if( 0 != dijk[2] )
897 : : {
898 : 44 : pto = ( pto + dijk[2] ) % np;
899 [ + - ][ - + ]: 44 : assert( pto >= 0 && pto < np );
900 [ + + ]: 44 : facedims[2] = facedims[5] = ( -1 == dijk[2] ? facedims[2] : facedims[5] );
901 : 44 : int dk = ( gdims[5] - gdims[2] ) / pijk[2];
902 [ + + ]: 44 : if( -1 == dijk[2] )
903 : : {
904 : 22 : facedims[5] = facedims[2];
905 : 22 : rdims[5] = ldims[2];
906 : 22 : rdims[2] = rdims[5] - dk; // never any kextra for alljkbal
907 : : }
908 : : else
909 : : {
910 : 22 : facedims[2] = facedims[5];
911 : 22 : rdims[2] = ldims[5];
912 : 44 : rdims[5] = rdims[2] + dk; // never any kextra for alljkbal
913 : : }
914 : : }
915 : :
916 [ + - ][ + - ]: 78 : assert( -1 == pto || ( rdims[0] >= gdims[0] && rdims[3] <= gdims[3] ) );
[ - + ]
917 [ + - ][ + - ]: 78 : assert( -1 == pto || ( rdims[1] >= gdims[1] && ( rdims[4] <= gdims[4] || ( across_bdy[1] && bot_j ) ) ) );
[ - + ][ # # ]
[ # # ]
918 [ + - ][ + - ]: 78 : assert( -1 == pto || ( rdims[2] >= gdims[2] && rdims[5] <= gdims[5] ) );
[ - + ]
919 [ - + ][ # # ]: 78 : assert( -1 == pto || ( ( facedims[0] >= rdims[0] ||
[ # # ][ # # ]
920 [ + - ]: 78 : ( gperiodic[0] && rdims[3] == gdims[3] + 1 && facedims[0] == gdims[0] ) ) ) );
921 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[3] <= rdims[3] ) );
922 [ - + ][ # # ]: 78 : assert( -1 == pto || ( ( facedims[1] >= rdims[1] ||
[ # # ][ # # ]
923 [ + - ]: 78 : ( gperiodic[1] && rdims[4] == gdims[4] + 1 && facedims[1] == gdims[1] ) ) ) );
924 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[4] <= rdims[4] ) );
925 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[2] >= rdims[2] ) );
926 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[5] <= rdims[5] ) );
927 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[0] >= ldims[0] ) );
928 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[3] <= ldims[3] ) );
929 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[1] >= ldims[1] ) );
930 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[4] <= ldims[4] ) );
931 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[2] >= ldims[2] ) );
932 [ + - ][ - + ]: 78 : assert( -1 == pto || ( facedims[5] <= ldims[5] ) );
933 : :
934 : 299 : return MB_SUCCESS;
935 : : }
936 : :
937 : 308 : ErrorCode ScdInterface::get_neighbor_sqij( int np, int pfrom, const int* const gdims, const int* const gperiodic,
938 : : const int* const dijk, int& pto, int* rdims, int* facedims, int* across_bdy )
939 : : {
940 [ + + ]: 308 : if( dijk[2] != 0 )
941 : : {
942 : : // for sqij, there is no k neighbor, ever
943 : 180 : pto = -1;
944 : 180 : return MB_SUCCESS;
945 : : }
946 : :
947 : 128 : pto = -1;
948 : 128 : across_bdy[0] = across_bdy[1] = across_bdy[2] = 0;
949 : : int lperiodic[3], pijk[3], ldims[6];
950 [ + - ]: 128 : ErrorCode rval = compute_partition_sqij( np, pfrom, gdims, gperiodic, ldims, lperiodic, pijk );
951 [ - + ]: 128 : if( MB_SUCCESS != rval ) return rval;
952 [ - + ]: 128 : assert( pijk[0] * pijk[1] == np );
953 : 128 : pto = -1;
954 : 128 : bool top_i = 0, top_j = 0, bot_i = 0, bot_j = 0;
955 : 128 : int ni = pfrom % pijk[0], nj = pfrom / pijk[0]; // row / column number of me
956 [ + + ]: 128 : if( ni == pijk[0] - 1 ) top_i = 1;
957 [ + + ]: 128 : if( nj == pijk[1] - 1 ) top_j = 1;
958 [ + + ]: 128 : if( !ni ) bot_i = 1;
959 [ + + ]: 128 : if( !nj ) bot_j = 1;
960 [ + - ][ + + ]: 128 : if( ( !gperiodic[0] && bot_i && -1 == dijk[0] ) || // left and not periodic
[ + + ][ + - ]
961 [ + + ][ + + ]: 119 : ( !gperiodic[0] && top_i && 1 == dijk[0] ) || // right and not periodic
[ + - ]
962 [ + + ][ + + ]: 110 : ( !gperiodic[1] && bot_j && -1 == dijk[1] ) || // bottom and not periodic
[ + - ]
963 [ + + ][ - + ]: 96 : ( !gperiodic[1] && top_j && 1 == dijk[1] ) ) // top and not periodic
964 : 32 : return MB_SUCCESS;
965 : :
966 [ + - ]: 96 : std::copy( ldims, ldims + 6, facedims );
967 [ + - ]: 96 : std::copy( ldims, ldims + 6, rdims );
968 : 96 : pto = pfrom;
969 : 96 : int j = gdims[4] - gdims[1], dj = j / pijk[1], jextra = ( gdims[4] - gdims[1] ) % dj, i = gdims[3] - gdims[0],
970 : 96 : di = i / pijk[0], iextra = ( gdims[3] - gdims[0] ) % di;
971 : :
972 [ + + ]: 96 : if( 0 != dijk[0] )
973 : : {
974 : 68 : pto = ( ni + dijk[0] + pijk[0] ) % pijk[0]; // get pto's ni value
975 : 68 : pto = nj * pijk[0] + pto; // then convert to pto
976 [ + - ][ - + ]: 68 : assert( pto >= 0 && pto < np );
977 [ + + ]: 68 : if( -1 == dijk[0] )
978 : : {
979 : 34 : facedims[3] = facedims[0];
980 [ - + ]: 34 : if( bot_i )
981 : : {
982 : : // going across lower periodic bdy in i
983 : 0 : across_bdy[0] = -1;
984 : 0 : rdims[3] = gdims[3] + 1; // +1 because ldims[3] on remote proc is gdims[3]+1
985 : 0 : rdims[0] = rdims[3] - di - 1; // -1 to account for rdims[3] being one larger
986 : : }
987 : : else
988 : : {
989 : 34 : rdims[3] = ldims[0];
990 : 34 : rdims[0] = rdims[3] - di;
991 : : }
992 : :
993 [ - + ]: 34 : if( pto % pijk[0] < iextra ) rdims[0]--;
994 : : }
995 : : else
996 : : {
997 [ - + ]: 34 : if( top_i )
998 : : {
999 : : // going across lower periodic bdy in i
1000 : 0 : facedims[3] = gdims[0];
1001 : 0 : across_bdy[0] = 1;
1002 : : }
1003 : 34 : facedims[0] = facedims[3];
1004 [ - + ]: 34 : rdims[0] = ( top_i ? gdims[0] : ldims[3] );
1005 : 34 : rdims[3] = rdims[0] + di;
1006 [ - + ]: 34 : if( pto % pijk[0] < iextra ) rdims[3]++;
1007 [ - + ][ # # ]: 34 : if( gperiodic[0] && ni == pijk[0] - 2 ) rdims[3]++; // remote proc is top_i and periodic
1008 : : }
1009 : : }
1010 [ + + ]: 96 : if( 0 != dijk[1] )
1011 : : {
1012 : 68 : pto = ( pto + dijk[1] * pijk[0] + np ) % np;
1013 [ + - ][ - + ]: 68 : assert( pto >= 0 && pto < np );
1014 [ + + ]: 68 : if( -1 == dijk[1] )
1015 : : {
1016 : 34 : facedims[4] = facedims[1];
1017 [ - + ]: 34 : if( bot_j )
1018 : : {
1019 : : // going across lower periodic bdy in j
1020 : 0 : rdims[4] = gdims[4] + 1; // +1 because ldims[4] on remote proc is gdims[4]+1
1021 : 0 : rdims[1] = rdims[4] - dj - 1; // -1 to account for gdims[4] being one larger
1022 : 0 : across_bdy[1] = -1;
1023 : : }
1024 : : else
1025 : : {
1026 : 34 : rdims[4] = ldims[1];
1027 : 34 : rdims[1] = rdims[4] - dj;
1028 : : }
1029 [ - + ]: 34 : if( pto / pijk[0] < jextra ) rdims[1]--;
1030 : : }
1031 : : else
1032 : : {
1033 [ - + ]: 34 : if( top_j )
1034 : : {
1035 : : // going across lower periodic bdy in j
1036 : 0 : facedims[4] = gdims[1];
1037 : 0 : rdims[1] = gdims[1];
1038 : 0 : across_bdy[1] = 1;
1039 : : }
1040 : : else
1041 : : {
1042 : 34 : rdims[1] = ldims[4];
1043 : : }
1044 : 34 : facedims[1] = facedims[4];
1045 : 34 : rdims[4] = rdims[1] + dj;
1046 [ - + ]: 34 : if( nj + 1 < jextra ) rdims[4]++;
1047 [ - + ][ # # ]: 34 : if( gperiodic[1] && nj == pijk[1] - 2 ) rdims[4]++; // remote proc is top_j and periodic
1048 : : }
1049 : : }
1050 : :
1051 : : // rdims within gdims
1052 [ + - ][ - + ]: 96 : assert( -1 == pto || ( rdims[0] >= gdims[0] &&
[ # # ][ - + ]
1053 [ + - ]: 96 : ( rdims[3] <= gdims[3] + ( gperiodic[0] && pto % pijk[0] == pijk[0] - 1 ? 1 : 0 ) ) ) );
1054 [ + - ][ - + ]: 96 : assert( -1 == pto || ( rdims[1] >= gdims[1] &&
[ # # ][ - + ]
1055 [ + - ]: 96 : ( rdims[4] <= gdims[4] + ( gperiodic[1] && pto / pijk[0] == pijk[1] - 1 ? 1 : 0 ) ) ) );
1056 [ + - ][ + - ]: 96 : assert( -1 == pto || ( rdims[2] >= gdims[2] && rdims[5] <= gdims[5] ) );
[ - + ]
1057 : : // facedims within rdims
1058 [ - + ][ # # ]: 96 : assert( -1 == pto || ( ( facedims[0] >= rdims[0] ||
[ # # ][ # # ]
1059 [ + - ]: 96 : ( gperiodic[0] && pto % pijk[0] == pijk[0] - 1 && facedims[0] == gdims[0] ) ) ) );
1060 [ + - ][ - + ]: 96 : assert( -1 == pto || ( facedims[3] <= rdims[3] ) );
1061 [ - + ][ # # ]: 96 : assert( -1 == pto || ( ( facedims[1] >= rdims[1] ||
[ # # ][ # # ]
1062 [ + - ]: 96 : ( gperiodic[1] && pto / pijk[0] == pijk[1] - 1 && facedims[1] == gdims[1] ) ) ) );
1063 [ + - ][ - + ]: 96 : assert( -1 == pto || ( facedims[4] <= rdims[4] ) );
1064 [ + - ][ + - ]: 96 : assert( -1 == pto || ( facedims[2] >= rdims[2] && facedims[5] <= rdims[5] ) );
[ - + ]
1065 : : // facedims within ldims
1066 [ + - ][ - + ]: 96 : assert( -1 == pto || ( ( facedims[0] >= ldims[0] || ( top_i && facedims[0] == gdims[0] ) ) ) );
[ # # ][ # # ]
1067 [ + - ][ - + ]: 96 : assert( -1 == pto || ( facedims[3] <= ldims[3] ) );
1068 [ + - ][ - + ]: 96 : assert( -1 == pto || ( ( facedims[1] >= ldims[1] || ( gperiodic[1] && top_j && facedims[1] == gdims[1] ) ) ) );
[ # # ][ # # ]
[ # # ]
1069 [ + - ][ - + ]: 96 : assert( -1 == pto || ( facedims[4] <= ldims[4] ) );
1070 [ + - ][ + - ]: 96 : assert( -1 == pto || ( facedims[2] >= ldims[2] && facedims[5] <= ldims[5] ) );
[ - + ]
1071 : :
1072 : 308 : return MB_SUCCESS;
1073 : : }
1074 : :
1075 : 278 : ErrorCode ScdInterface::get_neighbor_sqjk( int np, int pfrom, const int* const gdims, const int* const gperiodic,
1076 : : const int* const dijk, int& pto, int* rdims, int* facedims, int* across_bdy )
1077 : : {
1078 [ + + ]: 278 : if( dijk[0] != 0 )
1079 : : {
1080 : 180 : pto = -1;
1081 : 180 : return MB_SUCCESS;
1082 : : }
1083 : :
1084 : 98 : pto = -1;
1085 : 98 : across_bdy[0] = across_bdy[1] = across_bdy[2] = 0;
1086 : : int pijk[3], lperiodic[3], ldims[6];
1087 [ + - ]: 98 : ErrorCode rval = compute_partition_sqjk( np, pfrom, gdims, gperiodic, ldims, lperiodic, pijk );
1088 [ - + ]: 98 : if( MB_SUCCESS != rval ) return rval;
1089 [ - + ]: 98 : assert( pijk[1] * pijk[2] == np );
1090 : 98 : pto = -1;
1091 : 98 : bool top_j = 0, top_k = 0, bot_j = 0, bot_k = 0;
1092 : 98 : int nj = pfrom % pijk[1], nk = pfrom / pijk[1];
1093 [ + - ]: 98 : if( nj == pijk[1] - 1 ) top_j = 1;
1094 [ - + ]: 98 : if( nk == pijk[2] - 1 ) top_k = 1;
1095 [ + - ]: 98 : if( !nj ) bot_j = 1;
1096 [ + + ]: 98 : if( !nk ) bot_k = 1;
1097 [ + - ][ + - ]: 98 : if( ( !gperiodic[1] && bot_j && -1 == dijk[1] ) || // down and not periodic
[ + + ][ + - ]
1098 [ + - ][ + + ]: 68 : ( !gperiodic[1] && top_j && 1 == dijk[1] ) || // up and not periodic
[ + + ]
1099 [ + + ][ - + ]: 38 : ( bot_k && -1 == dijk[2] ) || // k- bdy
1100 [ # # ]: 0 : ( top_k && 1 == dijk[2] ) ) // k+ bdy
1101 : 62 : return MB_SUCCESS;
1102 : :
1103 [ + - ]: 36 : std::copy( ldims, ldims + 6, facedims );
1104 [ + - ]: 36 : std::copy( ldims, ldims + 6, rdims );
1105 : 36 : pto = pfrom;
1106 : 36 : int dj = ( gdims[4] - gdims[1] ) / pijk[1], jextra = ( gdims[4] - gdims[1] ) % dj,
1107 [ + - ]: 36 : dk = ( gdims[5] == gdims[2] ? 0 : ( gdims[5] - gdims[2] ) / pijk[2] ),
1108 : 36 : kextra = ( gdims[5] - gdims[2] ) - dk * pijk[2];
1109 [ - + ]: 36 : assert( ( dj * pijk[1] + jextra == ( gdims[4] - gdims[1] ) ) &&
1110 [ + - ]: 36 : ( dk * pijk[2] + kextra == ( gdims[5] - gdims[2] ) ) );
1111 [ - + ]: 36 : if( 0 != dijk[1] )
1112 : : {
1113 : 0 : pto = ( nj + dijk[1] + pijk[1] ) % pijk[1]; // get pto's ni value
1114 : 0 : pto = nk * pijk[1] + pto; // then convert to pto
1115 [ # # ][ # # ]: 0 : assert( pto >= 0 && pto < np );
1116 [ # # ]: 0 : if( -1 == dijk[1] )
1117 : : {
1118 : 0 : facedims[4] = facedims[1];
1119 [ # # ]: 0 : if( bot_j )
1120 : : {
1121 : : // going across lower periodic bdy in j
1122 : 0 : rdims[4] = gdims[4] + 1; // +1 because ldims[4] on remote proc is gdims[4]+1
1123 : 0 : across_bdy[1] = -1;
1124 : : }
1125 : : else
1126 : : {
1127 : 0 : rdims[4] = ldims[1];
1128 : : }
1129 : 0 : rdims[1] = rdims[4] - dj;
1130 [ # # ]: 0 : if( nj < jextra ) rdims[1]--;
1131 : : }
1132 : : else
1133 : : {
1134 [ # # ]: 0 : if( top_j )
1135 : : {
1136 : : // going across upper periodic bdy in j
1137 : 0 : rdims[1] = gdims[1];
1138 : 0 : facedims[4] = gdims[1];
1139 : 0 : across_bdy[1] = 1;
1140 : : }
1141 : : else
1142 : : {
1143 : 0 : rdims[1] = ldims[4];
1144 : : }
1145 : 0 : facedims[1] = facedims[4];
1146 : 0 : rdims[4] = rdims[1] + dj;
1147 [ # # ]: 0 : if( nj < jextra ) rdims[4]++;
1148 [ # # ][ # # ]: 0 : if( gperiodic[1] && nj == dijk[1] - 2 ) rdims[4]++; // +1 because next proc is on periodic bdy
1149 : : }
1150 : : }
1151 [ + - ]: 36 : if( 0 != dijk[2] )
1152 : : {
1153 : 36 : pto = ( pto + dijk[2] * pijk[1] + np ) % np;
1154 [ + - ][ - + ]: 36 : assert( pto >= 0 && pto < np );
1155 [ + + ]: 36 : if( -1 == dijk[2] )
1156 : : {
1157 : 18 : facedims[5] = facedims[2];
1158 : 18 : rdims[5] = ldims[2];
1159 : 18 : rdims[2] -= dk;
1160 [ + + ]: 18 : if( pto / pijk[1] < kextra ) rdims[2]--;
1161 : : }
1162 : : else
1163 : : {
1164 : 18 : facedims[2] = facedims[5];
1165 : 18 : rdims[2] = ldims[5];
1166 : 18 : rdims[5] += dk;
1167 [ + + ]: 18 : if( pto / pijk[1] < kextra ) rdims[5]++;
1168 : : }
1169 : : }
1170 : :
1171 [ + - ][ + - ]: 36 : assert( -1 == pto || ( rdims[0] >= gdims[0] && rdims[3] <= gdims[3] ) );
[ - + ]
1172 [ + - ][ + - ]: 36 : assert( -1 == pto || ( rdims[1] >= gdims[1] && ( rdims[4] <= gdims[4] || ( across_bdy[1] && bot_j ) ) ) );
[ - + ][ # # ]
[ # # ]
1173 [ + - ][ + - ]: 36 : assert( -1 == pto || ( rdims[2] >= gdims[2] && rdims[5] <= gdims[5] ) );
[ - + ]
1174 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[0] >= rdims[0] && facedims[3] <= rdims[3] ) );
[ - + ]
1175 [ - + ][ # # ]: 36 : assert( -1 == pto ||
[ # # ][ # # ]
1176 [ + - ]: 36 : ( ( facedims[1] >= rdims[1] || ( gperiodic[1] && rdims[4] == gdims[4] && facedims[1] == gdims[1] ) ) ) );
1177 [ + - ][ - + ]: 36 : assert( -1 == pto || ( facedims[4] <= rdims[4] ) );
1178 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[2] >= rdims[2] && facedims[5] <= rdims[5] ) );
[ - + ]
1179 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[0] >= ldims[0] && facedims[3] <= ldims[3] ) );
[ - + ]
1180 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[1] >= ldims[1] && facedims[4] <= ldims[4] ) );
[ - + ]
1181 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[2] >= ldims[2] && facedims[5] <= ldims[5] ) );
[ - + ]
1182 : :
1183 : 278 : return MB_SUCCESS;
1184 : : }
1185 : :
1186 : 308 : ErrorCode ScdInterface::get_neighbor_sqijk( int np, int pfrom, const int* const gdims, const int* const gperiodic,
1187 : : const int* const dijk, int& pto, int* rdims, int* facedims,
1188 : : int* across_bdy )
1189 : : {
1190 [ + - ][ + - ]: 308 : if( gperiodic[0] || gperiodic[1] || gperiodic[2] ) return MB_FAILURE;
[ - + ]
1191 : :
1192 : 308 : pto = -1;
1193 : 308 : across_bdy[0] = across_bdy[1] = across_bdy[2] = 0;
1194 : : int pijk[3], lperiodic[3], ldims[6];
1195 [ + - ]: 308 : ErrorCode rval = compute_partition_sqijk( np, pfrom, gdims, gperiodic, ldims, lperiodic, pijk );
1196 [ - + ]: 308 : if( MB_SUCCESS != rval ) return rval;
1197 [ - + ]: 308 : assert( pijk[0] * pijk[1] * pijk[2] == np );
1198 : 308 : pto = -1;
1199 : 308 : bool top[3] = { false, false, false }, bot[3] = { false, false, false };
1200 : : // nijk: rank in i/j/k direction
1201 : 308 : int nijk[3] = { pfrom % pijk[0], ( pfrom % ( pijk[0] * pijk[1] ) ) / pijk[0], pfrom / ( pijk[0] * pijk[1] ) };
1202 : :
1203 [ + + ]: 870 : for( int i = 0; i < 3; i++ )
1204 : : {
1205 [ + + ]: 774 : if( nijk[i] == pijk[i] - 1 ) top[i] = true;
1206 [ + + ]: 774 : if( !nijk[i] ) bot[i] = true;
1207 [ + - ][ + + ]: 774 : if( ( !gperiodic[i] && bot[i] && -1 == dijk[i] ) || // downward && not periodic
[ + + ][ + - ]
1208 [ + + ][ + + ]: 647 : ( !gperiodic[i] && top[i] && 1 == dijk[i] ) ) // upward && not periodic
1209 : 212 : return MB_SUCCESS;
1210 : : }
1211 : :
1212 [ + - ]: 96 : std::copy( ldims, ldims + 6, facedims );
1213 [ + - ]: 96 : std::copy( ldims, ldims + 6, rdims );
1214 : 96 : pto = pfrom;
1215 : : int delijk[3], extra[3];
1216 : : // nijk_to: rank of pto in i/j/k direction
1217 : : int nijk_to[3];
1218 [ + + ]: 384 : for( int i = 0; i < 3; i++ )
1219 : : {
1220 [ + - ]: 288 : delijk[i] = ( gdims[i + 3] == gdims[i] ? 0 : ( gdims[i + 3] - gdims[i] ) / pijk[i] );
1221 : 288 : extra[i] = ( gdims[i + 3] - gdims[i] ) % delijk[i];
1222 : 288 : nijk_to[i] = ( nijk[i] + dijk[i] + pijk[i] ) % pijk[i];
1223 : : }
1224 : 96 : pto = nijk_to[2] * pijk[0] * pijk[1] + nijk_to[1] * pijk[0] + nijk_to[0];
1225 [ + - ][ - + ]: 96 : assert( pto >= 0 && pto < np );
1226 [ + + ]: 384 : for( int i = 0; i < 3; i++ )
1227 : : {
1228 [ + + ]: 288 : if( 0 != dijk[i] )
1229 : : {
1230 [ + + ]: 136 : if( -1 == dijk[i] )
1231 : : {
1232 : 68 : facedims[i + 3] = facedims[i];
1233 [ - + ]: 68 : if( bot[i] )
1234 : : {
1235 : : // going across lower periodic bdy in i
1236 : 0 : rdims[i + 3] = gdims[i + 3] + 1; // +1 because ldims[4] on remote proc is gdims[4]+1
1237 : 0 : across_bdy[i] = -1;
1238 : : }
1239 : : else
1240 : : {
1241 : 68 : rdims[i + 3] = ldims[i];
1242 : : }
1243 : 68 : rdims[i] = rdims[i + 3] - delijk[i];
1244 [ - + ]: 68 : if( nijk[i] < extra[i] ) rdims[i]--;
1245 : : }
1246 : : else
1247 : : {
1248 [ - + ]: 68 : if( top[i] )
1249 : : {
1250 : : // going across upper periodic bdy in i
1251 : 0 : rdims[i] = gdims[i];
1252 : 0 : facedims[i + 3] = gdims[i];
1253 : 0 : across_bdy[i] = 1;
1254 : : }
1255 : : else
1256 : : {
1257 : 68 : rdims[i] = ldims[i + 3];
1258 : : }
1259 : 68 : facedims[i] = facedims[i + 3];
1260 : 68 : rdims[i + 3] = rdims[i] + delijk[i];
1261 [ - + ]: 68 : if( nijk[i] < extra[i] ) rdims[i + 3]++;
1262 [ - + ][ # # ]: 68 : if( gperiodic[i] && nijk[i] == dijk[i] - 2 ) rdims[i + 3]++; // +1 because next proc is on periodic bdy
1263 : : }
1264 : : }
1265 : : }
1266 : :
1267 [ - + ]: 96 : assert( -1 != pto );
1268 : : #ifndef NDEBUG
1269 [ + + ]: 384 : for( int i = 0; i < 3; i++ )
1270 : : {
1271 [ + - ][ - + ]: 288 : assert( ( rdims[i] >= gdims[i] && ( rdims[i + 3] <= gdims[i + 3] || ( across_bdy[i] && bot[i] ) ) ) );
[ # # ][ # # ]
1272 [ # # ][ # # ]: 0 : assert( ( ( facedims[i] >= rdims[i] ||
[ # # ]
1273 [ - + ]: 288 : ( gperiodic[i] && rdims[i + 3] == gdims[i + 3] && facedims[i] == gdims[i] ) ) ) );
1274 [ + - ][ - + ]: 288 : assert( ( facedims[i] >= ldims[i] && facedims[i + 3] <= ldims[i + 3] ) );
1275 : : }
1276 : : #endif
1277 : :
1278 : 308 : return MB_SUCCESS;
1279 : : }
1280 : :
1281 : 278 : ErrorCode ScdInterface::get_neighbor_alljorkori( int np, int pfrom, const int* const gdims, const int* const gperiodic,
1282 : : const int* const dijk, int& pto, int* rdims, int* facedims,
1283 : : int* across_bdy )
1284 : : {
1285 : 278 : ErrorCode rval = MB_SUCCESS;
1286 : 278 : pto = -1;
1287 [ - + ]: 278 : if( np == 1 ) return MB_SUCCESS;
1288 : :
1289 : : int pijk[3], lperiodic[3], ldims[6];
1290 [ + - ]: 278 : rval = compute_partition_alljorkori( np, pfrom, gdims, gperiodic, ldims, lperiodic, pijk );
1291 [ - + ]: 278 : if( MB_SUCCESS != rval ) return rval;
1292 : :
1293 : 278 : int ind = -1;
1294 : 278 : across_bdy[0] = across_bdy[1] = across_bdy[2] = 0;
1295 : :
1296 [ + - ]: 556 : for( int i = 0; i < 3; i++ )
1297 : : {
1298 [ + + ]: 556 : if( pijk[i] > 1 )
1299 : : {
1300 : 278 : ind = i;
1301 : 278 : break;
1302 : : }
1303 : : }
1304 : :
1305 [ - + ]: 278 : assert( -1 < ind );
1306 : :
1307 [ + + ]: 278 : if( !dijk[ind] )
1308 : : // no neighbor, pto is already -1, return
1309 : 80 : return MB_SUCCESS;
1310 : :
1311 [ - + ][ # # ]: 198 : bool is_periodic = ( ( gperiodic[0] && ind == 0 ) || ( gperiodic[1] && ind == 1 ) );
[ - + ][ # # ]
1312 [ + + ][ + + ]: 198 : if( dijk[( ind + 1 ) % 3] || dijk[( ind + 2 ) % 3] || // stepping in either other two directions
[ + - ]
1313 [ + + ][ + + ]: 38 : ( !is_periodic && ldims[ind] == gdims[ind] && dijk[ind] == -1 ) || // lower side and going lower
[ + - ]
1314 [ - + ][ # # ]: 36 : ( !is_periodic && ldims[3 + ind] >= gdims[3 + ind] &&
1315 : 0 : dijk[ind] == 1 ) ) // not >= because ldims is only > gdims when periodic;
1316 : : // higher side and going higher
1317 : 162 : return MB_SUCCESS;
1318 : :
1319 [ + - ]: 36 : std::copy( ldims, ldims + 6, facedims );
1320 [ + - ]: 36 : std::copy( ldims, ldims + 6, rdims );
1321 : :
1322 : 36 : int dind = ( gdims[ind + 3] - gdims[ind] ) / np;
1323 : 36 : int extra = ( gdims[ind + 3] - gdims[ind] ) % np;
1324 [ + + ][ + - ]: 36 : if( -1 == dijk[ind] && pfrom )
1325 : : {
1326 : : // actual left neighbor
1327 : 18 : pto = pfrom - 1; // no need for %np, because pfrom > 0
1328 : 18 : facedims[ind + 3] = facedims[ind];
1329 : 18 : rdims[ind + 3] = ldims[ind];
1330 [ + + ]: 18 : rdims[ind] = rdims[ind + 3] - dind - ( pto < extra ? 1 : 0 );
1331 : : }
1332 [ + - ][ + - ]: 18 : else if( 1 == dijk[ind] && pfrom < np - 1 )
1333 : : {
1334 : : // actual right neighbor
1335 : 18 : pto = pfrom + 1;
1336 : 18 : facedims[ind] = facedims[ind + 3];
1337 : 18 : rdims[ind] = ldims[ind + 3];
1338 [ + + ]: 18 : rdims[ind + 3] = rdims[ind] + dind + ( pto < extra ? 1 : 0 );
1339 [ - + ][ # # ]: 18 : if( is_periodic && pfrom == np - 2 ) rdims[ind + 3]++; // neighbor is on periodic bdy
1340 : : }
1341 [ # # ][ # # ]: 0 : else if( -1 == dijk[ind] && !pfrom && gperiodic[ind] )
[ # # ]
1342 : : {
1343 : : // downward across periodic bdy
1344 : 0 : pto = np - 1;
1345 : 0 : facedims[ind + 3] = facedims[ind] = gdims[ind]; // by convention, facedims is within gdims, so lower value
1346 : 0 : rdims[ind + 3] =
1347 : 0 : gdims[ind + 3] + 1; // by convention, local dims one greater than gdims to indicate global lower value
1348 : 0 : rdims[ind] = rdims[ind + 3] - dind - 1;
1349 : 0 : across_bdy[ind] = -1;
1350 : : }
1351 [ # # ][ # # ]: 0 : else if( 1 == dijk[ind] && pfrom == np - 1 && is_periodic )
[ # # ]
1352 : : {
1353 : : // right across periodic bdy
1354 : 0 : pto = 0;
1355 : 0 : facedims[ind + 3] = facedims[ind] = gdims[ind]; // by convention, facedims is within gdims, so lowest value
1356 : 0 : rdims[ind] = gdims[ind];
1357 [ # # ]: 0 : rdims[ind + 3] = rdims[ind] + dind + ( pto < extra ? 1 : 0 );
1358 : 0 : across_bdy[ind] = 1;
1359 : : }
1360 : :
1361 [ + - ][ + - ]: 36 : assert( -1 == pto || ( rdims[0] >= gdims[0] && ( rdims[3] <= gdims[3] || ( across_bdy[0] && !pfrom ) ) ) );
[ - + ][ # # ]
[ # # ]
1362 [ + - ][ + - ]: 36 : assert( -1 == pto || ( rdims[1] >= gdims[1] && ( rdims[4] <= gdims[4] || ( across_bdy[1] && !pfrom ) ) ) );
[ - + ][ # # ]
[ # # ]
1363 [ + - ][ + - ]: 36 : assert( -1 == pto || ( rdims[2] >= gdims[2] && rdims[5] <= gdims[5] ) );
[ - + ]
1364 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[0] >= rdims[0] && facedims[3] <= rdims[3] ) );
[ - + ]
1365 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[1] >= rdims[1] && facedims[4] <= rdims[4] ) );
[ - + ]
1366 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[2] >= rdims[2] && facedims[5] <= rdims[5] ) );
[ - + ]
1367 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[0] >= ldims[0] && facedims[3] <= ldims[3] ) );
[ - + ]
1368 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[1] >= ldims[1] && facedims[4] <= ldims[4] ) );
[ - + ]
1369 [ + - ][ + - ]: 36 : assert( -1 == pto || ( facedims[2] >= ldims[2] && facedims[5] <= ldims[5] ) );
[ - + ]
1370 : :
1371 : 278 : return rval;
1372 : : }
1373 : :
1374 : : //! get shared vertices for alljorkori partition scheme
1375 : : #ifndef MOAB_HAVE_MPI
1376 : : ErrorCode ScdInterface::get_shared_vertices( ParallelComm*, ScdBox*, std::vector< int >&, std::vector< int >&,
1377 : : std::vector< int >& )
1378 : : {
1379 : : return MB_FAILURE;
1380 : : #else
1381 : 0 : ErrorCode ScdInterface::get_shared_vertices( ParallelComm* pcomm, ScdBox* box, std::vector< int >& procs,
1382 : : std::vector< int >& offsets, std::vector< int >& shared_indices )
1383 : : {
1384 : : // get index of partitioned dimension
1385 [ # # ]: 0 : const int* ldims = box->box_dims();
1386 : : ErrorCode rval;
1387 : : int ijkrem[6], ijkface[6], across_bdy[3];
1388 : :
1389 [ # # ]: 0 : for( int k = -1; k <= 1; k++ )
1390 : : {
1391 [ # # ]: 0 : for( int j = -1; j <= 1; j++ )
1392 : : {
1393 [ # # ]: 0 : for( int i = -1; i <= 1; i++ )
1394 : : {
1395 [ # # ][ # # ]: 0 : if( !i && !j && !k ) continue;
[ # # ]
1396 : : int pto;
1397 : 0 : int dijk[] = { i, j, k };
1398 [ # # ][ # # ]: 0 : rval = get_neighbor( pcomm->proc_config().proc_size(), pcomm->proc_config().proc_rank(),
[ # # ][ # # ]
1399 [ # # ][ # # ]: 0 : box->par_data(), dijk, pto, ijkrem, ijkface, across_bdy );
1400 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
1401 [ # # ]: 0 : if( -1 != pto )
1402 : : {
1403 [ # # ][ # # ]: 0 : if( procs.empty() || pto != *procs.rbegin() )
[ # # ][ # # ]
[ # # # # ]
1404 : : {
1405 [ # # ]: 0 : procs.push_back( pto );
1406 [ # # ]: 0 : offsets.push_back( shared_indices.size() );
1407 : : }
1408 [ # # ]: 0 : rval = get_indices( ldims, ijkrem, across_bdy, ijkface, shared_indices );
1409 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
1410 : :
1411 : : // check indices against known #verts on local and remote
1412 : : // begin of this block is shared_indices[*offsets.rbegin()], end is
1413 : : // shared_indices.end(), halfway is
1414 : : // (shared_indices.size()-*offsets.rbegin())/2
1415 : : #ifndef NDEBUG
1416 [ # # ]: 0 : int start_idx = *offsets.rbegin(), end_idx = shared_indices.size(),
1417 : 0 : mid_idx = ( start_idx + end_idx ) / 2;
1418 : :
1419 [ # # ]: 0 : int num_local_verts = ( ldims[3] - ldims[0] + 1 ) * ( ldims[4] - ldims[1] + 1 ) *
1420 [ # # ]: 0 : ( -1 == ldims[2] && -1 == ldims[5] ? 1 : ( ldims[5] - ldims[2] + 1 ) ),
1421 [ # # ]: 0 : num_remote_verts = ( ijkrem[3] - ijkrem[0] + 1 ) * ( ijkrem[4] - ijkrem[1] + 1 ) *
1422 [ # # ]: 0 : ( -1 == ijkrem[2] && -1 == ijkrem[5] ? 1 : ( ijkrem[5] - ijkrem[2] + 1 ) );
1423 : :
1424 [ # # ][ # # ]: 0 : assert(
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
1425 : : *std::min_element( &shared_indices[start_idx], &shared_indices[mid_idx] ) >= 0 &&
1426 : : *std::max_element( &shared_indices[start_idx], &shared_indices[mid_idx] ) < num_local_verts &&
1427 : : *std::min_element( &shared_indices[mid_idx], &shared_indices[end_idx] ) >= 0 &&
1428 [ # # ]: 0 : *std::max_element( &shared_indices[mid_idx], &shared_indices[end_idx] ) < num_remote_verts );
1429 : : #endif
1430 : : }
1431 : : }
1432 : : }
1433 : : }
1434 : :
1435 [ # # ]: 0 : offsets.push_back( shared_indices.size() );
1436 : :
1437 : 0 : return MB_SUCCESS;
1438 : : #endif
1439 : : }
1440 : :
1441 : 0 : std::ostream& operator<<( std::ostream& str, const ScdParData& pd )
1442 : : {
1443 : 0 : str << "Partition method = " << ScdParData::PartitionMethodNames[pd.partMethod] << ", gDims = (" << pd.gDims[0]
1444 : 0 : << "," << pd.gDims[1] << "," << pd.gDims[2] << ")-(" << pd.gDims[3] << "," << pd.gDims[4] << "," << pd.gDims[5]
1445 : 0 : << "), gPeriodic = (" << pd.gPeriodic[0] << "," << pd.gPeriodic[1] << "," << pd.gPeriodic[2] << "), pDims = ("
1446 : 0 : << pd.pDims[0] << "," << pd.pDims[1] << "," << pd.pDims[2] << ")" << std::endl;
1447 : 0 : return str;
1448 : : }
1449 : :
1450 [ + - ][ + - ]: 228 : } // namespace moab
|
