MOAB: Mesh Oriented datABase  (version 5.5.0)
scdtest.cpp
Go to the documentation of this file.
1 #include <string>
2 #include <iomanip>
3 #include <iostream>
4 #include <cassert>
5 
6 #include "moab/Core.hpp"
7 #include "moab/ParallelComm.hpp"
8 #include "moab/HomXform.hpp"
9 #include "MBParallelConventions.h"
10 #include "MBTagConventions.hpp"
11 
12 using namespace moab;
13 
14 // Number of cells in each direction:
15 // const int NC = 2;
16 
17 const int NI = 2;
18 const int NJ = 2;
19 const int NK = 1;
20 
21 // Number of processes:
22 // const int NPROCS = 4;
23 
24 // Domain size:
25 // const double DSIZE = 10.0;
26 
27 // MOAB objects:
28 Interface* mbint = NULL;
29 ParallelComm* mbpc = NULL;
30 
31 // Local domain starting and ending hex indexes:
32 int is, js, ks;
33 int ie, je, ke;
34 
35 // Obvious:
36 int rank;
37 int size;
38 
39 void set_local_domain_bounds();
40 void create_hexes_and_verts();
41 void resolve_and_exchange();
42 void error( ErrorCode err );
43 
44 int main( int argc, char* argv[] )
45 {
46  MPI_Init( &argc, &argv );
47  MPI_Comm_rank( MPI_COMM_WORLD, &rank );
48  MPI_Comm_size( MPI_COMM_WORLD, &size );
49  if( size != 4 && size != 2 )
50  {
51  std::cerr << "Run this with 2 or 4 processes\n";
52  MPI_Finalize();
53  exit( 1 );
54  }
55 
56  mbint = new Core();
57  mbpc = new ParallelComm( mbint, MPI_COMM_WORLD );
58 
59  set_local_domain_bounds();
60  create_hexes_and_verts();
61  resolve_and_exchange();
62 
63  delete mbpc; // ParallelComm instance should be deleted before MOAB instance is deleted
64  delete mbint;
65 
66  MPI_Finalize();
67  return 0;
68 }
69 
70 void set_local_domain_bounds()
71 {
72  switch( size )
73  {
74  case 2:
75  switch( rank )
76  {
77  case 0:
78  is = 0;
79  ie = NI / 2;
80  js = 0;
81  je = NJ;
82  ks = 0;
83  ke = NK;
84  break;
85  case 1:
86  is = NI / 2;
87  ie = NI;
88  js = 0;
89  je = NJ;
90  ks = 0;
91  ke = NK;
92  break;
93  }
94  break;
95 
96  case 4:
97  switch( rank )
98  {
99  case 0:
100  is = 0;
101  ie = NI / 2;
102  js = 0;
103  je = NJ / 2;
104  ks = 0;
105  ke = NK;
106  break;
107  case 1:
108  is = NI / 2;
109  ie = NI;
110  js = 0;
111  je = NJ / 2;
112  ks = 0;
113  ke = NK;
114  break;
115  case 2:
116  is = 0;
117  ie = NI / 2;
118  js = NJ / 2;
119  je = NJ;
120  ks = 0;
121  ke = NK;
122  break;
123  case 3:
124  is = NI / 2;
125  ie = NI;
126  js = NJ / 2;
127  je = NJ;
128  ks = 0;
129  ke = NK;
130  break;
131  }
132  break;
133 
134  default:
135  std::cerr << "Run this with 4 processes\n";
136  exit( 1 );
137  }
138 }
139 
140 void create_hexes_and_verts()
141 {
142  Core* mbcore = dynamic_cast< Core* >( mbint );
143  HomCoord coord_min( 0, 0, 0 );
144  HomCoord coord_max( ie - is, je - js, ke - ks );
145  EntitySequence* vertex_seq = NULL;
146  EntitySequence* cell_seq = NULL;
147  EntityHandle vs, cs;
148 
149  error( mbcore->create_scd_sequence( coord_min, coord_max, MBVERTEX, 1, vs, vertex_seq ) );
150  error( mbcore->create_scd_sequence( coord_min, coord_max, MBHEX, 1, cs, cell_seq ) );
151 
152  HomCoord p1( 0, 0, 0 );
153  HomCoord p2( 1, 0, 0 );
154  HomCoord p3( 0, 1, 0 );
155 
156  error( mbcore->add_vsequence( vertex_seq, cell_seq, p1, p1, p2, p2, p3, p3 ) );
157 
158  // Set global id's:
159  int gid;
160  Tag global_id_tag = mbint->globalId_tag();
161  EntityHandle handle = vs;
162  int i, j, k;
163 
164  ErrorCode err;
165 
166  for( k = ks; k < ke + 1; k++ )
167  for( j = js; j < je + 1; j++ )
168  for( i = is; i < ie + 1; i++ )
169  {
170  gid = 1 + i + j * ( NI + 1 ) + k * ( NI + 1 ) * ( NJ + 1 );
171  err = mbint->tag_set_data( global_id_tag, &handle, 1, &gid );
172  if( err != MB_SUCCESS )
173  {
174  exit( 1 );
175  }
176  handle++;
177  }
178 
179  handle = cs;
180  for( k = ks; k < ke; k++ )
181  for( j = js; j < je; j++ )
182  for( i = is; i < ie; i++ )
183  {
184  gid = 1 + i + j * NI + k * NI * NJ;
185  error( mbint->tag_set_data( global_id_tag, &handle, 1, &gid ) );
186  handle++;
187  }
188 }
189 
190 void resolve_and_exchange()
191 {
192  EntityHandle entity_set;
193 
194  // Create the entity set:
195  error( mbint->create_meshset( MESHSET_SET, entity_set ) );
196 
197  // Get a list of hexes:
198  Range range;
199  error( mbint->get_entities_by_type( 0, MBHEX, range ) );
200 
201  // Add entities to the entity set:
202  error( mbint->add_entities( entity_set, range ) );
203 
204  // Add the MATERIAL_SET tag to the entity set:
205  Tag tag;
206  error( mbint->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, tag ) );
207  error( mbint->tag_set_data( tag, &entity_set, 1, &rank ) );
208 
209  // Set up partition sets. This is where MOAB is actually told what
210  // entities each process owns:
211  error( mbint->get_entities_by_type_and_tag( 0, MBENTITYSET, &tag, NULL, 1, mbpc->partition_sets() ) );
212 
213  // Finally, determine which entites are shared and exchange the
214  // ghosted entities:
215  error( mbpc->resolve_shared_ents( 0, -1, -1 ) );
216  error( mbpc->exchange_ghost_cells( -1, 0, 1, 0, true ) );
217 }
218 
219 void error( ErrorCode err )
220 {
221  if( err != MB_SUCCESS )
222  {
223  std::cerr << "Error: MOAB function failed\n";
224  assert( 0 );
225  }
226 }