Child helper class for Domain grid. More...
#include <NCHelperDomain.hpp>
Public Member Functions | |
NCHelperDomain (ReadNC *readNC, int fileId, const FileOptions &opts, EntityHandle fileSet) | |
ErrorCode | create_mesh (Range &faces) |
Implementation of NCHelper::create_mesh() More... | |
![]() | |
ScdNCHelper (ReadNC *readNC, int fileId, const FileOptions &opts, EntityHandle fileSet) | |
virtual | ~ScdNCHelper () |
![]() | |
NCHelper (ReadNC *readNC, int fileId, const FileOptions &opts, EntityHandle fileSet) | |
virtual | ~NCHelper () |
ErrorCode | create_conventional_tags (const std::vector< int > &tstep_nums) |
Create NC conventional tags. More... | |
ErrorCode | update_time_tag_vals () |
Update time tag values if timesteps spread across files. More... | |
Static Public Member Functions | |
static bool | can_read_file (ReadNC *readNC, int fileId) |
![]() | |
static ReadNC::NCFormatType | get_nc_format (ReadNC *readNC, int fileId) |
Get appropriate format to read the file. More... | |
static std::string | get_default_ncformat_options (ReadNC::NCFormatType format) |
Get appropriate format to read the file. More... | |
static NCHelper * | get_nc_helper (ReadNC *readNC, int fileId, const FileOptions &opts, EntityHandle fileSet) |
Get appropriate helper instance for ReadNC class. More... | |
Private Member Functions | |
virtual ErrorCode | init_mesh_vals () |
Interfaces to be implemented in child classes. More... | |
virtual std::string | get_mesh_type_name () |
Private Attributes | |
int | nv |
int | nvDim |
Additional Inherited Members | |
![]() | |
ErrorCode | read_variables_setup (std::vector< std::string > &var_names, std::vector< int > &tstep_nums, std::vector< ReadNC::VarData > &vdatas, std::vector< ReadNC::VarData > &vsetdatas) |
Separate set and non-set variables (common to scd mesh and ucd mesh) More... | |
ErrorCode | read_variables_to_set (std::vector< ReadNC::VarData > &vdatas, std::vector< int > &tstep_nums) |
Read set variables (common to scd mesh and ucd mesh) More... | |
ErrorCode | read_coordinate (const char *var_name, int lmin, int lmax, std::vector< double > &cvals) |
ErrorCode | get_tag_to_set (ReadNC::VarData &var_data, int tstep_num, Tag &tagh) |
ErrorCode | get_tag_to_nonset (ReadNC::VarData &var_data, int tstep_num, Tag &tagh, int num_lev) |
ErrorCode | create_attrib_string (const std::map< std::string, ReadNC::AttData > &attMap, std::string &attString, std::vector< int > &attLen) |
Create a character string attString of attMap. with '\0' terminating each attribute name, ';' separating the data type and value, and ';' separating one name/data type/value from the next'. attLen stores the end position for each name/data type/ value. More... | |
ErrorCode | create_dummy_variables () |
For a dimension that does not have a corresponding coordinate variable (e.g. ncol for HOMME), create a dummy variable with a sparse tag to store the dimension length. More... | |
![]() | |
int | gDims [6] |
Dimensions of global grid in file. More... | |
int | lDims [6] |
Dimensions of my local part of grid. More... | |
int | gCDims [6] |
Center dimensions of global grid in file. More... | |
int | lCDims [6] |
Center dimensions of my local part of grid. More... | |
std::vector< double > | ilVals |
Values for i/j. More... | |
std::vector< double > | jlVals |
std::vector< double > | ilCVals |
Center values for i/j. More... | |
std::vector< double > | jlCVals |
int | iDim |
Dimension numbers for i/j. More... | |
int | jDim |
int | iCDim |
Center dimension numbers for i/j. More... | |
int | jCDim |
int | locallyPeriodic [3] |
Whether mesh is locally periodic in i or j or k. More... | |
int | globallyPeriodic [3] |
Whether mesh is globally periodic in i or j or k. More... | |
![]() | |
ReadNC * | _readNC |
Allow NCHelper to directly access members of ReadNC. More... | |
int | _fileId |
Cache some information from ReadNC. More... | |
const FileOptions & | _opts |
EntityHandle | _fileSet |
int | nTimeSteps |
Dimensions of time and level. More... | |
int | nLevels |
std::vector< double > | tVals |
Values for time and level. More... | |
std::vector< double > | levVals |
int | tDim |
Dimension numbers for time and level. More... | |
int | levDim |
std::set< std::string > | ignoredVarNames |
Ignored variables. More... | |
std::set< std::string > | dummyVarNames |
Dummy variables. More... | |
Child helper class for Domain grid.
Definition at line 18 of file NCHelperDomain.hpp.
|
inline |
Definition at line 21 of file NCHelperDomain.hpp.
22 : ScdNCHelper( readNC, fileId, opts, fileSet ) 23 { 24 }
|
static |
Definition at line 20 of file NCHelperDomain.cpp.
21 {
22 std::vector< std::string >& dimNames = readNC->dimNames;
23
24 // If dimension names "n" AND "ni" AND "nj" AND "nv" exist then it should be the Domain grid
25 // some files do not have n, which should be just ni*nj
26 if( ( std::find( dimNames.begin(), dimNames.end(), std::string( "ni" ) ) != dimNames.end() ) &&
27 ( std::find( dimNames.begin(), dimNames.end(), std::string( "nj" ) ) != dimNames.end() ) &&
28 ( std::find( dimNames.begin(), dimNames.end(), std::string( "nv" ) ) != dimNames.end() ) )
29 {
30 // we do not test anymore if it is an actual CAM file
31 return true;
32 }
33
34 return false;
35 }
References moab::ReadNC::dimNames.
Referenced by moab::NCHelper::get_nc_format().
Implementation of NCHelper::create_mesh()
Reimplemented from moab::ScdNCHelper.
Definition at line 230 of file NCHelperDomain.cpp.
231 {
232 Interface*& mbImpl = _readNC->mbImpl;
233 // std::string& fileName = _readNC->fileName;
234 Tag& mGlobalIdTag = _readNC->mGlobalIdTag;
235 // const Tag*& mpFileIdTag = _readNC->mpFileIdTag;
236 DebugOutput& dbgOut = _readNC->dbgOut;
237
238 bool& culling = _readNC->culling;
239 bool& repartition = _readNC->repartition;
240
241 ErrorCode rval;
242 int success = 0;
243
244 int local_elems = ( lDims[4] - lDims[1] ) * ( lDims[3] - lDims[0] );
245 dbgOut.tprintf( 1, "local cells: %d \n", local_elems );
246
247 // count how many will be with mask 1 here
248 // basically, read the mask variable on the local elements;
249 std::string maskstr( "mask" );
250 ReadNC::VarData& vmask = _readNC->varInfo[maskstr];
251
252 // mask is (nj, ni)
253 vmask.readStarts.push_back( lDims[1] );
254 vmask.readStarts.push_back( lDims[0] );
255 vmask.readCounts.push_back( lDims[4] - lDims[1] );
256 vmask.readCounts.push_back( lDims[3] - lDims[0] );
257 std::vector< int > mask( local_elems );
258 success = NCFUNCAG( _vara_int )( _fileId, vmask.varId, &vmask.readStarts[0], &vmask.readCounts[0], &mask[0] );
259 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read int data for mask variable " );
260
261 std::vector< int > gids( local_elems );
262 int elem_index = 0;
263 int global_row_size = gDims[3] - gDims[0]; // this is along first dimension in global decomposition
264 // create global id array for cells, for all cells, including those with 0 mask; which will be not used eventually
265 for( int j = lCDims[1]; j < lCDims[4]; j++ )
266 for( int i = lCDims[0]; i < lCDims[3]; i++ )
267 {
268 gids[elem_index] = j * global_row_size + i + 1;
269 elem_index++;
270 }
271 std::vector< NCDF_SIZE > startsv( 3 );
272 startsv[0] = vmask.readStarts[0];
273 startsv[1] = vmask.readStarts[1];
274 startsv[2] = 0;
275 std::vector< NCDF_SIZE > countsv( 3 );
276 countsv[0] = vmask.readCounts[0];
277 countsv[1] = vmask.readCounts[1];
278 countsv[2] = nv; // number of vertices per element
279
280 // read xv and yv coords for vertices, and create elements;
281 std::string xvstr( "xv" );
282 ReadNC::VarData& var_xv = _readNC->varInfo[xvstr];
283 std::vector< double > xv( local_elems * nv );
284
285 std::vector< NCDF_SIZE > starts_vv, counts_vv;
286 starts_vv.resize( 3 );
287 counts_vv.resize( 3 );
288 bool nv_last = true;
289 if( _readNC->dimNames[var_xv.varDims[0]] == std::string( "nv" ) ) // it means that nv is the first dimension
290 {
291 starts_vv[0] = 0;
292 starts_vv[1] = startsv[0];
293 starts_vv[2] = startsv[1];
294 counts_vv[0] = nv;
295 counts_vv[1] = countsv[0];
296 counts_vv[2] = countsv[1];
297 nv_last = false;
298 }
299 else
300 {
301 starts_vv = startsv;
302 counts_vv = countsv;
303 }
304 dbgOut.tprintf( 1, " nv is the last dimension in xv array (0 or 1) %d \n", (int)nv_last );
305 success = NCFUNCAG( _vara_double )( _fileId, var_xv.varId, &starts_vv[0], &counts_vv[0], &xv[0] );
306 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read double data for xv variable " );
307
308 std::string yvstr( "yv" );
309 ReadNC::VarData& var_yv = _readNC->varInfo[yvstr];
310 std::vector< double > yv( local_elems * nv );
311 success = NCFUNCAG( _vara_double )( _fileId, var_yv.varId, &starts_vv[0], &counts_vv[0], &yv[0] );
312 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read double data for yv variable " );
313
314 // read other variables, like xc, yc, frac, area
315 std::string xcstr( "xc" );
316 ReadNC::VarData& var_xc = _readNC->varInfo[xcstr];
317 std::vector< double > xc( local_elems );
318 success = NCFUNCAG( _vara_double )( _fileId, var_xc.varId, &vmask.readStarts[0], &vmask.readCounts[0], &xc[0] );
319 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read double data for xc variable " );
320
321 std::string ycstr( "yc" );
322 ReadNC::VarData& var_yc = _readNC->varInfo[ycstr];
323 std::vector< double > yc( local_elems );
324 success = NCFUNCAG( _vara_double )( _fileId, var_yc.varId, &vmask.readStarts[0], &vmask.readCounts[0], &yc[0] );
325 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read double data for yc variable " );
326
327 std::string fracstr( "frac" );
328 ReadNC::VarData& var_frac = _readNC->varInfo[fracstr];
329 std::vector< double > frac( local_elems );
330 if( var_frac.varId >= 0 )
331 {
332 success =
333 NCFUNCAG( _vara_double )( _fileId, var_frac.varId, &vmask.readStarts[0], &vmask.readCounts[0], &frac[0] );
334 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read double data for frac variable " );
335 }
336 std::string areastr( "area" );
337 std::vector< double > area( local_elems );
338 ReadNC::VarData& var_area = _readNC->varInfo[areastr];
339 success = NCFUNCAG( _vara_double )( _fileId, var_area.varId, &vmask.readStarts[0], &vmask.readCounts[0], &area[0] );
340 if( success ) MB_SET_ERR( MB_FAILURE, "Failed to read double data for area variable " );
341 // create tags for them
342 Tag areaTag, fracTag, xcTag, ycTag;
343 rval = mbImpl->tag_get_handle( "area", 1, MB_TYPE_DOUBLE, areaTag, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Trouble creating area tag" );
344 rval = mbImpl->tag_get_handle( "frac", 1, MB_TYPE_DOUBLE, fracTag, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Trouble creating frac tag" );
345 rval = mbImpl->tag_get_handle( "xc", 1, MB_TYPE_DOUBLE, xcTag, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Trouble creating xc tag" );
346 rval = mbImpl->tag_get_handle( "yc", 1, MB_TYPE_DOUBLE, ycTag, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Trouble creating yc tag" );
347
348 // create tags for GRID_IMASK, which will be the same name as the Scrip helper tag that holds the mask
349 // create the maskTag GRID_IMASK, with default value of 1
350 Tag maskTag;
351 int def_val = 1;
352 rval = mbImpl->tag_get_handle( "GRID_IMASK", 1, MB_TYPE_INTEGER, maskTag, MB_TAG_DENSE | MB_TAG_CREAT, &def_val );MB_CHK_SET_ERR( rval, "Trouble creating GRID_IMASK tag" );
353
354 // will now look to repartition the cells, using zoltan, looking at the xc and yc coordinates in 2d, convert to 3d,
355 // and decide based on those partitioning info to what task to send each cell, along with its vertices, and global id
356 #ifdef MOAB_HAVE_MPI
357 int procs = 1;
358 bool& isParallel = _readNC->isParallel;
359 ParallelComm* myPcomm = NULL;
360 if( isParallel )
361 {
362 myPcomm = _readNC->myPcomm;
363 procs = myPcomm->proc_config().proc_size();
364 }
365
366 if( procs >= 2 && repartition )
367 {
368 // Redistribute local cells after trivial partition (e.g. apply Zoltan partition)
369 rval = redistribute_cells( myPcomm, xc, yc, xv, yv, frac, mask, area, gids, nv, nv_last );MB_CHK_SET_ERR( rval, "Failed to redistribute local cells" );
370 local_elems = (int)xc.size();
371 dbgOut.tprintf( 1, "local cells after repartition: %d \n", local_elems );
372 }
373
374 #endif
375
376 int nb_with_mask1 = 0;
377 for( int i = 0; i < local_elems; i++ )
378 if( 1 == mask[i] ) nb_with_mask1++;
379 dbgOut.tprintf( 1, "local cells with mask 1: %d \n", nb_with_mask1 );
380
381 EntityHandle* conn_arr;
382 EntityHandle vtx_handle;
383 Range tmp_range;
384
385 // set connectivity into that space
386
387 EntityHandle start_cell;
388 EntityType mdb_type = MBVERTEX;
389 if( nv == 3 )
390 mdb_type = MBTRI;
391 else if( nv == 4 )
392 mdb_type = MBQUAD;
393 else if( nv > 4 ) // (nv > 4)
394 mdb_type = MBPOLYGON;
395 // for nv = 1 , type is vertex
396
397 int num_actual_cells = local_elems;
398 if( culling ) num_actual_cells = nb_with_mask1;
399
400 if( nv > 1 && num_actual_cells > 0 )
401 {
402 rval = _readNC->readMeshIface->get_element_connect( num_actual_cells, nv, mdb_type, 0, start_cell, conn_arr );MB_CHK_SET_ERR( rval, "Failed to create local cells" );
403 tmp_range.insert( start_cell, start_cell + num_actual_cells - 1 );
404 }
405
406 // Create vertices; first identify different ones, with a tolerance
407 std::map< Node3D, EntityHandle > vertex_map;
408
409 if( num_actual_cells > 0 )
410 {
411 // Set vertex coordinates
412 // will read all xv, yv, but use only those with correct mask on
413
414 // total index in netcdf arrays
415 const double pideg = acos( -1.0 ) / 180.0;
416
417 for( elem_index = 0; elem_index < local_elems; elem_index++ )
418 {
419 if( culling && 0 == mask[elem_index] )
420 continue; // nothing to do, do not advance elem_index in actual moab arrays
421 // set area and fraction on those elements too
422 dbgOut.tprintf( 3, "elem index %d \n", elem_index );
423 for( int k = 0; k < nv; k++ )
424 {
425 int index_v_arr = nv * elem_index + k;
426 if( !nv_last ) index_v_arr = k * local_elems + elem_index;
427 double x, y;
428 if( nv > 1 )
429 {
430 x = xv[index_v_arr];
431 y = yv[index_v_arr];
432 double cosphi = cos( pideg * y );
433 double zmult = sin( pideg * y );
434 double xmult = cosphi * cos( x * pideg );
435 double ymult = cosphi * sin( x * pideg );
436 Node3D pt( xmult, ymult, zmult );
437 vertex_map[pt] = 0;
438 dbgOut.tprintf( 3, " %d x=%5.1f y=%5.1f pt: %f %f %f \n", k, x, y, pt.coords[0], pt.coords[1],
439 pt.coords[2] );
440 }
441 else
442 {
443 x = xc[elem_index];
444 y = yc[elem_index];
445 Node3D pt( x, y, 0 );
446 vertex_map[pt] = 0;
447 }
448 }
449 }
450 int nLocalVertices = (int)vertex_map.size();
451 std::vector< double* > arrays;
452 EntityHandle start_vertex;
453 rval = _readNC->readMeshIface->get_node_coords( 3, nLocalVertices, 0, start_vertex, arrays );MB_CHK_SET_ERR( rval, "Failed to create local vertices" );
454
455 vtx_handle = start_vertex;
456 // Copy vertex coordinates into entity sequence coordinate arrays
457 // and copy handle into vertex_map.
458 double *x = arrays[0], *y = arrays[1], *z = arrays[2];
459 for( auto i = vertex_map.begin(); i != vertex_map.end(); ++i )
460 {
461 i->second = vtx_handle;
462 ++vtx_handle;
463 *x = i->first.coords[0];
464 ++x;
465 *y = i->first.coords[1];
466 ++y;
467 *z = i->first.coords[2];
468 ++z;
469 }
470
471 // int nj = gDims[4]-gDims[1]; // is it about 1 in irregular cases
472
473 // int local_row_size = lCDims[3] - lCDims[0];
474 int index = 0; // consider the mask for advancing in moab arrays;
475
476 // create now vertex arrays, size vertex_map.size()
477 for( elem_index = 0; elem_index < local_elems; elem_index++ )
478 {
479 if( culling && 0 == mask[elem_index] )
480 continue; // nothing to do, do not advance elem_index in actual moab arrays
481 // set area and fraction on those elements too
482 for( int k = 0; k < nv; k++ )
483 {
484 int index_v_arr = nv * elem_index + k;
485 if( !nv_last ) index_v_arr = k * local_elems + elem_index;
486 if( nv > 1 )
487 {
488 double x = xv[index_v_arr];
489 double y = yv[index_v_arr];
490 double cosphi = cos( pideg * y );
491 double zmult = sin( pideg * y );
492 double xmult = cosphi * cos( x * pideg );
493 double ymult = cosphi * sin( x * pideg );
494 Node3D pt( xmult, ymult, zmult );
495 conn_arr[index * nv + k] = vertex_map[pt];
496 }
497 }
498 EntityHandle cell = start_vertex + index;
499 if( nv > 1 ) cell = start_cell + index;
500 // set other tags, like xc, yc, frac, area
501 rval = mbImpl->tag_set_data( xcTag, &cell, 1, &xc[elem_index] );MB_CHK_SET_ERR( rval, "Failed to set xc tag" );
502 rval = mbImpl->tag_set_data( ycTag, &cell, 1, &yc[elem_index] );MB_CHK_SET_ERR( rval, "Failed to set yc tag" );
503 rval = mbImpl->tag_set_data( areaTag, &cell, 1, &area[elem_index] );MB_CHK_SET_ERR( rval, "Failed to set area tag" );
504 rval = mbImpl->tag_set_data( fracTag, &cell, 1, &frac[elem_index] );MB_CHK_SET_ERR( rval, "Failed to set frac tag" );
505 rval = mbImpl->tag_set_data( maskTag, &cell, 1, &mask[elem_index] );MB_CHK_SET_ERR( rval, "Failed to set mask tag" );
506
507 // set the global id too:
508 int globalId = gids[elem_index];
509 rval = mbImpl->tag_set_data( mGlobalIdTag, &cell, 1, &globalId );MB_CHK_SET_ERR( rval, "Failed to set global id tag" );
510 index++;
511 }
512
513 rval = mbImpl->add_entities( _fileSet, tmp_range );MB_CHK_SET_ERR( rval, "Failed to add new cells to current file set" );
514
515 // modify local file set, to merge coincident vertices, and to correct repeated vertices in elements
516 std::vector< Tag > tagList;
517 tagList.push_back( mGlobalIdTag );
518 tagList.push_back( xcTag );
519 tagList.push_back( ycTag );
520 tagList.push_back( areaTag );
521 tagList.push_back( fracTag );
522 tagList.push_back( maskTag ); // not sure this is needed though ? on cells or on vertices?
523 rval = IntxUtils::remove_padded_vertices( mbImpl, _fileSet, tagList );MB_CHK_SET_ERR( rval, "Failed to remove duplicate vertices" );
524
525 rval = mbImpl->get_entities_by_dimension( _fileSet, 2, faces );MB_CHK_ERR( rval );
526 Range all_verts;
527 rval = mbImpl->get_connectivity( faces, all_verts );MB_CHK_ERR( rval );
528 //printf(" range vert size :%ld \n", all_verts.size());
529 rval = mbImpl->add_entities( _fileSet, all_verts );MB_CHK_ERR( rval );
530
531 // need to add adjacencies; TODO: fix this for all nc readers
532 // copy this logic from migrate mesh in par comm graph
533 Core* mb = (Core*)mbImpl;
534 AEntityFactory* adj_fact = mb->a_entity_factory();
535 if( !adj_fact->vert_elem_adjacencies() )
536 adj_fact->create_vert_elem_adjacencies();
537 else
538 {
539 for( Range::iterator it = faces.begin(); it != faces.end(); ++it )
540 {
541 EntityHandle eh = *it;
542 const EntityHandle* conn = NULL;
543 int num_nodes = 0;
544 rval = mb->get_connectivity( eh, conn, num_nodes );MB_CHK_ERR( rval );
545 adj_fact->notify_create_entity( eh, conn, num_nodes );
546 }
547 }
548 }
549
550 #ifdef MOAB_HAVE_MPI
551 myPcomm = _readNC->myPcomm; // we will have to set the global id on vertices in any case, even
552 if( myPcomm && procs >= 2 )
553 {
554 double tol = 1.e-12; // this is the same as static tolerance in NCHelper
555 ParallelMergeMesh pmm( myPcomm, tol );
556 rval = pmm.merge( _fileSet,
557 /* do not do local merge*/ false,
558 /* 2d cells*/ 2 );MB_CHK_SET_ERR( rval, "Failed to merge vertices in parallel" );
559 // assign global ids only for vertices, cells have them fine
560 rval = myPcomm->assign_global_ids( _fileSet, /*dim*/ 0 );MB_CHK_ERR( rval );
561 }
562 #endif
563
564 return MB_SUCCESS;
565 }
References moab::NCHelper::_fileId, moab::NCHelper::_fileSet, moab::NCHelper::_readNC, moab::Interface::add_entities(), moab::ParallelComm::assign_global_ids(), moab::Range::begin(), moab::NCHelper::Node3D::coords, moab::AEntityFactory::create_vert_elem_adjacencies(), moab::ReadNC::culling, moab::ReadNC::dbgOut, moab::ReadNC::dimNames, moab::Range::end(), ErrorCode, moab::ScdNCHelper::gDims, moab::Interface::get_connectivity(), moab::ReadUtilIface::get_element_connect(), moab::Interface::get_entities_by_dimension(), moab::ReadUtilIface::get_node_coords(), moab::Range::insert(), moab::ReadNC::isParallel, moab::ScdNCHelper::lCDims, moab::ScdNCHelper::lDims, mb, MB_CHK_ERR, MB_CHK_SET_ERR, MB_SET_ERR, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_DENSE, MB_TYPE_DOUBLE, MB_TYPE_INTEGER, moab::ReadNC::mbImpl, MBPOLYGON, MBQUAD, MBTRI, MBVERTEX, moab::ParallelMergeMesh::merge(), moab::ReadNC::mGlobalIdTag, NCFUNCAG, moab::AEntityFactory::notify_create_entity(), nv, moab::pideg, moab::ParallelComm::proc_config(), moab::ProcConfig::proc_size(), moab::ReadNC::VarData::readCounts, moab::ReadNC::readMeshIface, moab::ReadNC::VarData::readStarts, moab::IntxUtils::remove_padded_vertices(), moab::ReadNC::repartition, moab::Interface::tag_get_handle(), moab::Interface::tag_set_data(), moab::DebugOutput::tprintf(), moab::ReadNC::VarData::varDims, moab::ReadNC::VarData::varId, moab::ReadNC::varInfo, and moab::AEntityFactory::vert_elem_adjacencies().
|
inlineprivatevirtual |
Implements moab::NCHelper.
Definition at line 31 of file NCHelperDomain.hpp.
32 {
33 return "DOMAIN";
34 }
|
privatevirtual |
Interfaces to be implemented in child classes.
Implements moab::NCHelper.
Definition at line 37 of file NCHelperDomain.cpp.
38 {
39 Interface*& mbImpl = _readNC->mbImpl;
40 std::vector< std::string >& dimNames = _readNC->dimNames;
41 std::vector< int >& dimLens = _readNC->dimLens;
42 std::map< std::string, ReadNC::VarData >& varInfo = _readNC->varInfo;
43 DebugOutput& dbgOut = _readNC->dbgOut;
44 #ifdef MOAB_HAVE_MPI
45 bool& isParallel = _readNC->isParallel;
46 #endif
47 int& partMethod = _readNC->partMethod;
48 ScdParData& parData = _readNC->parData;
49
50 ErrorCode rval;
51
52 // Look for names of i/j dimensions
53 // First i
54 std::vector< std::string >::iterator vit;
55 unsigned int idx;
56 if( ( vit = std::find( dimNames.begin(), dimNames.end(), "ni" ) ) != dimNames.end() )
57 idx = vit - dimNames.begin();
58 else
59 {
60 MB_SET_ERR( MB_FAILURE, "Couldn't find 'ni' variable" );
61 }
62 iDim = idx;
63 gDims[0] = 0;
64 gDims[3] = dimLens[idx];
65
66 // Then j
67 if( ( vit = std::find( dimNames.begin(), dimNames.end(), "nj" ) ) != dimNames.end() )
68 idx = vit - dimNames.begin();
69 else
70 {
71 MB_SET_ERR( MB_FAILURE, "Couldn't find 'nj' variable" );
72 }
73 jDim = idx;
74 gDims[1] = 0;
75 gDims[4] = dimLens[idx]; // Add 2 for the pole points ? not needed
76
77 // do not use gcdims ? or use only gcdims?
78
79 // Try a truly 2D mesh
80 gDims[2] = -1;
81 gDims[5] = -1;
82
83 // Get number of vertices per cell
84 if( ( vit = std::find( dimNames.begin(), dimNames.end(), "nv" ) ) != dimNames.end() )
85 idx = vit - dimNames.begin();
86 else
87 {
88 MB_SET_ERR( MB_FAILURE, "Couldn't find 'nv' dimension" );
89 }
90 nvDim = idx;
91 nv = dimLens[idx];
92
93 // Parse options to get subset
94 int rank = 0, procs = 1;
95 #ifdef MOAB_HAVE_MPI
96 if( isParallel )
97 {
98 ParallelComm*& myPcomm = _readNC->myPcomm;
99 rank = myPcomm->proc_config().proc_rank();
100 procs = myPcomm->proc_config().proc_size();
101 }
102 #endif
103 if( procs > 1 )
104 {
105 for( int i = 0; i < 6; i++ )
106 parData.gDims[i] = gDims[i];
107 parData.partMethod = partMethod;
108 int pdims[3];
109
110 locallyPeriodic[0] = locallyPeriodic[1] = locallyPeriodic[2] = 0;
111 rval = ScdInterface::compute_partition( procs, rank, parData, lDims, locallyPeriodic, pdims );MB_CHK_ERR( rval );
112 for( int i = 0; i < 3; i++ )
113 parData.pDims[i] = pdims[i];
114
115 dbgOut.tprintf( 1, "Partition: %dx%d (out of %dx%d)\n", lDims[3] - lDims[0], lDims[4] - lDims[1],
116 gDims[3] - gDims[0], gDims[4] - gDims[1] );
117 if( 0 == rank )
118 dbgOut.tprintf( 1, "Contiguous chunks of size %d bytes.\n",
119 8 * ( lDims[3] - lDims[0] ) * ( lDims[4] - lDims[1] ) );
120 }
121 else
122 {
123 for( int i = 0; i < 6; i++ )
124 lDims[i] = gDims[i];
125 locallyPeriodic[0] = globallyPeriodic[0];
126 }
127
128 // Now get actual coordinate values for vertices and cell centers
129 lCDims[0] = lDims[0];
130
131 lCDims[3] = lDims[3];
132
133 // For FV models, will always be non-periodic in j
134 lCDims[1] = lDims[1];
135 lCDims[4] = lDims[4];
136
137 #if 0
138 // Resize vectors to store values later
139 if (-1 != lDims[0])
140 ilVals.resize(lDims[3] - lDims[0] + 1);
141 if (-1 != lCDims[0])
142 ilCVals.resize(lCDims[3] - lCDims[0] + 1);
143 if (-1 != lDims[1])
144 jlVals.resize(lDims[4] - lDims[1] + 1);
145 if (-1 != lCDims[1])
146 jlCVals.resize(lCDims[4] - lCDims[1] + 1);
147 if (nTimeSteps > 0)
148 tVals.resize(nTimeSteps);
149 #endif
150
151 dbgOut.tprintf( 1, "I=%d-%d, J=%d-%d\n", lDims[0], lDims[3], lDims[1], lDims[4] );
152 dbgOut.tprintf( 1, "%d elements, %d vertices\n", ( lDims[3] - lDims[0] ) * ( lDims[4] - lDims[1] ),
153 ( lDims[3] - lDims[0] ) * ( lDims[4] - lDims[1] ) * nv );
154
155 // For each variable, determine the entity location type and number of levels
156 std::map< std::string, ReadNC::VarData >::iterator mit;
157 for( mit = varInfo.begin(); mit != varInfo.end(); ++mit )
158 {
159 ReadNC::VarData& vd = ( *mit ).second;
160
161 // Default entLoc is ENTLOCSET
162 if( std::find( vd.varDims.begin(), vd.varDims.end(), tDim ) != vd.varDims.end() )
163 {
164 if( ( std::find( vd.varDims.begin(), vd.varDims.end(), iCDim ) != vd.varDims.end() ) &&
165 ( std::find( vd.varDims.begin(), vd.varDims.end(), jCDim ) != vd.varDims.end() ) )
166 vd.entLoc = ReadNC::ENTLOCFACE;
167 else if( ( std::find( vd.varDims.begin(), vd.varDims.end(), jDim ) != vd.varDims.end() ) &&
168 ( std::find( vd.varDims.begin(), vd.varDims.end(), iCDim ) != vd.varDims.end() ) )
169 vd.entLoc = ReadNC::ENTLOCNSEDGE;
170 else if( ( std::find( vd.varDims.begin(), vd.varDims.end(), jCDim ) != vd.varDims.end() ) &&
171 ( std::find( vd.varDims.begin(), vd.varDims.end(), iDim ) != vd.varDims.end() ) )
172 vd.entLoc = ReadNC::ENTLOCEWEDGE;
173 }
174
175 // Default numLev is 0
176 if( std::find( vd.varDims.begin(), vd.varDims.end(), levDim ) != vd.varDims.end() ) vd.numLev = nLevels;
177 }
178
179 std::vector< std::string > ijdimNames( 2 );
180 ijdimNames[0] = "__ni";
181 ijdimNames[1] = "__nj";
182 std::string tag_name;
183 Tag tagh;
184
185 // __<dim_name>_LOC_MINMAX (for slon, slat, lon and lat)
186 for( unsigned int i = 0; i != ijdimNames.size(); i++ )
187 {
188 std::vector< int > val( 2, 0 );
189 if( ijdimNames[i] == "__ni" )
190 {
191 val[0] = lDims[0];
192 val[1] = lDims[3];
193 }
194 else if( ijdimNames[i] == "__nj" )
195 {
196 val[0] = lDims[1];
197 val[1] = lDims[4];
198 }
199
200 std::stringstream ss_tag_name;
201 ss_tag_name << ijdimNames[i] << "_LOC_MINMAX";
202 tag_name = ss_tag_name.str();
203 rval = mbImpl->tag_get_handle( tag_name.c_str(), 2, MB_TYPE_INTEGER, tagh, MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Trouble creating conventional tag " << tag_name );
204 rval = mbImpl->tag_set_data( tagh, &_fileSet, 1, &val[0] );MB_CHK_SET_ERR( rval, "Trouble setting data to conventional tag " << tag_name );
205 if( MB_SUCCESS == rval ) dbgOut.tprintf( 2, "Conventional tag %s is created.\n", tag_name.c_str() );
206 }
207
208 // __<dim_name>_LOC_VALS (for slon, slat, lon and lat)
209 // Assume all have the same data type as lon (expected type is float or double)
210 switch( varInfo["xc"].varDataType )
211 {
212 case NC_FLOAT:
213 case NC_DOUBLE:
214 break;
215 default:
216 MB_SET_ERR( MB_FAILURE, "Unexpected variable data type for 'lon'" );
217 }
218
219 // do not need conventional tags
220 Tag convTagsCreated = 0;
221 int def_val = 0;
222 rval = mbImpl->tag_get_handle( "__CONV_TAGS_CREATED", 1, MB_TYPE_INTEGER, convTagsCreated,
223 MB_TAG_SPARSE | MB_TAG_CREAT, &def_val );MB_CHK_SET_ERR( rval, "Trouble getting _CONV_TAGS_CREATED tag" );
224 int create_conv_tags_flag = 1;
225 rval = mbImpl->tag_set_data( convTagsCreated, &_fileSet, 1, &create_conv_tags_flag );MB_CHK_SET_ERR( rval, "Trouble setting _CONV_TAGS_CREATED tag" );
226
227 return MB_SUCCESS;
228 }
References moab::NCHelper::_fileSet, moab::NCHelper::_readNC, moab::ScdInterface::compute_partition(), moab::ReadNC::dbgOut, moab::ReadNC::dimLens, moab::ReadNC::dimNames, moab::ReadNC::VarData::entLoc, moab::ReadNC::ENTLOCEWEDGE, moab::ReadNC::ENTLOCFACE, moab::ReadNC::ENTLOCNSEDGE, ErrorCode, moab::ScdNCHelper::gDims, moab::ScdParData::gDims, moab::ScdNCHelper::globallyPeriodic, moab::ScdNCHelper::iCDim, moab::ScdNCHelper::iDim, moab::ScdNCHelper::ilCVals, moab::ScdNCHelper::ilVals, moab::ReadNC::isParallel, moab::ScdNCHelper::jCDim, moab::ScdNCHelper::jDim, moab::ScdNCHelper::jlCVals, moab::ScdNCHelper::jlVals, moab::ScdNCHelper::lCDims, moab::ScdNCHelper::lDims, moab::NCHelper::levDim, moab::ScdNCHelper::locallyPeriodic, MB_CHK_ERR, MB_CHK_SET_ERR, MB_SET_ERR, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_SPARSE, MB_TYPE_INTEGER, moab::ReadNC::mbImpl, moab::NCHelper::nLevels, moab::NCHelper::nTimeSteps, moab::ReadNC::VarData::numLev, nv, nvDim, moab::ReadNC::parData, moab::ReadNC::partMethod, moab::ScdParData::partMethod, moab::ScdParData::pDims, moab::ParallelComm::proc_config(), moab::ProcConfig::proc_rank(), moab::ProcConfig::proc_size(), moab::Interface::tag_get_handle(), moab::Interface::tag_set_data(), moab::NCHelper::tDim, moab::DebugOutput::tprintf(), moab::NCHelper::tVals, moab::ReadNC::VarData::varDims, and moab::ReadNC::varInfo.
|
private |
Definition at line 48 of file NCHelperDomain.hpp.
Referenced by create_mesh(), and init_mesh_vals().
|
private |
Definition at line 49 of file NCHelperDomain.hpp.
Referenced by init_mesh_vals().