 |
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
|
#include <iostream>#include <vector>#include <string>#include <iomanip>#include <fstream>#include <ctime>#include <cmath>#include <cassert>#include <cfloat>#include "moab/Core.hpp"#include "moab/Interface.hpp"#include "moab/verdict/VerdictWrapper.hpp"#include "moab/NestedRefine.hpp"
Include dependency graph for umr.cpp:Go to the source code of this file.
Macros | |
| #define | SUCCESS 0 |
| #define | USAGE_ERROR 1 |
| #define | NOT_IMPLEMENTED 2 |
Functions | |
| static void | print_usage (const char *name, std::ostream &stream) |
| static void | usage_error (const char *name) |
| bool | parse_id_list (const char *string, int dim, int nval, std::vector< int > &results) |
| bool | make_opts_string (std::vector< std::string > options, std::string &opts) |
| ErrorCode | get_degree_seq (Core &mb, EntityHandle fileset, int dim, double desired_vol, int &num_levels, std::vector< int > &level_degs) |
| ErrorCode | get_max_volume (Core &mb, EntityHandle fileset, int dim, double &vmax) |
| int | main (int argc, char *argv[]) |
| ErrorCode get_degree_seq | ( | Core & | mb, |
| EntityHandle | fileset, | ||
| int | dim, | ||
| double | desired_vol, | ||
| int & | num_levels, | ||
| std::vector< int > & | level_degs | ||
| ) |
Definition at line 411 of file umr.cpp.
417 {
418 // Find max volume
419 double vmax_global;
420 ErrorCode error = get_max_volume( mb, fileset, dim, vmax_global );MB_CHK_ERR( error );
421
422 int init_nl = num_levels;
423 num_levels = 0;
424 level_degs.clear();
425
426 // Find a sequence that reduces the maximum volume to desired.
427 // desired_vol should be a fraction or absolute value ?
428
429 // double remV = vmax_global*desired_vol/dim;
430 double Vdesired = desired_vol;
431 double try_x;
432 double remV = vmax_global;
433 int degs[3][3] = { { 5, 3, 2 }, { 25, 9, 4 }, { 0, 27, 8 } };
434
435 if( dim == 1 || dim == 2 )
436 {
437 while( remV - Vdesired >= 0 )
438 {
439 try_x = degs[dim - 1][0];
440 if( ( remV / try_x - Vdesired ) >= 0 )
441 {
442 level_degs.push_back( 5 );
443 num_levels += 1;
444 remV /= try_x;
445 }
446 else
447 {
448 try_x = degs[dim - 1][1];
449 if( ( remV / try_x - Vdesired ) >= 0 )
450 {
451 level_degs.push_back( 3 );
452 num_levels += 1;
453 remV /= try_x;
454 }
455 else
456 {
457 try_x = degs[dim - 1][2];
458 if( ( remV / try_x - Vdesired ) >= 0 )
459 {
460 level_degs.push_back( 2 );
461 num_levels += 1;
462 remV /= try_x;
463 }
464 else
465 break;
466 }
467 }
468 }
469 }
470 else
471 {
472 while( remV - Vdesired >= 0 )
473 {
474 try_x = degs[dim - 1][1];
475 if( ( remV / try_x - Vdesired ) >= 0 )
476 {
477 level_degs.push_back( 3 );
478 num_levels += 1;
479 remV /= try_x;
480 }
481 else
482 {
483 try_x = degs[dim - 1][2];
484 if( ( remV / try_x - Vdesired ) >= 0 )
485 {
486 level_degs.push_back( 2 );
487 num_levels += 1;
488 remV /= try_x;
489 }
490 else
491 break;
492 }
493 }
494 }
495
496 if( init_nl != 0 && init_nl < num_levels )
497 {
498 for( int i = level_degs.size(); i >= init_nl; i-- )
499 level_degs.pop_back();
500 num_levels = init_nl;
501 }
502
503 return MB_SUCCESS;
504 }
References dim, moab::error(), ErrorCode, get_max_volume(), mb, MB_CHK_ERR, and MB_SUCCESS.
Referenced by main().
| ErrorCode get_max_volume | ( | Core & | mb, |
| EntityHandle | fileset, | ||
| int | dim, | ||
| double & | vmax | ||
| ) |
Definition at line 506 of file umr.cpp.
507 {
508 ErrorCode error;
509 VerdictWrapper vw( &mb );
510 QualityType q;
511
512 switch( dim )
513 {
514 case 1:
515 q = MB_LENGTH;
516 break;
517 case 2:
518 q = MB_AREA;
519 break;
520 case 3:
521 q = MB_VOLUME;
522 break;
523 default:
524 return MB_FAILURE;
525 break;
526 }
527
528 // Get all entities of the highest dimension which is passed as a command line argument.
529 Range allents, owned;
530 error = mb.get_entities_by_handle( fileset, allents );MB_CHK_ERR( error );
531 owned = allents.subset_by_dimension( dim );MB_CHK_ERR( error );
532
533 // Get all owned entities
534 #ifdef MOAB_HAVE_MPI
535 int size = 1;
536 MPI_Comm_size( MPI_COMM_WORLD, &size );
537 int mpi_err;
538 if( size > 1 )
539 {
540 // filter the entities not owned, so we do not process them more than once
541 ParallelComm* pcomm = moab::ParallelComm::get_pcomm( &mb, 0 );
542 Range current = owned;
543 owned.clear();
544 error = pcomm->filter_pstatus( current, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1, &owned );
545 if( error != MB_SUCCESS )
546 {
547 MPI_Finalize();
548 return MB_FAILURE;
549 }
550 }
551 #endif
552
553 double vmax_local = 0;
554 // Find the maximum volume of an entity in the owned mesh
555 for( Range::iterator it = owned.begin(); it != owned.end(); it++ )
556 {
557 double volume;
558 error = vw.quality_measure( *it, q, volume );MB_CHK_ERR( error );
559 if( volume > vmax_local ) vmax_local = volume;
560 }
561
562 // Get the global maximum
563 double vmax_global = vmax_local;
564 #ifdef MOAB_HAVE_MPI
565 mpi_err = MPI_Reduce( &vmax_local, &vmax_global, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD );
566 if( mpi_err )
567 {
568 MPI_Finalize();
569 return MB_FAILURE;
570 }
571 #endif
572
573 vmax = vmax_global;
574
575 return MB_SUCCESS;
576 }
References moab::Range::begin(), moab::Range::clear(), dim, moab::Range::end(), moab::error(), ErrorCode, moab::ParallelComm::filter_pstatus(), moab::Core::get_entities_by_handle(), moab::ParallelComm::get_pcomm(), mb, moab::MB_AREA, MB_CHK_ERR, moab::MB_LENGTH, MB_SUCCESS, moab::MB_VOLUME, PSTATUS_NOT, PSTATUS_NOT_OWNED, moab::VerdictWrapper::quality_measure(), size, and moab::Range::subset_by_dimension().
Referenced by get_degree_seq(), and main().
| int main | ( | int | argc, |
| char * | argv[] | ||
| ) |
Definition at line 89 of file umr.cpp.
90 {
91 int proc_id = 0, size = 1;
92 #ifdef MOAB_HAVE_MPI
93 MPI_Init( &argc, &argv );
94 MPI_Comm_rank( MPI_COMM_WORLD, &proc_id );
95 MPI_Comm_size( MPI_COMM_WORLD, &size );
96 #endif
97
98 int num_levels = 0, dim = 0;
99 std::vector< int > level_degrees;
100 bool optimize = false;
101 bool do_flag = true;
102 bool print_times = false, print_size = false, output = false;
103 bool parallel = false, resolve_shared = false, exchange_ghosts = false;
104 bool printusage = false, parselevels = true;
105 bool qc_vol = false, only_quality = false;
106 double cvol = 0;
107 std::string infile;
108
109 int i = 1;
110 while( i < argc )
111 {
112 if( !argv[i][0] && 0 == proc_id )
113 {
114 usage_error( argv[0] );
115 #ifdef MOAB_HAVE_MPI
116 MPI_Finalize();
117 #endif
118 }
119
120 if( do_flag && argv[i][0] == '-' )
121 {
122 switch( argv[i][1] )
123 {
124 // do flag arguments:
125 case '-':
126 do_flag = false;
127 break;
128 case 't':
129 print_times = true;
130 break;
131 case 's':
132 print_size = true;
133 break;
134 case 'h':
135 case 'H':
136 print_usage( argv[0], std::cerr );
137 printusage = true;
138 break;
139 case 'd':
140 dim = atol( argv[i + 1] );
141 ++i;
142 break;
143 case 'n':
144 num_levels = atol( argv[i + 1] );
145 ++i;
146 break;
147 case 'L':
148 if( dim != 0 && num_levels != 0 )
149 {
150 parselevels = parse_id_list( argv[i + 1], dim, num_levels, level_degrees );
151 ++i;
152 }
153 else
154 {
155 print_usage( argv[0], std::cerr );
156 printusage = true;
157 }
158 break;
159 case 'V':
160 qc_vol = true;
161 cvol = strtod( argv[i + 1], NULL );
162 ++i;
163 break;
164 case 'q':
165 only_quality = true;
166 break;
167 case 'o':
168 output = true;
169 break;
170 case 'O':
171 optimize = true;
172 break;
173 #ifdef MOAB_HAVE_MPI
174 case 'p':
175 parallel = true;
176 resolve_shared = true;
177 if( argv[i][1] == '1' ) exchange_ghosts = true;
178 break;
179 #endif
180 default:
181 ++i;
182 switch( argv[i - 1][1] )
183 {
184 // case 'O': read_opts.push_back(argv[i]); break;
185 default:
186 std::cerr << "Invalid option: " << argv[i] << std::endl;
187 }
188 }
189 ++i;
190 }
191 // do file names
192 else
193 {
194 infile = argv[i];
195 ++i;
196 }
197 }
198
199 if( infile.empty() && !printusage ) print_usage( argv[0], std::cerr );
200
201 if( !parselevels ) exit( USAGE_ERROR );
202
203 #ifdef MOAB_HAVE_MPI
204 parallel = true;
205 resolve_shared = true;
206 #endif
207
208 ErrorCode error;
209
210 Core* moab = new Core;
211 Interface* mb = (Interface*)moab;
212 EntityHandle fileset;
213
214 // Create a fileset
215 error = mb->create_meshset( MESHSET_SET, fileset );MB_CHK_ERR( error );
216
217 // Set the read options for parallel file loading
218 std::vector< std::string > read_opts, write_opts;
219 std::string read_options, write_options;
220
221 if( parallel && size > 1 )
222 {
223 read_opts.push_back( "PARALLEL=READ_PART" );
224 read_opts.push_back( "PARTITION=PARALLEL_PARTITION" );
225 if( resolve_shared ) read_opts.push_back( "PARALLEL_RESOLVE_SHARED_ENTS" );
226 if( exchange_ghosts ) read_opts.push_back( "PARALLEL_GHOSTS=3.0.1" );
227
228 /* if (output)
229 {
230 write_opts.push_back(";;PARALLEL=WRITE_PART");
231 std::cout<<"Here"<<std::endl;
232 }*/
233 }
234
235 if( !make_opts_string( read_opts, read_options ) || !make_opts_string( write_opts, write_options ) )
236 {
237 #ifdef MOAB_HAVE_MPI
238 MPI_Finalize();
239 #endif
240 return USAGE_ERROR;
241 }
242
243 // Load file
244 std::cout << "READ OPTS=" << (char*)read_options.c_str() << std::endl;
245 error = mb->load_file( infile.c_str(), &fileset, read_options.c_str() );MB_CHK_ERR( error );
246
247 // Create the nestedrefine instance
248
249 #ifdef MOAB_HAVE_MPI
250 ParallelComm* pc = new ParallelComm( mb, MPI_COMM_WORLD );
251 NestedRefine* uref = new NestedRefine( moab, pc, fileset );
252 #else
253 NestedRefine* uref = new NestedRefine( moab );
254 #endif
255
256 std::vector< EntityHandle > lsets;
257
258 // std::cout<<"Read input file"<<std::endl;
259
260 if( only_quality )
261 {
262 double vmax;
263 error = get_max_volume( *moab, fileset, dim, vmax );MB_CHK_ERR( error );
264 #ifdef MOAB_HAVE_MPI
265 int rank = 0;
266 MPI_Comm_rank( MPI_COMM_WORLD, &rank );
267 if( rank == 0 ) std::cout << "Maximum global volume = " << vmax << std::endl;
268 #else
269 std::cout << "Maximum volume = " << vmax << std::endl;
270 #endif
271 exit( SUCCESS );
272 }
273
274 // If a volume constraint is given, find an optimal degree sequence to reach the desired volume
275 // constraint.
276 if( qc_vol )
277 {
278 error = get_degree_seq( *moab, fileset, dim, cvol, num_levels, level_degrees );MB_CHK_ERR( error );
279
280 if( dim == 0 ) print_usage( argv[0], std::cerr );
281 }
282
283 if( num_levels == 0 ) num_levels = 1;
284
285 int* ldeg = new int[num_levels];
286 // std::cout<<"level_degrees.size = "<<level_degrees.size()<<std::endl;
287 if( level_degrees.empty() )
288 {
289 for( int l = 0; l < num_levels; l++ )
290 ldeg[l] = 2;
291 }
292 else
293 {
294 for( int l = 0; l < num_levels; l++ )
295 ldeg[l] = level_degrees[l];
296 }
297
298 std::cout << "Starting hierarchy generation" << std::endl;
299
300 std::cout << "opt = " << optimize << std::endl;
301
302 error = uref->generate_mesh_hierarchy( num_levels, ldeg, lsets, optimize );MB_CHK_ERR( error );
303
304 if( print_times )
305 {
306 std::cout << "Finished hierarchy generation in " << uref->timeall.tm_total << " secs" << std::endl;
307 if( parallel )
308 {
309 std::cout << "Time taken for refinement " << uref->timeall.tm_refine << " secs" << std::endl;
310 std::cout << "Time taken for resolving shared interface " << uref->timeall.tm_resolve << " secs"
311 << std::endl;
312 }
313 }
314 else
315 std::cout << "Finished hierarchy generation " << std::endl;
316
317 if( print_size )
318 {
319 Range all_ents, ents[4];
320 error = mb->get_entities_by_handle( fileset, all_ents );MB_CHK_ERR( error );
321
322 for( int k = 0; k < 4; k++ )
323 ents[k] = all_ents.subset_by_dimension( k );
324
325 std::cout << std::endl;
326 if( qc_vol )
327 {
328 double volume;
329 error = get_max_volume( *moab, fileset, dim, volume );MB_CHK_ERR( error );
330 std::cout << "Mesh size for level 0"
331 << " :: nverts = " << ents[0].size() << ", nedges = " << ents[1].size()
332 << ", nfaces = " << ents[2].size() << ", ncells = " << ents[3].size() << " :: Vmax = " << volume
333 << std::endl;
334 }
335 else
336 std::cout << "Mesh size for level 0"
337 << " :: nverts = " << ents[0].size() << ", nedges = " << ents[1].size()
338 << ", nfaces = " << ents[2].size() << ", ncells = " << ents[3].size() << std::endl;
339
340 for( int l = 0; l < num_levels; l++ )
341 {
342 all_ents.clear();
343 ents[0].clear();
344 ents[1].clear();
345 ents[2].clear();
346 ents[3].clear();
347 error = mb->get_entities_by_handle( lsets[l + 1], all_ents );MB_CHK_ERR( error );
348
349 for( int k = 0; k < 4; k++ )
350 ents[k] = all_ents.subset_by_dimension( k );
351
352 // std::cout<<std::endl;
353
354 if( qc_vol )
355 {
356 double volume;
357 error = get_max_volume( *moab, lsets[l + 1], dim, volume );MB_CHK_ERR( error );
358 std::cout << "Mesh size for level " << l + 1 << " :: nverts = " << ents[0].size()
359 << ", nedges = " << ents[1].size() << ", nfaces = " << ents[2].size()
360 << ", ncells = " << ents[3].size() << " :: Vmax = " << volume << std::endl;
361 }
362 else
363 std::cout << "Mesh size for level " << l + 1 << " :: nverts = " << ents[0].size()
364 << ", nedges = " << ents[1].size() << ", nfaces = " << ents[2].size()
365 << ", ncells = " << ents[3].size() << std::endl;
366 }
367 }
368
369 if( output )
370 {
371 for( int l = 0; l < num_levels; l++ )
372 {
373 std::string::size_type idx1 = infile.find_last_of( "\\/" );
374 std::string::size_type idx2 = infile.find_last_of( "." );
375 std::string file = infile.substr( idx1 + 1, idx2 - idx1 - 1 );
376 std::stringstream out;
377 if( parallel )
378 out << "_ML_" << l + 1 << ".h5m";
379 else
380 out << "_ML_" << l + 1 << ".vtk";
381 file = file + out.str();
382 const char* output_file = file.c_str();
383 #ifdef MOAB_HAVE_MPI
384 error = mb->write_file( output_file, 0, ";;PARALLEL=WRITE_PART", &lsets[l + 1], 1 );MB_CHK_ERR( error );
385 #else
386 error = mb->write_file( output_file, 0, NULL, &lsets[l + 1], 1 );MB_CHK_ERR( error );
387 #endif
388 // const char* output_file = file.c_str();
389 // error = mb->write_file(output_file, 0, write_options.c_str(), &lsets[l+1],
390 // 1);MB_CHK_ERR(error);
391 // mb->list_entity(lsets[l+1]);
392 // mb->write_file(output_file, 0, "PARALLEL=WRITE_PART", &lsets[l+1], 1);
393 }
394 }
395
396 delete uref;
397 #ifdef MOAB_HAVE_MPI
398 delete pc;
399 #endif
400
401 delete[] ldeg;
402 delete moab;
403
404 #ifdef MOAB_HAVE_MPI
405 MPI_Finalize();
406 #endif
407
408 exit( SUCCESS );
409 }
References moab::Range::clear(), moab::Core::create_meshset(), dim, moab::error(), ErrorCode, moab::NestedRefine::generate_mesh_hierarchy(), get_degree_seq(), moab::Core::get_entities_by_handle(), get_max_volume(), moab::Core::load_file(), make_opts_string(), mb, MB_CHK_ERR, MESHSET_SET, output, parse_id_list(), print_usage(), size, moab::Range::size(), moab::Range::subset_by_dimension(), SUCCESS, moab::NestedRefine::timeall, moab::NestedRefine::codeperf::tm_refine, moab::NestedRefine::codeperf::tm_resolve, moab::NestedRefine::codeperf::tm_total, USAGE_ERROR, usage_error(), and moab::Core::write_file().
| bool make_opts_string | ( | std::vector< std::string > | options, |
| std::string & | opts | ||
| ) |
Definition at line 631 of file umr.cpp.
632 {
633 opts.clear();
634 if( options.empty() ) return true;
635
636 // choose a separator character
637 std::vector< std::string >::const_iterator i;
638 char separator = '\0';
639 const char* alt_separators = ";+,:\t\n";
640 for( const char* sep_ptr = alt_separators; *sep_ptr; ++sep_ptr )
641 {
642 bool seen = false;
643 for( i = options.begin(); i != options.end(); ++i )
644 if( i->find( *sep_ptr, 0 ) != std::string::npos )
645 {
646 seen = true;
647 break;
648 }
649 if( !seen )
650 {
651 separator = *sep_ptr;
652 break;
653 }
654 }
655 if( !separator )
656 {
657 std::cerr << "Error: cannot find separator character for options string" << std::endl;
658 return false;
659 }
660 if( separator != ';' )
661 {
662 opts = ";";
663 opts += separator;
664 }
665
666 // concatenate options
667 i = options.begin();
668 opts += *i;
669 for( ++i; i != options.end(); ++i )
670 {
671 opts += separator;
672 opts += *i;
673 }
674
675 return true;
676 }
Referenced by main().
| bool parse_id_list | ( | const char * | string, |
| int | dim, | ||
| int | nval, | ||
| std::vector< int > & | results | ||
| ) |
Definition at line 578 of file umr.cpp.
579 {
580 bool okay = true;
581 char* mystr = strdup( string );
582 for( const char* ptr = strtok( mystr, "," ); ptr; ptr = strtok( 0, "," ) )
583 {
584 char* endptr;
585 int val = strtol( ptr, &endptr, 0 );
586
587 if( dim == 1 || dim == 2 )
588 {
589 if( val != 2 && val != 3 && val != 5 )
590 {
591 std::cerr << "Not a valid degree for the passed dimension" << std::endl;
592 okay = false;
593 break;
594 }
595 }
596 else if( dim == 3 )
597 {
598 if( val != 2 && val != 3 )
599 {
600 std::cerr << "Not a valid degree for the passed dimension" << std::endl;
601 okay = false;
602 break;
603 }
604 }
605
606 if( endptr == ptr || val <= 0 )
607 {
608 std::cerr << "Not a valid id: " << ptr << std::endl;
609 okay = false;
610 break;
611 }
612
613 results.push_back( val );
614 }
615
616 if( (int)results.size() < nval )
617 {
618 for( int i = results.size(); i <= nval - 1; i++ )
619 results.push_back( results[0] );
620 }
621 else if( (int)results.size() > nval )
622 {
623 for( int i = results.size(); i > nval; i-- )
624 results.pop_back();
625 }
626
627 free( mystr );
628 return okay;
629 }
References dim.
Referenced by main().
|
static |
Definition at line 35 of file umr.cpp.
36 {
37 stream << "Usage: " << name << " <options> <input_file> [<input_file2> ...]" << std::endl
38 << "Options: " << std::endl
39 << "\t-h - Print this help text and exit." << std::endl
40 << "\t-d - Dimension of the mesh." << std::endl
41 << "\t-n - Exact or a maximum number of levels for the hierarchy. Default 1." << std::endl
42 << "\t-L - Degree of refinement for each level. Pass an array or a number. It "
43 "is mandatory to pass dimension and num_levels before to use this option. If this flag "
44 "is not used, a default of deg 2 refinement is used. "
45 << std::endl
46 << "\t-V - Pass a desired volume (absolute) . This will generate a hierarchy "
47 "such that the maximum volume is reduced to the given amount approximately. The length "
48 "of the hierarchy can be constrained further if a maximum number of levels is passed. "
49 "It is mandatory to pass the dimension for this option. "
50 << std::endl
51 << "\t-q - Prints out the maximum volume of the mesh and exits the program. "
52 "This option can be used as a guide to volume constrainted mesh hierarchy later. "
53 << std::endl
54 << "\t-t - Print out the time taken to generate hierarchy." << std::endl
55 << "\t-s - Print out the mesh sizes of each level of the generated hierarchy." << std::endl
56 << "\t-o - Specify true for output files for the mesh levels of the hierarchy." << std::endl
57 //<< "\t-O option - Specify read option." << std::endl
58 #ifdef MOAB_HAVE_MPI
59 << "\t-p[0|1|2] - Read in parallel[0], optionally also doing resolve_shared_ents (1) "
60 "and exchange_ghosts (2)"
61 << std::endl
62 #endif
63 ;
64 exit( USAGE_ERROR );
65 }
References USAGE_ERROR.
Referenced by main(), and usage_error().
|
static |
Definition at line 67 of file umr.cpp.
68 {
69 print_usage( name, std::cerr );
70 #ifdef MOAB_HAVE_MPI
71 MPI_Finalize();
72 #endif
73 exit( USAGE_ERROR );
74 }
References print_usage(), and USAGE_ERROR.
Referenced by main().
1.9.1.