 |
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
|
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
1 /*! \file NestedRefine.hpp
2 * This class implements the generation of a hierarchy of meshes via uniform refinement from an
3 * input mesh. The internal upper bound on the number of levels is set to 20.
4 */
5
6 #ifndef NESTED_REFINE_HPP
7 #define NESTED_REFINE_HPP
8
9 #include "moab/MOABConfig.h"
10 #include "moab/Range.hpp"
11 #include "moab/CN.hpp"
12 #include <map>
13 #include <set>
14
15 namespace moab
16 {
17
18 #define MAX_DEGREE 3
19 #define MAX_VERTS 64
20 #define MAX_CHILDRENS 27
21 #define MAX_HE 12
22 #define MAX_HF 6
23 #define MAX_CONN 8
24 #define MAX_VHF 20
25 #define MAX_LEVELS 20
26 #define MAX_PROCS 64
27
28 // Forward Declarations
29 class Core;
30 class HalfFacetRep;
31 class ParallelComm;
32 class CpuTimer;
33
34 class NestedRefine
35 {
36
37 public:
38 NestedRefine( Core* impl, ParallelComm* comm = 0, EntityHandle rset = 0 );
39
40 ~NestedRefine();
41
42 ErrorCode initialize();
43
44 //! \brief Generate a mesh hierarchy.
45 /** Given a mesh, generate a sequence of meshes via uniform refinement. The inputs are: a) an
46 * array(level_degrees) storing the degrees which will be used to refine the previous level mesh
47 * to generate a new level and b) the total number of levels(should be same length as that of
48 * the array in a). Each mesh level in the hierarchy are stored in different meshsets whose
49 * handles are returned after the hierarchy generation. These handles can be used to work with a
50 * specific mesh level. \param level_degrees Integer array storing the degrees used in each
51 * level. \param num_level The total number of levels in the hierarchy. \param hm_set
52 * EntityHandle STL vector that returns the handles of the sets created for each mesh level.
53 */
54
55 ErrorCode generate_mesh_hierarchy( int num_level,
56 int* level_degrees,
57 std::vector< EntityHandle >& level_sets,
58 bool optimize = false );
59
60 //! Given an entity and its level, return its connectivity.
61 /** Given an entity at a certain level, it finds the connectivity via direct access to a stored
62 * internal pointer to the memory to connectivity sequence for the given level. \param ent
63 * EntityHandle of the entity \param level Integer level of the entity for which connectivity is
64 * requested \param conn std::vector returning the connectivity of the entity
65 */
66
67 ErrorCode get_connectivity( EntityHandle ent, int level, std::vector< EntityHandle >& conn );
68
69 //! Given a vector of vertices and their level, return its coordinates.
70 /** Given a vector of vertices at a certain level, it finds the coordinates via direct access to
71 * a stored internal pointer to the memory to coordinate sequence for the given level. \param
72 * verts std::vector of the entity handles of the vertices \param num_verts The number of
73 * vertices \param level Integer level of the entity for which connectivity is requested \param
74 * coords double pointer returning the coordinates of the vertices
75 */
76
77 ErrorCode get_coordinates( EntityHandle* verts, int num_verts, int level, double* coords );
78
79 //! Get the adjacencies associated with an entity.
80 /** Given an entity of dimension <em>d</em>, gather all the adjacent <em>D</em> dimensional
81 * entities where <em>D >, = , < d </em>.
82 *
83 * \param source_entity EntityHandle to which adjacent entities have to be found.
84 * \param target_dimension Int Dimension of the desired adjacent entities.
85 * \param target_entities Vector in which the adjacent EntityHandle are returned.
86 */
87
88 ErrorCode get_adjacencies( const EntityHandle source_entity,
89 const unsigned int target_dimension,
90 std::vector< EntityHandle >& target_entities );
91
92 // Interlevel parent-child or vice-versa queries
93 /** Given an entity from a certain level, it returns a pointer to its parent at the requested
94 * parent level. NOTE: This query does not support vertices.
95 *
96 * \param child EntityHandle of the entity whose parent is requested
97 * \param child_level Mesh level where the child exists
98 * \param parent_level Mesh level from which parent is requested
99 * \param parent Pointer to the parent in the requested parent_level
100 */
101
102 ErrorCode child_to_parent( EntityHandle child, int child_level, int parent_level, EntityHandle* parent );
103
104 /** Given an entity from a certain level, it returns a std::vector of all its children from the
105 * requested child level. NOTE: This query does not support vertices.
106 *
107 * \param parent EntityHandle of the entity whose children in subsequent level is requested
108 * \param parent_level Mesh level where the parent exists
109 * \param child_level Mesh level from which its children are requested
110 * \param children Vector containing all childrens from the requested child_level
111 */
112
113 ErrorCode parent_to_child( EntityHandle parent,
114 int parent_level,
115 int child_level,
116 std::vector< EntityHandle >& children );
117
118 /** Given a vertex from a certain level, it returns a std::vector of all entities from any
119 * previous levels that contains it.
120 *
121 * \param vertex EntityHandle of the vertex
122 * \param vert_level Mesh level of the vertex
123 * \param parent_level Mesh level from which entities containing vertex is requested
124 * \param incident_entities Vector containing entities from the parent level incident on the
125 * vertex
126 */
127
128 ErrorCode vertex_to_entities_up( EntityHandle vertex,
129 int vert_level,
130 int parent_level,
131 std::vector< EntityHandle >& incident_entities );
132
133 /** Given a vertex from a certain level, it returns a std::vector of all children entities of
134 * incident entities to vertex from any subsequent levels
135 *
136 * \param vertex EntityHandle of the vertex
137 * \param vert_level Mesh level of the vertex
138 * \param child_level Mesh level from which child entities are requested
139 * \param incident_entities Vector containing entities from the child level
140 */
141
142 ErrorCode vertex_to_entities_down( EntityHandle vertex,
143 int vert_level,
144 int child_level,
145 std::vector< EntityHandle >& incident_entities );
146
147 ErrorCode get_vertex_duplicates( EntityHandle vertex, int level, EntityHandle& dupvertex );
148
149 /** Given an entity at a certain level, it returns a boolean value true if it lies on the domain
150 * boundary. \param entity
151 */
152
153 bool is_entity_on_boundary( const EntityHandle& entity );
154
155 ErrorCode exchange_ghosts( std::vector< EntityHandle >& lsets, int num_glayers );
156
157 ErrorCode update_special_tags( int level, EntityHandle& lset );
158
159 struct codeperf
160 {
161 double tm_total;
162 double tm_refine;
163 double tm_resolve;
164 };
165
166 codeperf timeall;
167
168 protected:
169 Core* mbImpl;
170 ParallelComm* pcomm;
171 HalfFacetRep* ahf;
172 CpuTimer* tm;
173 EntityHandle _rset;
174
175 Range _inverts, _inedges, _infaces, _incells;
176
177 EntityType elementype;
178 int meshdim, nlevels;
179 int level_dsequence[MAX_LEVELS];
180 std::map< int, int > deg_index;
181 bool hasghost;
182
183 /*! \struct refPatterns
184 * Refinement patterns w.r.t the reference element. It consists of a locally indexed vertex list
185 * along with their natural coordinates, the connectivity of the subdivided entities with local
186 * indices, their local AHF maps along with other helper fields to aid in general book keeping
187 * such as avoiding vertex duplication during refinement. The entity and degree specific values
188 * are stored in the Templates.hpp. \sa Templates.hpp
189 */
190
191 //! refPatterns
192 struct refPatterns
193 {
194 //! Number of new vertices on edge
195 short int nv_edge;
196 //! Number of new vertices on face, does not include those on edge
197 short int nv_face;
198 //! Number of new vertices in cell
199 short int nv_cell;
200 //! Total number of new vertices per entity
201 short int total_new_verts;
202 //! Total number of new child entities
203 short int total_new_ents;
204 //! Lower and upper indices of the new vertices
205 int vert_index_bnds[2];
206 //! Natural coordinates of the new vertices w.r.t reference
207 double vert_nat_coord[MAX_VERTS][3];
208 //! Connectivity of the new entities
209 int ents_conn[MAX_CHILDRENS][MAX_CONN];
210 //! Vertex to half-facet map of the new vertices
211 int v2hf[MAX_VERTS][2];
212 //! Opposite half-facet map of the new entities
213 int ents_opphfs[MAX_CHILDRENS][2 * MAX_CONN];
214 //! Helper: storing the local ids of vertices on each local edge
215 int vert_on_edges[MAX_HE][MAX_VHF];
216 //! Helper: storing local ids of verts on each local face, doesnt include those on edges of
217 //! the face
218 int vert_on_faces[MAX_HF][MAX_VHF];
219 //! Helper: stores child half-facets incident on parent half-facet. First column contain the
220 //! number of such children
221 int ents_on_pent[MAX_HF][MAX_CHILDRENS];
222 //! Helper: stores child ents incident on new verts on edge.
223 // Each triad in the column consists of :
224 // 1) a local child incident on the vertex on the edge
225 // 2) the local face id from the child
226 // 3) the local vertex id wrt the child connectivity
227 int ents_on_vedge[MAX_HE][MAX_VHF * 3];
228 };
229 //! refPatterns
230
231 static const refPatterns refTemplates[9][MAX_DEGREE];
232
233 //! Helper
234 struct intFEdge
235 {
236 //! Number of edges interior to a face
237 short int nie;
238 //! Local connectivity of the interior edges
239 int ieconn[12][2];
240 };
241 static const intFEdge intFacEdg[2][2];
242
243 int get_index_from_degree( int degree );
244
245 // HM Storage Helper
246 struct level_memory
247 {
248 int num_verts, num_edges, num_faces, num_cells;
249 EntityHandle start_vertex, start_edge, start_face, start_cell;
250 std::vector< double* > coordinates;
251 EntityHandle *edge_conn, *face_conn, *cell_conn;
252 Range verts, edges, faces, cells;
253 };
254
255 level_memory level_mesh[MAX_LEVELS];
256
257 // Basic Functions
258
259 // Estimate and create storage for the levels
260 ErrorCode estimate_hm_storage( EntityHandle set, int level_degree, int cur_level, int hmest[4] );
261 ErrorCode create_hm_storage_single_level( EntityHandle* set, int cur_level, int estL[4] );
262
263 // Generate HM : Construct the hierarchical mesh: 1D, 2D, 3D
264 ErrorCode generate_hm( int* level_degrees, int num_level, EntityHandle* hm_set, bool optimize );
265 ErrorCode construct_hm_entities( int cur_level, int deg );
266 ErrorCode construct_hm_1D( int cur_level, int deg );
267 ErrorCode construct_hm_1D( int cur_level, int deg, EntityType type, std::vector< EntityHandle >& trackverts );
268 ErrorCode construct_hm_2D( int cur_level, int deg );
269 ErrorCode construct_hm_2D( int cur_level,
270 int deg,
271 EntityType type,
272 std::vector< EntityHandle >& trackvertsC_edg,
273 std::vector< EntityHandle >& trackvertsF );
274 ErrorCode construct_hm_3D( int cur_level, int deg );
275
276 ErrorCode subdivide_cells( EntityType type, int cur_level, int deg );
277 ErrorCode subdivide_tets( int cur_level, int deg );
278
279 // General helper functions
280 ErrorCode copy_vertices_from_prev_level( int cur_level );
281 ErrorCode count_subentities( EntityHandle set, int cur_level, int* nedges, int* nfaces );
282 ErrorCode get_octahedron_corner_coords( int cur_level, int deg, EntityHandle* vbuffer, double* ocoords );
283 int find_shortest_diagonal_octahedron( int cur_level, int deg, EntityHandle* vbuffer );
284 int get_local_vid( EntityHandle vid, EntityHandle ent, int level );
285
286 // Book-keeping functions
287 ErrorCode update_tracking_verts( EntityHandle cid,
288 int cur_level,
289 int deg,
290 std::vector< EntityHandle >& trackvertsC_edg,
291 std::vector< EntityHandle >& trackvertsC_face,
292 EntityHandle* vbuffer );
293 ErrorCode reorder_indices( int cur_level,
294 int deg,
295 EntityHandle cell,
296 int lfid,
297 EntityHandle sib_cell,
298 int sib_lfid,
299 int index,
300 int* id_sib );
301 ErrorCode reorder_indices( int deg,
302 EntityHandle* face1_conn,
303 EntityHandle* face2_conn,
304 int nvF,
305 std::vector< int >& lemap,
306 std::vector< int >& vidx,
307 int* leorient = NULL );
308 ErrorCode reorder_indices( int deg, int nvF, int comb, int* childfid_map );
309 ErrorCode reorder_indices( EntityHandle* face1_conn,
310 EntityHandle* face2_conn,
311 int nvF,
312 int* conn_map,
313 int& comb,
314 int* orient = NULL );
315 ErrorCode get_lid_inci_child( EntityType type,
316 int deg,
317 int lfid,
318 int leid,
319 std::vector< int >& child_ids,
320 std::vector< int >& child_lvids );
321
322 // Permutation matrices
323 struct pmat
324 {
325 short int num_comb; // Number of combinations
326 int comb[MAX_HE][MAX_HE]; // Combinations
327 int lemap[MAX_HE][MAX_HE]; // Local edge map
328 int orient[MAX_HE]; // Orientation
329 int porder2[MAX_HE][MAX_HE]; // Permuted order degree 2
330 int porder3[MAX_HE][MAX_HE]; // Permuted order degree 3
331 };
332
333 static const pmat permutation[2];
334
335 // Print functions
336 ErrorCode print_maps_1D( int level );
337 ErrorCode print_maps_2D( int level, EntityType type );
338 ErrorCode print_maps_3D( int level, EntityType type );
339
340 // Coordinates
341 ErrorCode compute_coordinates( int cur_level,
342 int deg,
343 EntityType type,
344 EntityHandle* vbuffer,
345 int vtotal,
346 double* corner_coords,
347 std::vector< int >& vflag,
348 int nverts_prev );
349
350 // Update the ahf maps
351
352 ErrorCode update_local_ahf( int deg, EntityType type, EntityHandle* vbuffer, EntityHandle* ent_buffer, int etotal );
353
354 ErrorCode update_local_ahf( int deg,
355 EntityType type,
356 int pat_id,
357 EntityHandle* vbuffer,
358 EntityHandle* ent_buffer,
359 int etotal );
360
361 // ErrorCode update_global_ahf(EntityType type, int cur_level, int deg);
362
363 // ErrorCode update_global_ahf(int cur_level, int deg, std::vector<int> &pattern_ids);
364
365 ErrorCode update_global_ahf( EntityType type, int cur_level, int deg, std::vector< int >* pattern_ids = NULL );
366
367 ErrorCode update_global_ahf_1D( int cur_level, int deg );
368
369 ErrorCode update_global_ahf_1D_sub( int cur_level, int deg );
370
371 ErrorCode update_ahf_1D( int cur_level );
372
373 ErrorCode update_global_ahf_2D( int cur_level, int deg );
374
375 ErrorCode update_global_ahf_2D_sub( int cur_level, int deg );
376
377 ErrorCode update_global_ahf_3D( int cur_level, int deg, std::vector< int >* pattern_ids = NULL );
378
379 // ErrorCode update_global_ahf_3D(int cur_level, int deg);
380
381 // ErrorCode update_global_ahf_3D(int cur_level, int deg, std::vector<int> &pattern_ids);
382
383 /** Boundary extraction functions
384 * Given a vertex at a certain level, it returns a boolean value true if it lies on the domain
385 * boundary. Note: This is a specialization of the NestedRefine::is_entity_on_boundary function
386 * and applies only to vertex queries. \param entity
387 */
388 bool is_vertex_on_boundary( const EntityHandle& entity );
389 bool is_edge_on_boundary( const EntityHandle& entity );
390 bool is_face_on_boundary( const EntityHandle& entity );
391 bool is_cell_on_boundary( const EntityHandle& entity );
392
393 // ErrorCode find_skin_faces(EntityHandle set, int level, int nskinF);
394
395 /** Parallel communication routines
396 * We implement two strategies to resolve the shared entities of the newly created entities.
397 * The first strategy is to use the existing parallel merge capability which essentially uses
398 * a coordinate-based matching of vertices and subsequently the entity handles through
399 * their connectivities. The second strategy is an optimized and a new algorithm. It uses
400 * the existing shared information from the coarse entities and propagates the parallel
401 * information appropriately.
402 */
403
404 // Send/Recv Buffer storage
405 /* struct pbuffer{
406 int rank;
407 std::vector<int> msgsize;
408 std::vector<EntityHandle> localBuff;
409 std::vector<EntityHandle> remoteBuff;
410 };
411
412 pbuffer commBuffers[MAX_PROCS];
413
414 //Parallel tag values
415 struct parinfo{
416 std::multimap<EntityHandle, int> remoteps;
417 std::multimap<EntityHandle, EntityHandle> remotehs;
418 };
419
420 parinfo parallelInfo[MAX_LEVELS];*/
421
422 #ifdef MOAB_HAVE_MPI
423
424 ErrorCode resolve_shared_ents_parmerge( int level, EntityHandle levelset );
425 ErrorCode resolve_shared_ents_opt( EntityHandle* hm_set, int num_levels );
426
427 ErrorCode collect_shared_entities_by_dimension( Range sharedEnts, Range& allEnts );
428 ErrorCode collect_FList( int to_proc,
429 Range faces,
430 std::vector< EntityHandle >& FList,
431 std::vector< EntityHandle >& RList );
432 ErrorCode collect_EList( int to_proc,
433 Range edges,
434 std::vector< EntityHandle >& EList,
435 std::vector< EntityHandle >& RList );
436 ErrorCode collect_VList( int to_proc,
437 Range verts,
438 std::vector< EntityHandle >& VList,
439 std::vector< EntityHandle >& RList );
440
441 ErrorCode decipher_remote_handles( std::vector< int >& sharedprocs,
442 std::vector< std::vector< int > >& auxinfo,
443 std::vector< std::vector< EntityHandle > >& localbuffers,
444 std::vector< std::vector< EntityHandle > >& remotebuffers,
445 std::multimap< EntityHandle, int >& remProcs,
446 std::multimap< EntityHandle, EntityHandle >& remHandles );
447
448 ErrorCode decipher_remote_handles_face( int shared_proc,
449 int numfaces,
450 std::vector< EntityHandle >& localFaceList,
451 std::vector< EntityHandle >& remFaceList,
452 std::multimap< EntityHandle, int >& remProcs,
453 std::multimap< EntityHandle, EntityHandle >& remHandles );
454
455 ErrorCode decipher_remote_handles_edge( int shared_proc,
456 int numedges,
457 std::vector< EntityHandle >& localEdgeList,
458 std::vector< EntityHandle >& remEdgeList,
459 std::multimap< EntityHandle, int >& remProcs,
460 std::multimap< EntityHandle, EntityHandle >& remHandles );
461
462 ErrorCode decipher_remote_handles_vertex( int shared_proc,
463 int numverts,
464 std::vector< EntityHandle >& localVertexList,
465 std::vector< EntityHandle >& remVertexList,
466 std::multimap< EntityHandle, int >& remProcs,
467 std::multimap< EntityHandle, EntityHandle >& remHandles );
468
469 ErrorCode update_parallel_tags( std::multimap< EntityHandle, int >& remProcs,
470 std::multimap< EntityHandle, EntityHandle >& remHandles );
471
472 ErrorCode get_data_from_buff( int dim,
473 int type,
474 int level,
475 int entityidx,
476 int nentities,
477 std::vector< EntityHandle >& buffer,
478 std::vector< EntityHandle >& data );
479
480 bool check_for_parallelinfo( EntityHandle entity, int proc, std::multimap< EntityHandle, int >& remProcs );
481
482 ErrorCode check_for_parallelinfo( EntityHandle entity,
483 int proc,
484 std::multimap< EntityHandle, EntityHandle >& remHandles,
485 std::multimap< EntityHandle, int >& remProcs,
486 EntityHandle& rhandle );
487
488 #endif
489 };
490
491 } // namespace moab
492 #endif
1.9.1.