164 MPI_Init( &argc, &argv );
165 MPI_Comm_rank( MPI_COMM_WORLD, &proc_id );
168 #ifdef MOAB_HAVE_TEMPESTREMAP
169 bool tempestin =
false, tempestout =
false;
177 #if( defined( MOAB_HAVE_MPI ) && defined( MOAB_HAVE_TEMPESTREMAP ) )
181 bool append_rank =
false;
182 bool percent_rank_subst =
false;
183 bool file_written =
false;
185 std::list< std::string >::iterator j;
186 bool dims[4] = {
false,
false,
false,
false };
187 const char* format = NULL;
188 std::list< std::string > in;
191 std::set< int > geom[4], mesh[4];
192 std::vector< EntityHandle > set_list;
193 std::vector< std::string > write_opts, read_opts;
194 std::string metis_partition_file;
195 #ifdef MOAB_HAVE_TEMPESTREMAP
196 std::string globalid_tag_name;
197 int spectral_order = 1;
198 bool unitscaling =
false;
203 const char*
const geom_names[] = {
"VERTEX",
"CURVE",
"SURFACE",
"VOLUME" };
207 bool print_times =
false;
208 bool generate[] = {
false,
false,
false };
210 bool parallel =
false, resolve_shared =
false, exchange_ghosts =
false;
211 for( i = 1; i < argc; i++ )
215 if( do_flag && argv[i][0] ==
'-' )
217 if( !argv[i][1] || ( argv[i][1] !=
'M' && argv[i][2] ) )
usage_error( argv[0] );
245 percent_rank_subst =
true;
249 if( argv[i][2] ==
'1' || argv[i][2] ==
'2' ) resolve_shared =
true;
250 if( argv[i][2] ==
'2' ) exchange_ghosts =
true;
253 #ifdef MOAB_HAVE_TEMPESTREMAP
267 dims[argv[i][1] -
'0'] =
true;
272 if( i == argc || argv[i][0] ==
'-' )
274 std::cerr <<
"Expected argument following " << argv[i - 1] << std::endl;
277 if( argv[i - 1][1] ==
'I' )
279 dim = atoi( argv[i] );
280 if( dim < 1 || dim > 2 )
282 std::cerr <<
"Invalid dimension value following -I" << std::endl;
285 generate[
dim] =
true;
289 switch( argv[i - 1][1] )
292 read_opts.push_back( std::string(
"SAT_FILE=" ) + argv[i] );
300 write_opts.push_back( argv[i] );
304 read_opts.push_back( argv[i] );
307 #ifdef MOAB_HAVE_TEMPESTREMAP
309 globalid_tag_name = std::string( argv[i] );
313 spectral_order = atoi( argv[i] );
342 metis_partition_file = argv[i];
346 std::cerr <<
"Invalid option: " << argv[i] << std::endl;
351 std::cerr <<
"Invalid flag or flag value: " << argv[i - 1] <<
" " << argv[i] << std::endl;
359 in.push_back( argv[i] );
364 std::cerr <<
"No output file name specified." << std::endl;
373 std::ostringstream mod;
374 mod << out <<
"." << proc_id;
378 if( percent_rank_subst )
380 for( j = in.begin(); j != in.end(); ++j )
386 std::string read_options, write_options;
389 read_opts.push_back(
"PARALLEL=READ_PART" );
390 read_opts.push_back(
"PARTITION=PARALLEL_PARTITION" );
391 if( !append_rank && !percent_rank_subst ) write_opts.push_back(
"PARALLEL=WRITE_PART" );
393 if( resolve_shared ) read_opts.push_back(
"PARALLEL_RESOLVE_SHARED_ENTS" );
394 if( exchange_ghosts ) read_opts.push_back(
"PARALLEL_GHOSTS=3.0.1" );
404 if( !metis_partition_file.empty() )
406 if( ( in.size() != 1 ) || ( proc_id != 0 ) )
408 std::cerr <<
" mpas partition allows only one input file, in serial conversion\n";
419 #ifdef MOAB_HAVE_TEMPESTREMAP
420 if( tempestin && in.size() > 1 )
422 std::cerr <<
" we can read only one tempest files at a time\n";
429 if( tempestin or tempestout )
438 bool use_overlap_context =
false;
439 Tag srcParentTag, tgtParentTag;
442 for( j = in.begin(); j != in.end(); ++j )
444 std::string inFileName = *j;
448 #ifdef MOAB_HAVE_TEMPESTREMAP
463 tempestMesh->RemoveZeroEdges();
464 tempestMesh->RemoveCoincidentNodes();
475 NcFile ncInput( inFileName.c_str(), NcFile::ReadOnly );
477 NcError error_temp( NcError::silent_nonfatal );
479 NcAtt* attRectilinear = ncInput.get_att(
"rectilinear" );
480 NcVar* varGridDims = ncInput.get_var(
"grid_dims" );
483 std::vector< int > vecDimSizes( 3, 0 );
491 if( attRectilinear !=
nullptr )
494 NcAtt* attRectilinearDim0Size = ncInput.get_att(
"rectilinear_dim0_size" );
495 NcAtt* attRectilinearDim1Size = ncInput.get_att(
"rectilinear_dim1_size" );
497 if( attRectilinearDim0Size ==
nullptr )
499 _EXCEPTIONT(
"Missing attribute \"rectilinear_dim0_size\"" );
501 if( attRectilinearDim1Size ==
nullptr )
503 _EXCEPTIONT(
"Missing attribute \"rectilinear_dim1_size\"" );
506 int nDim0Size = attRectilinearDim0Size->as_int( 0 );
507 int nDim1Size = attRectilinearDim1Size->as_int( 0 );
510 NcAtt* attRectilinearDim0Name = ncInput.get_att(
"rectilinear_dim0_name" );
511 NcAtt* attRectilinearDim1Name = ncInput.get_att(
"rectilinear_dim1_name" );
513 if( attRectilinearDim0Name ==
nullptr )
515 _EXCEPTIONT(
"Missing attribute \"rectilinear_dim0_name\"" );
517 if( attRectilinearDim1Name ==
nullptr )
519 _EXCEPTIONT(
"Missing attribute \"rectilinear_dim1_name\"" );
522 std::string strDim0Name = attRectilinearDim0Name->as_string( 0 );
523 std::string strDim1Name = attRectilinearDim1Name->as_string( 0 );
525 std::map< std::string, int > vecDimNameSizes;
527 vecDimNameSizes[strDim0Name] = nDim0Size;
528 vecDimNameSizes[strDim1Name] = nDim1Size;
530 vecDimSizes[1] = vecDimNameSizes[
"lat"];
531 vecDimSizes[2] = vecDimNameSizes[
"lon"];
533 else if( varGridDims !=
nullptr )
536 NcDim* dimGridRank = varGridDims->get_dim( 0 );
537 if( dimGridRank == NULL )
539 _EXCEPTIONT(
"Variable \"grid_dims\" has no dimensions" );
542 int gridrank = dimGridRank->size();
545 varGridDims->get( &( gridsizes[0] ), dimGridRank->size() );
556 vecDimSizes[1] = elems.
size();
562 vecDimSizes[1] = gridsizes[0];
563 vecDimSizes[2] = gridsizes[1];
575 vecDimSizes[1] = elems.
size();
582 switch( vecDimSizes[0] )
585 printf(
"Cubed-Sphere mesh: %d (elements), %d (vertices)\n", vecDimSizes[1], vecDimSizes[2] );
588 printf(
"Rectilinear RLL mesh: (lon) %d X (lat) %d\n", vecDimSizes[2], vecDimSizes[1] );
591 printf(
"Icosahedral (triangular) mesh: %d (elements), %d (vertices)\n", vecDimSizes[1],
596 printf(
"Polygonal mesh: %d (elements), %d (vertices)\n", vecDimSizes[1], vecDimSizes[2] );
602 const size_t nOverlapFaces = tempestMesh->faces.size();
603 if( tempestMesh->vecSourceFaceIx.size() == nOverlapFaces &&
604 tempestMesh->vecSourceFaceIx.size() == nOverlapFaces )
607 use_overlap_context =
true;
618 std::vector< int > gids( faces.
size() ), srcpar( faces.
size() ), tgtpar( faces.
size() );
621 for(
unsigned ii = 0; ii < faces.
size(); ++ii )
623 srcpar[ii] = tempestMesh->vecSourceFaceIx[gids[ii] - 1];
624 tgtpar[ii] = tempestMesh->vecTargetFaceIx[gids[ii] - 1];
635 else if( tempestout )
641 std::vector< int > metadata( 2 );
655 use_overlap_context =
true;
664 result = reorder_tool.handle_order_from_int_tag( srcParentTag, -1, order );
MB_CHK_ERR( result );
665 result = reorder_tool.reorder_entities( order );
MB_CHK_ERR( result );
673 assert( metadata.size() );
674 std::cout <<
"Converting a RLL mesh with rectilinear dimension: " << metadata[0] <<
" X "
675 << metadata[1] << std::endl;
684 result =
gMB->
load_file( j->c_str(), 0, read_options.c_str() );
686 result =
gMB->
load_file( j->c_str(), 0, read_options.c_str() );
690 std::cerr <<
"Failed to load \"" << *j <<
"\"." << std::endl;
691 std::cerr <<
"Error code: " <<
gMB->
get_error_string( result ) <<
" (" << result <<
")" << std::endl;
694 std::cerr <<
"Error message: " << message << std::endl;
700 if( !proc_id ) std::cerr <<
"Read \"" << *j <<
"\"" << std::endl;
701 if( print_times && !proc_id )
write_times( std::cout );
705 bool have_geom =
false;
708 if( !geom[
dim].empty() ) have_geom =
true;
713 bool have_sets = have_geom;
721 std::cerr <<
"No ID tag defined." << std::endl;
727 std::cerr <<
"No geometry tag defined." << std::endl;
736 Tag tags[] = { id_tag, dim_tag };
737 const void* vals[] = { &id_val, &
dim };
740 int init_count = set_list.size();
741 for( std::set< int >::iterator iter = geom[
dim].begin(); iter != geom[
dim].end(); ++iter )
749 std::cerr << geom_names[
dim] <<
" " << id_val <<
" not found.\n";
751 std::copy( range.
begin(), range.
end(), std::back_inserter( set_list ) );
755 std::cout <<
"Found " << ( set_list.size() - init_count ) <<
' ' << geom_names[
dim] <<
" sets"
761 for( i = 0; i < 4; ++i )
765 if( mesh[i].empty() )
continue;
773 std::cerr <<
"Tag not found: " << mesh_tag_names[i] << std::endl;
778 int init_count = set_list.size();
779 for( std::set< int >::iterator iter = mesh[i].begin(); iter != mesh[i].end(); ++iter )
782 const void* vals[] = { &*iter };
787 std::cerr << mesh_tag_names[i] <<
" " << *iter <<
" not found.\n";
789 std::copy( range.
begin(), range.
end(), std::back_inserter( set_list ) );
793 std::cout <<
"Found " << ( set_list.size() - init_count ) <<
' ' << mesh_tag_names[i] <<
" sets"
800 if( dims[
dim] ) bydim =
true;
805 if( generate[1] && !dims[1] )
807 std::cerr <<
"Warning: Request to generate 1D internal entities but not export them." << std::endl;
810 if( generate[2] && !dims[2] )
812 std::cerr <<
"Warning: Request to generate 2D internal entities but not export them." << std::endl;
818 if( generate[1] || generate[2] )
825 num_sets = set_list.size();
828 for( i = 0; i < num_sets; ++i )
830 Range dim3, dim2, adj;
850 Range dead_entities, tmp_range;
853 if( dims[
dim] )
continue;
855 dead_entities.
merge( tmp_range );
861 while( !empty_sets.
empty() )
864 dead_entities.
merge( empty_sets );
867 empty_sets =
subtract( tmp_range, dead_entities );
874 if( have_sets && set_list.empty() )
876 std::cerr <<
"Nothing to write." << std::endl;
884 if( !metis_partition_file.empty() )
889 std::cerr <<
"Failed to process partition file \"" << metis_partition_file <<
"\"." << std::endl;
899 std::cout <<
"Found " << set_list.size() <<
" specified sets to write (total)." << std::endl;
901 std::cout <<
"No sets specifed. Writing entire mesh." << std::endl;
906 #ifdef MOAB_HAVE_TEMPESTREMAP
915 int ntot_elements = 0, nelements = faces.
size();
917 int ierr = MPI_Allreduce( &nelements, &ntot_elements, 1, MPI_INT, MPI_SUM, pcomm->
comm() );
918 if( ierr != 0 )
MB_CHK_SET_ERR( MB_FAILURE,
"MPI_Allreduce failed to get total source elements" );
920 ntot_elements = nelements;
924 std::vector< int > gids( faces.
size() );
927 if( faces.
size() > 1 && gids[0] == gids[1] && !use_overlap_context )
940 if( spectral_order > 1 && globalid_tag_name.size() > 1 )
942 result = remapper->
GenerateMeshMetadata( *tempestMesh, ntot_elements, faces, NULL, globalid_tag_name,
950 if( use_overlap_context &&
false )
952 const int nOverlapFaces = faces.
size();
954 tempestMesh->vecSourceFaceIx.resize( nOverlapFaces );
955 tempestMesh->vecTargetFaceIx.resize( nOverlapFaces );
960 tempestMesh->Write( out, NcFile::Netcdf4 );
970 result =
gMB->
write_file( out.c_str(), format, write_options.c_str(), &set_list[0], set_list.size() );
972 result =
gMB->
write_file( out.c_str(), format, write_options.c_str() );
975 std::cerr <<
"Failed to write \"" << out <<
"\"." << std::endl;
976 std::cerr <<
"Error code: " <<
gMB->
get_error_string( result ) <<
" (" << result <<
")" << std::endl;
979 std::cerr <<
"Error message: " << message << std::endl;
987 if( !proc_id ) std::cerr <<
"Wrote \"" << out <<
"\"" << std::endl;
988 if( print_times && !proc_id )
write_times( std::cout );