moab
|
Child helper class for MPAS grid. More...
#include <NCHelperMPAS.hpp>
Public Member Functions | |
NCHelperMPAS (ReadNC *readNC, int fileId, const FileOptions &opts, EntityHandle fileSet) | |
Static Public Member Functions | |
static bool | can_read_file (ReadNC *readNC) |
Private Member Functions | |
virtual ErrorCode | init_mesh_vals () |
Implementation of NCHelper::init_mesh_vals() | |
virtual ErrorCode | check_existing_mesh () |
Implementation of NCHelper::check_existing_mesh() | |
virtual ErrorCode | create_mesh (Range &faces) |
Implementation of NCHelper::create_mesh() | |
virtual std::string | get_mesh_type_name () |
Implementation of NCHelper::get_mesh_type_name() | |
virtual ErrorCode | read_ucd_variable_to_nonset_allocate (std::vector< ReadNC::VarData > &vdatas, std::vector< int > &tstep_nums) |
Implementation of UcdNCHelper::read_ucd_variable_to_nonset_allocate() | |
virtual ErrorCode | read_ucd_variable_to_nonset (std::vector< ReadNC::VarData > &vdatas, std::vector< int > &tstep_nums) |
Implementation of UcdNCHelper::read_ucd_variable_to_nonset() | |
ErrorCode | redistribute_local_cells (int start_cell_index) |
Redistribute local cells after trivial partition (e.g. Zoltan partition, if applicable) | |
ErrorCode | create_local_vertices (const std::vector< int > &vertices_on_local_cells, EntityHandle &start_vertex) |
Create local vertices. | |
ErrorCode | create_local_edges (EntityHandle start_vertex, const std::vector< int > &num_edges_on_local_cells) |
Create local edges (optional) | |
ErrorCode | create_local_cells (const std::vector< int > &vertices_on_local_cells, const std::vector< int > &num_edges_on_local_cells, EntityHandle start_vertex, Range &faces) |
Create local cells without padding (cells are divided into groups based on the number of edges) | |
ErrorCode | create_padded_local_cells (const std::vector< int > &vertices_on_local_cells, EntityHandle start_vertex, Range &faces) |
Create local cells with padding (padded cells will have the same number of edges) | |
ErrorCode | create_gather_set_vertices (EntityHandle gather_set, EntityHandle &gather_set_start_vertex) |
Create gather set vertices. | |
ErrorCode | create_gather_set_edges (EntityHandle gather_set, EntityHandle gather_set_start_vertex) |
Create gather set edges (optional) | |
ErrorCode | create_gather_set_cells (EntityHandle gather_set, EntityHandle gather_set_start_vertex) |
Create gather set cells without padding (cells are divided into groups based on the number of edges) | |
ErrorCode | create_padded_gather_set_cells (EntityHandle gather_set, EntityHandle gather_set_start_vertex) |
Create gather set cells with padding (padded cells will have the same number of edges) | |
Private Attributes | |
int | maxEdgesPerCell |
int | numCellGroups |
bool | createGatherSet |
std::map< EntityHandle, int > | cellHandleToGlobalID |
Range | facesOwned |
Child helper class for MPAS grid.
Definition at line 17 of file NCHelperMPAS.hpp.
moab::NCHelperMPAS::NCHelperMPAS | ( | ReadNC * | readNC, |
int | fileId, | ||
const FileOptions & | opts, | ||
EntityHandle | fileSet | ||
) |
Definition at line 23 of file NCHelperMPAS.cpp.
: UcdNCHelper(readNC, fileId, opts, fileSet) , maxEdgesPerCell(DEFAULT_MAX_EDGES_PER_CELL) , numCellGroups(0) , createGatherSet(false) { // Ignore variables containing topological information ignoredVarNames.insert("nEdgesOnEdge"); ignoredVarNames.insert("nEdgesOnCell"); ignoredVarNames.insert("edgesOnVertex"); ignoredVarNames.insert("cellsOnVertex"); ignoredVarNames.insert("verticesOnEdge"); ignoredVarNames.insert("edgesOnEdge"); ignoredVarNames.insert("cellsOnEdge"); ignoredVarNames.insert("verticesOnCell"); ignoredVarNames.insert("edgesOnCell"); ignoredVarNames.insert("cellsOnCell"); // Ignore variables for index conversion ignoredVarNames.insert("indexToVertexID"); ignoredVarNames.insert("indexToEdgeID"); ignoredVarNames.insert("indexToCellID"); }
bool moab::NCHelperMPAS::can_read_file | ( | ReadNC * | readNC | ) | [static] |
Definition at line 47 of file NCHelperMPAS.cpp.
{ std::vector<std::string>& dimNames = readNC->dimNames; // If dimension name "vertexDegree" exists then it should be the MPAS grid if (std::find(dimNames.begin(), dimNames.end(), std::string("vertexDegree")) != dimNames.end()) return true; return false; }
ErrorCode moab::NCHelperMPAS::check_existing_mesh | ( | ) | [private, virtual] |
Implementation of NCHelper::check_existing_mesh()
Implements moab::NCHelper.
Definition at line 206 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; Tag& mGlobalIdTag = _readNC->mGlobalIdTag; bool& noMesh = _readNC->noMesh; if (noMesh) { ErrorCode rval; // Restore numCellGroups if (0 == numCellGroups) { Tag numCellGroupsTag; rval = mbImpl->tag_get_handle("__NUM_CELL_GROUPS", 1, MB_TYPE_INTEGER, numCellGroupsTag); if (MB_SUCCESS == rval) rval = mbImpl->tag_get_data(numCellGroupsTag, &_fileSet, 1, &numCellGroups); } if (localGidVerts.empty()) { // Get all vertices from tmp_set (it is the input set in no_mesh scenario) Range local_verts; rval = mbImpl->get_entities_by_dimension(_fileSet, 0, local_verts); if (MB_FAILURE == rval) return rval; if (!local_verts.empty()) { std::vector<int> gids(local_verts.size()); // !IMPORTANT : this has to be the GLOBAL_ID tag rval = mbImpl->tag_get_data(mGlobalIdTag, local_verts, &gids[0]); if (MB_FAILURE == rval) return rval; // Restore localGidVerts std::copy(gids.rbegin(), gids.rend(), range_inserter(localGidVerts)); nLocalVertices = localGidVerts.size(); } } if (localGidEdges.empty()) { // Get all edges from _fileSet (it is the input set in no_mesh scenario) Range local_edges; rval = mbImpl->get_entities_by_dimension(_fileSet, 1, local_edges); if (MB_FAILURE == rval) return rval; if (!local_edges.empty()) { std::vector<int> gids(local_edges.size()); // !IMPORTANT : this has to be the GLOBAL_ID tag rval = mbImpl->tag_get_data(mGlobalIdTag, local_edges, &gids[0]); if (MB_FAILURE == rval) return rval; // Restore localGidEdges std::copy(gids.rbegin(), gids.rend(), range_inserter(localGidEdges)); nLocalEdges = localGidEdges.size(); } } if (localGidCells.empty()) { // Get all cells from tmp_set (it is the input set in no_mesh scenario) Range local_cells; rval = mbImpl->get_entities_by_dimension(_fileSet, 2, local_cells); if (MB_FAILURE == rval) return rval; if (!local_cells.empty()) { std::vector<int> gids(local_cells.size()); // !IMPORTANT : this has to be the GLOBAL_ID tag rval = mbImpl->tag_get_data(mGlobalIdTag, local_cells, &gids[0]); if (MB_FAILURE == rval) return rval; // Restore localGidCells std::copy(gids.rbegin(), gids.rend(), range_inserter(localGidCells)); nLocalCells = localGidCells.size(); if (numCellGroups > 1) { // Restore cellHandleToGlobalID map Range::const_iterator rit; int i; for (rit = local_cells.begin(), i = 0; rit != local_cells.end(); ++rit, i++) cellHandleToGlobalID[*rit] = gids[i]; } } } } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_gather_set_cells | ( | EntityHandle | gather_set, |
EntityHandle | gather_set_start_vertex | ||
) | [private] |
Create gather set cells without padding (cells are divided into groups based on the number of edges)
Definition at line 1582 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; // Read number of edges on each gather set cell int nEdgesOnCellVarId; int success = NCFUNC(inq_varid)(_fileId, "nEdgesOnCell", &nEdgesOnCellVarId); ERRORS(success, "Failed to get variable id of nEdgesOnCell."); std::vector<int> num_edges_on_gather_set_cells(nCells); NCDF_SIZE read_start = 0; NCDF_SIZE read_count = static_cast<NCDF_SIZE>(nCells); #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_int)(_fileId, nEdgesOnCellVarId, &read_start, &read_count, &num_edges_on_gather_set_cells[0]); ERRORS(success, "Failed to read nEdgesOnCell data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_int)(_fileId, nEdgesOnCellVarId, &read_start, &read_count, &num_edges_on_gather_set_cells[0]); ERRORS(success, "Failed to read nEdgesOnCell data."); #endif // Read vertices on each gather set cell (connectivity) int verticesOnCellVarId; success = NCFUNC(inq_varid)(_fileId, "verticesOnCell", &verticesOnCellVarId); ERRORS(success, "Failed to get variable id of verticesOnCell."); std::vector<int> vertices_on_gather_set_cells(nCells * maxEdgesPerCell); NCDF_SIZE read_starts[2] = {0, 0}; NCDF_SIZE read_counts[2] = {static_cast<NCDF_SIZE>(nCells), static_cast<NCDF_SIZE>(maxEdgesPerCell)}; #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_int)(_fileId, verticesOnCellVarId, read_starts, read_counts, &vertices_on_gather_set_cells[0]); ERRORS(success, "Failed to read verticesOnCell data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_int)(_fileId, verticesOnCellVarId, read_starts, read_counts, &vertices_on_gather_set_cells[0]); ERRORS(success, "Failed to read verticesOnCell data."); #endif // Divide gather set cells into groups based on the number of edges Range gather_set_cells_with_n_edges[DEFAULT_MAX_EDGES_PER_CELL + 1]; // Insert larger values before smaller values to increase efficiency for (int i = nCells - 1; i >= 0; i--) { int num_edges = num_edges_on_gather_set_cells[i]; gather_set_cells_with_n_edges[num_edges].insert(i + 1); // 0 based -> 1 based } // Create gather set cells EntityHandle* conn_arr_gather_set_cells_with_n_edges[DEFAULT_MAX_EDGES_PER_CELL + 1]; for (int num_edges_per_cell = 3; num_edges_per_cell <= maxEdgesPerCell; num_edges_per_cell++) { int num_group_cells = gather_set_cells_with_n_edges[num_edges_per_cell].size(); if (num_group_cells > 0) { EntityHandle start_element; ErrorCode rval = _readNC->readMeshIface->get_element_connect(num_group_cells, num_edges_per_cell, MBPOLYGON, 0, start_element, conn_arr_gather_set_cells_with_n_edges[num_edges_per_cell], num_group_cells); ERRORR(rval, "Failed to create cells."); // Add cells to the gather set Range gather_set_cells_range(start_element, start_element + num_group_cells - 1); rval = mbImpl->add_entities(gather_set, gather_set_cells_range); ERRORR(rval, "Failed to add cells to the gather set."); for (int j = 0; j < num_group_cells; j++) { int gather_set_cell_idx = gather_set_cells_with_n_edges[num_edges_per_cell][j]; // Global cell index, 1 based gather_set_cell_idx--; // 1 based -> 0 based for (int k = 0; k < num_edges_per_cell; k++) { EntityHandle gather_set_vert_idx = vertices_on_gather_set_cells[gather_set_cell_idx * maxEdgesPerCell + k]; // Global vertex index, 1 based gather_set_vert_idx--; // 1 based -> 0 based // Connectivity array is shifted by where the gather set vertices start conn_arr_gather_set_cells_with_n_edges[num_edges_per_cell][j * num_edges_per_cell + k] = gather_set_start_vertex + gather_set_vert_idx; } } } } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_gather_set_edges | ( | EntityHandle | gather_set, |
EntityHandle | gather_set_start_vertex | ||
) | [private] |
Create gather set edges (optional)
Definition at line 1532 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; // Create gather set edges EntityHandle start_edge; EntityHandle* conn_arr_gather_set_edges = NULL; // Don't need to specify allocation number here, because we know enough edges were created before ErrorCode rval = _readNC->readMeshIface->get_element_connect(nEdges, 2, MBEDGE, 0, start_edge, conn_arr_gather_set_edges); ERRORR(rval, "Failed to create edges."); // Add edges to the gather set Range gather_set_edges_range(start_edge, start_edge + nEdges - 1); rval = mbImpl->add_entities(gather_set, gather_set_edges_range); ERRORR(rval, "Failed to add edges to the gather set."); // Read vertices on each edge int verticesOnEdgeVarId; int success = NCFUNC(inq_varid)(_fileId, "verticesOnEdge", &verticesOnEdgeVarId); ERRORS(success, "Failed to get variable id of verticesOnEdge."); // Utilize the memory storage pointed by conn_arr_gather_set_edges int* vertices_on_gather_set_edges = (int*) conn_arr_gather_set_edges; NCDF_SIZE read_starts[2] = {0, 0}; NCDF_SIZE read_counts[2] = {static_cast<NCDF_SIZE>(nEdges), 2}; #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_int)(_fileId, verticesOnEdgeVarId, read_starts, read_counts, vertices_on_gather_set_edges); ERRORS(success, "Failed to read verticesOnEdge data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_int)(_fileId, verticesOnEdgeVarId, read_starts, read_counts, vertices_on_gather_set_edges); ERRORS(success, "Failed to read verticesOnEdge data."); #endif // Populate connectivity data for gather set edges // Convert in-place from int (stored in the first half) to EntityHandle // Reading backward is the trick for (int edge_vert = nEdges * 2 - 1; edge_vert >= 0; edge_vert--) { int gather_set_vert_idx = vertices_on_gather_set_edges[edge_vert]; // Global vertex index, 1 based gather_set_vert_idx--; // 1 based -> 0 based // Connectivity array is shifted by where the gather set vertices start conn_arr_gather_set_edges[edge_vert] = gather_set_start_vertex + gather_set_vert_idx; } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_gather_set_vertices | ( | EntityHandle | gather_set, |
EntityHandle & | gather_set_start_vertex | ||
) | [private] |
Create gather set vertices.
Definition at line 1436 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; Tag& mGlobalIdTag = _readNC->mGlobalIdTag; const Tag*& mpFileIdTag = _readNC->mpFileIdTag; // Create gather set vertices std::vector<double*> arrays; // Don't need to specify allocation number here, because we know enough vertices were created before ErrorCode rval = _readNC->readMeshIface->get_node_coords(3, nVertices, 0, gather_set_start_vertex, arrays); ERRORR(rval, "Failed to create vertices."); // Add vertices to the gather set Range gather_set_verts_range(gather_set_start_vertex, gather_set_start_vertex + nVertices - 1); rval = mbImpl->add_entities(gather_set, gather_set_verts_range); ERRORR(rval, "Failed to add vertices to the gather set."); // Read x coordinates for gather set vertices double* xptr = arrays[0]; int xVertexVarId; int success = NCFUNC(inq_varid)(_fileId, "xVertex", &xVertexVarId); ERRORS(success, "Failed to get variable id of xVertex."); NCDF_SIZE read_start = 0; NCDF_SIZE read_count = static_cast<NCDF_SIZE>(nVertices); #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_double)(_fileId, xVertexVarId, &read_start, &read_count, xptr); ERRORS(success, "Failed to read xVertex data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_double)(_fileId, xVertexVarId, &read_start, &read_count, xptr); ERRORS(success, "Failed to read xVertex data."); #endif // Read y coordinates for gather set vertices double* yptr = arrays[1]; int yVertexVarId; success = NCFUNC(inq_varid)(_fileId, "yVertex", &yVertexVarId); ERRORS(success, "Failed to get variable id of yVertex."); #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_double)(_fileId, yVertexVarId, &read_start, &read_count, yptr); ERRORS(success, "Failed to read yVertex data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_double)(_fileId, yVertexVarId, &read_start, &read_count, yptr); ERRORS(success, "Failed to read yVertex data."); #endif // Read z coordinates for gather set vertices double* zptr = arrays[2]; int zVertexVarId; success = NCFUNC(inq_varid)(_fileId, "zVertex", &zVertexVarId); ERRORS(success, "Failed to get variable id of zVertex."); #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_double)(_fileId, zVertexVarId, &read_start, &read_count, zptr); ERRORS(success, "Failed to read zVertex data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_double)(_fileId, zVertexVarId, &read_start, &read_count, zptr); ERRORS(success, "Failed to read zVertex data."); #endif // Get ptr to GID memory for gather set vertices int count = 0; void* data = NULL; rval = mbImpl->tag_iterate(mGlobalIdTag, gather_set_verts_range.begin(), gather_set_verts_range.end(), count, data); ERRORR(rval, "Failed to iterate global id tag on gather set vertices."); assert(count == nVertices); int* gid_data = (int*) data; for (int j = 1; j <= nVertices; j++) gid_data[j - 1] = j; // Set the file id tag too, it should be bigger something not interfering with global id if (mpFileIdTag) { rval = mbImpl->tag_iterate(*mpFileIdTag, gather_set_verts_range.begin(), gather_set_verts_range.end(), count, data); ERRORR(rval, "Failed to iterate file id tag on gather set vertices."); assert(count == nVertices); gid_data = (int*) data; for (int j = 1; j <= nVertices; j++) gid_data[j - 1] = nVertices + j; // Bigger than global id tag } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_local_cells | ( | const std::vector< int > & | vertices_on_local_cells, |
const std::vector< int > & | num_edges_on_local_cells, | ||
EntityHandle | start_vertex, | ||
Range & | faces | ||
) | [private] |
Create local cells without padding (cells are divided into groups based on the number of edges)
Definition at line 1317 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; Tag& mGlobalIdTag = _readNC->mGlobalIdTag; // Divide local cells into groups based on the number of edges Range local_cells_with_n_edges[DEFAULT_MAX_EDGES_PER_CELL + 1]; // Insert larger values before smaller ones to increase efficiency for (int i = nLocalCells - 1; i >= 0; i--) { int num_edges = num_edges_on_local_cells[i]; local_cells_with_n_edges[num_edges].insert(localGidCells[i]); // Global cell index } std::vector<int> num_edges_on_cell_groups; for (int i = 3; i <= maxEdgesPerCell; i++) { if (local_cells_with_n_edges[i].size() > 0) num_edges_on_cell_groups.push_back(i); } numCellGroups = num_edges_on_cell_groups.size(); EntityHandle* conn_arr_local_cells_with_n_edges[DEFAULT_MAX_EDGES_PER_CELL + 1]; for (int i = 0; i < numCellGroups; i++) { int num_edges_per_cell = num_edges_on_cell_groups[i]; int num_group_cells = (int)local_cells_with_n_edges[num_edges_per_cell].size(); // Create local cells for each non-empty cell group EntityHandle start_element; ErrorCode rval = _readNC->readMeshIface->get_element_connect(num_group_cells, num_edges_per_cell, MBPOLYGON, 0, start_element, conn_arr_local_cells_with_n_edges[num_edges_per_cell], num_group_cells); ERRORR(rval, "Failed to create cells"); faces.insert(start_element, start_element + num_group_cells - 1); // Add local cells to the file set Range local_cells_range(start_element, start_element + num_group_cells - 1); rval = _readNC->mbImpl->add_entities(_fileSet, local_cells_range); ERRORR(rval, "Failed to add local cells to the file set."); // Get ptr to gid memory for local cells int count = 0; void* data = NULL; rval = mbImpl->tag_iterate(mGlobalIdTag, local_cells_range.begin(), local_cells_range.end(), count, data); ERRORR(rval, "Failed to iterate global id tag on local cells."); assert(count == num_group_cells); int* gid_data = (int*) data; std::copy(local_cells_with_n_edges[num_edges_per_cell].begin(), local_cells_with_n_edges[num_edges_per_cell].end(), gid_data); // Set connectivity array with proper local vertices handles for (int j = 0; j < num_group_cells; j++) { EntityHandle global_cell_idx = local_cells_with_n_edges[num_edges_per_cell][j]; // Global cell index, 1 based int local_cell_idx = localGidCells.index(global_cell_idx); // Local cell index, 0 based assert(local_cell_idx != -1); if (numCellGroups > 1) { // Populate cellHandleToGlobalID map to read cell variables cellHandleToGlobalID[start_element + j] = global_cell_idx; } for (int k = 0; k < num_edges_per_cell; k++) { EntityHandle global_vert_idx = vertices_on_local_cells[local_cell_idx * maxEdgesPerCell + k]; // Global vertex index, 1 based int local_vert_idx = localGidVerts.index(global_vert_idx); // Local vertex index, 0 based assert(local_vert_idx != -1); conn_arr_local_cells_with_n_edges[num_edges_per_cell][j * num_edges_per_cell + k] = start_vertex + local_vert_idx; } } } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_local_edges | ( | EntityHandle | start_vertex, |
const std::vector< int > & | num_edges_on_local_cells | ||
) | [private] |
Create local edges (optional)
Definition at line 1173 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; Tag& mGlobalIdTag = _readNC->mGlobalIdTag; DebugOutput& dbgOut = _readNC->dbgOut; // Read edges on each local cell, to get localGidEdges int edgesOnCellVarId; int success = NCFUNC(inq_varid)(_fileId, "edgesOnCell", &edgesOnCellVarId); ERRORS(success, "Failed to get variable id of edgesOnCell."); std::vector<int> edges_on_local_cells(nLocalCells * maxEdgesPerCell); dbgOut.tprintf(1, " edges_on_local_cells.size() = %d\n", (int)edges_on_local_cells.size()); #ifdef PNETCDF_FILE size_t nb_reads = localGidCells.psize(); std::vector<int> requests(nb_reads); std::vector<int> statuss(nb_reads); size_t idxReq = 0; #endif size_t indexInArray = 0; for (Range::pair_iterator pair_iter = localGidCells.pair_begin(); pair_iter != localGidCells.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_starts[2] = {static_cast<NCDF_SIZE>(starth - 1), 0}; NCDF_SIZE read_counts[2] = {static_cast<NCDF_SIZE>(endh - starth + 1), static_cast<NCDF_SIZE>(maxEdgesPerCell)}; // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_int)(_fileId, edgesOnCellVarId, read_starts, read_counts, &(edges_on_local_cells[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_int)(_fileId, edgesOnCellVarId, read_starts, read_counts, &(edges_on_local_cells[indexInArray])); #endif ERRORS(success, "Failed to read edgesOnCell data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1) * maxEdgesPerCell; } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif // Correct local cell edges array in the same way as local cell vertices array, replace the // padded edges with the last edges in the corresponding cells for (int local_cell_idx = 0; local_cell_idx < nLocalCells; local_cell_idx++) { int num_edges = num_edges_on_local_cells[local_cell_idx]; int idx_in_local_edge_arr = local_cell_idx * maxEdgesPerCell; int last_edge_idx = edges_on_local_cells[idx_in_local_edge_arr + num_edges - 1]; for (int i = num_edges; i < maxEdgesPerCell; i++) edges_on_local_cells[idx_in_local_edge_arr + i] = last_edge_idx; } // Collect local edges std::sort(edges_on_local_cells.begin(), edges_on_local_cells.end()); std::copy(edges_on_local_cells.rbegin(), edges_on_local_cells.rend(), range_inserter(localGidEdges)); nLocalEdges = localGidEdges.size(); dbgOut.tprintf(1, " localGidEdges.psize() = %d\n", (int)localGidEdges.psize()); dbgOut.tprintf(1, " localGidEdges.size() = %d\n", (int)localGidEdges.size()); // Create local edges EntityHandle start_edge; EntityHandle* conn_arr_edges = NULL; ErrorCode rval = _readNC->readMeshIface->get_element_connect(nLocalEdges, 2, MBEDGE, 0, start_edge, conn_arr_edges, // Might have to create gather mesh later (createGatherSet ? nLocalEdges + nEdges : nLocalEdges)); ERRORR(rval, "Failed to create edges."); // Add local edges to the file set Range local_edges_range(start_edge, start_edge + nLocalEdges - 1); rval = _readNC->mbImpl->add_entities(_fileSet, local_edges_range); ERRORR(rval, "Failed to add local edges to the file set."); // Get ptr to GID memory for edges int count = 0; void* data = NULL; rval = mbImpl->tag_iterate(mGlobalIdTag, local_edges_range.begin(), local_edges_range.end(), count, data); ERRORR(rval, "Failed to iterate global id tag on local edges."); assert(count == nLocalEdges); int* gid_data = (int*) data; std::copy(localGidEdges.begin(), localGidEdges.end(), gid_data); int verticesOnEdgeVarId; // Read vertices on each local edge, to get edge connectivity success = NCFUNC(inq_varid)(_fileId, "verticesOnEdge", &verticesOnEdgeVarId); ERRORS(success, "Failed to get variable id of verticesOnEdge."); // Utilize the memory storage pointed by conn_arr_edges int* vertices_on_local_edges = (int*) conn_arr_edges; #ifdef PNETCDF_FILE nb_reads = localGidEdges.psize(); requests.resize(nb_reads); statuss.resize(nb_reads); idxReq = 0; #endif indexInArray = 0; for (Range::pair_iterator pair_iter = localGidEdges.pair_begin(); pair_iter != localGidEdges.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_starts[2] = {static_cast<NCDF_SIZE>(starth - 1), 0}; NCDF_SIZE read_counts[2] = {static_cast<NCDF_SIZE>(endh - starth + 1), 2}; // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_int)(_fileId, verticesOnEdgeVarId, read_starts, read_counts, &(vertices_on_local_edges[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_int)(_fileId, verticesOnEdgeVarId, read_starts, read_counts, &(vertices_on_local_edges[indexInArray])); #endif ERRORS(success, "Failed to read verticesOnEdge data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1) * 2; } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif // Populate connectivity data for local edges // Convert in-place from int (stored in the first half) to EntityHandle // Reading backward is the trick for (int edge_vert = nLocalEdges * 2 - 1; edge_vert >= 0; edge_vert--) { int global_vert_idx = vertices_on_local_edges[edge_vert]; // Global vertex index, 1 based int local_vert_idx = localGidVerts.index(global_vert_idx); // Local vertex index, 0 based assert(local_vert_idx != -1); conn_arr_edges[edge_vert] = start_vertex + local_vert_idx; } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_local_vertices | ( | const std::vector< int > & | vertices_on_local_cells, |
EntityHandle & | start_vertex | ||
) | [private] |
Create local vertices.
Definition at line 1009 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; Tag& mGlobalIdTag = _readNC->mGlobalIdTag; const Tag*& mpFileIdTag = _readNC->mpFileIdTag; DebugOutput& dbgOut = _readNC->dbgOut; // Make a copy of vertices_on_local_cells for sorting (keep original one to set cell connectivity later) std::vector<int> vertices_on_local_cells_sorted(vertices_on_local_cells); std::sort(vertices_on_local_cells_sorted.begin(), vertices_on_local_cells_sorted.end()); std::copy(vertices_on_local_cells_sorted.rbegin(), vertices_on_local_cells_sorted.rend(), range_inserter(localGidVerts)); nLocalVertices = localGidVerts.size(); dbgOut.tprintf(1, " localGidVerts.psize() = %d\n", (int)localGidVerts.psize()); dbgOut.tprintf(1, " localGidVerts.size() = %d\n", (int)localGidVerts.size()); // Create local vertices std::vector<double*> arrays; ErrorCode rval = _readNC->readMeshIface->get_node_coords(3, nLocalVertices, 0, start_vertex, arrays, // Might have to create gather mesh later (createGatherSet ? nLocalVertices + nVertices : nLocalVertices)); ERRORR(rval, "Failed to create local vertices."); // Add local vertices to the file set Range local_verts_range(start_vertex, start_vertex + nLocalVertices - 1); rval = _readNC->mbImpl->add_entities(_fileSet, local_verts_range); ERRORR(rval, "Failed to add local vertices to the file set."); // Get ptr to GID memory for local vertices int count = 0; void* data = NULL; rval = mbImpl->tag_iterate(mGlobalIdTag, local_verts_range.begin(), local_verts_range.end(), count, data); ERRORR(rval, "Failed to iterate global id tag on local vertices."); assert(count == nLocalVertices); int* gid_data = (int*) data; std::copy(localGidVerts.begin(), localGidVerts.end(), gid_data); // Duplicate GID data, which will be used to resolve sharing if (mpFileIdTag) { rval = mbImpl->tag_iterate(*mpFileIdTag, local_verts_range.begin(), local_verts_range.end(), count, data); ERRORR(rval, "Failed to iterate file id tag on local vertices."); assert(count == nLocalVertices); gid_data = (int*) data; std::copy(localGidVerts.begin(), localGidVerts.end(), gid_data); } #ifdef PNETCDF_FILE size_t nb_reads = localGidVerts.psize(); std::vector<int> requests(nb_reads); std::vector<int> statuss(nb_reads); size_t idxReq = 0; #endif // Read x coordinates for local vertices double* xptr = arrays[0]; int xVertexVarId; int success = NCFUNC(inq_varid)(_fileId, "xVertex", &xVertexVarId); ERRORS(success, "Failed to get variable id of xVertex."); size_t indexInArray = 0; for (Range::pair_iterator pair_iter = localGidVerts.pair_begin(); pair_iter != localGidVerts.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_start = (NCDF_SIZE) (starth - 1); NCDF_SIZE read_count = (NCDF_SIZE) (endh - starth + 1); // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_double)(_fileId, xVertexVarId, &read_start, &read_count, &(xptr[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_double)(_fileId, xVertexVarId, &read_start, &read_count, &(xptr[indexInArray])); #endif ERRORS(success, "Failed to read xVertex data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1); } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif // Read y coordinates for local vertices double* yptr = arrays[1]; int yVertexVarId; success = NCFUNC(inq_varid)(_fileId, "yVertex", &yVertexVarId); ERRORS(success, "Failed to get variable id of yVertex."); #ifdef PNETCDF_FILE idxReq = 0; #endif indexInArray = 0; for (Range::pair_iterator pair_iter = localGidVerts.pair_begin(); pair_iter != localGidVerts.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_start = (NCDF_SIZE) (starth - 1); NCDF_SIZE read_count = (NCDF_SIZE) (endh - starth + 1); // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_double)(_fileId, yVertexVarId, &read_start, &read_count, &(yptr[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_double)(_fileId, yVertexVarId, &read_start, &read_count, &(yptr[indexInArray])); #endif ERRORS(success, "Failed to read yVertex data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1); } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif // Read z coordinates for local vertices double* zptr = arrays[2]; int zVertexVarId; success = NCFUNC(inq_varid)(_fileId, "zVertex", &zVertexVarId); ERRORS(success, "Failed to get variable id of zVertex."); #ifdef PNETCDF_FILE idxReq = 0; #endif indexInArray = 0; for (Range::pair_iterator pair_iter = localGidVerts.pair_begin(); pair_iter != localGidVerts.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_start = (NCDF_SIZE) (starth - 1); NCDF_SIZE read_count = (NCDF_SIZE) (endh - starth + 1); // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_double)(_fileId, zVertexVarId, &read_start, &read_count, &(zptr[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_double)(_fileId, zVertexVarId, &read_start, &read_count, &(zptr[indexInArray])); #endif ERRORS(success, "Failed to read zVertex data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1); } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_mesh | ( | Range & | faces | ) | [private, virtual] |
Implementation of NCHelper::create_mesh()
Implements moab::NCHelper.
Definition at line 298 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; int& gatherSetRank = _readNC->gatherSetRank; bool& noMixedElements = _readNC->noMixedElements; bool& noEdges = _readNC->noEdges; int rank = 0; int procs = 1; #ifdef USE_MPI bool& isParallel = _readNC->isParallel; if (isParallel) { ParallelComm*& myPcomm = _readNC->myPcomm; rank = myPcomm->proc_config().proc_rank(); procs = myPcomm->proc_config().proc_size(); } #endif // Need to know whether we'll be creating gather mesh if (rank == gatherSetRank) createGatherSet = true; if (procs >= 2) { // Compute the number of local cells on this proc nLocalCells = int(std::floor(1.0 * nCells / procs)); // The starting global cell index in the MPAS file for this proc int start_cell_idx = rank * nLocalCells; // Number of extra cells after equal split over procs int iextra = nCells % procs; // Allocate extra cells over procs if (rank < iextra) nLocalCells++; start_cell_idx += std::min(rank, iextra); start_cell_idx++; // 0 based -> 1 based // Redistribute local cells after trivial partition (e.g. apply Zoltan partition) ErrorCode rval = redistribute_local_cells(start_cell_idx); ERRORR(rval, "Failed to redistribute local cells after trivial partition."); } else { nLocalCells = nCells; localGidCells.insert(1, nLocalCells); } // Read number of edges on each local cell, to calculate actual maxEdgesPerCell int nEdgesOnCellVarId; int success = NCFUNC(inq_varid)(_fileId, "nEdgesOnCell", &nEdgesOnCellVarId); ERRORS(success, "Failed to get variable id of nEdgesOnCell."); std::vector<int> num_edges_on_local_cells(nLocalCells); #ifdef PNETCDF_FILE size_t nb_reads = localGidCells.psize(); std::vector<int> requests(nb_reads); std::vector<int> statuss(nb_reads); size_t idxReq = 0; #endif size_t indexInArray = 0; for (Range::pair_iterator pair_iter = localGidCells.pair_begin(); pair_iter != localGidCells.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_start = (NCDF_SIZE) (starth - 1); NCDF_SIZE read_count = (NCDF_SIZE) (endh - starth + 1); // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_int)(_fileId, nEdgesOnCellVarId, &read_start, &read_count, &(num_edges_on_local_cells[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_int)(_fileId, nEdgesOnCellVarId, &read_start, &read_count, &(num_edges_on_local_cells[indexInArray])); #endif ERRORS(success, "Failed to read nEdgesOnCell data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1); } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif // Get local maxEdgesPerCell on this proc int local_max_edges_per_cell = *(std::max_element(num_edges_on_local_cells.begin(), num_edges_on_local_cells.end())); maxEdgesPerCell = local_max_edges_per_cell; // If parallel, do a MPI_Allreduce to get global maxEdgesPerCell across all procs #ifdef USE_MPI if (procs > 1) { int global_max_edges_per_cell; ParallelComm*& myPcomm = _readNC->myPcomm; MPI_Allreduce(&local_max_edges_per_cell, &global_max_edges_per_cell, 1, MPI_INTEGER, MPI_MAX, myPcomm->proc_config().proc_comm()); assert(local_max_edges_per_cell <= global_max_edges_per_cell); maxEdgesPerCell = global_max_edges_per_cell; if (0 == rank) { DebugOutput& dbgOut = _readNC->dbgOut; dbgOut.tprintf(1, " global_max_edges_per_cell = %d\n", global_max_edges_per_cell); } } #endif // Read vertices on each local cell, to get localGidVerts and cell connectivity later int verticesOnCellVarId; success = NCFUNC(inq_varid)(_fileId, "verticesOnCell", &verticesOnCellVarId); ERRORS(success, "Failed to get variable id of verticesOnCell."); std::vector<int> vertices_on_local_cells(nLocalCells * maxEdgesPerCell); #ifdef PNETCDF_FILE idxReq = 0; #endif indexInArray = 0; for (Range::pair_iterator pair_iter = localGidCells.pair_begin(); pair_iter != localGidCells.pair_end(); pair_iter++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; NCDF_SIZE read_starts[2] = {static_cast<NCDF_SIZE>(starth - 1), 0}; NCDF_SIZE read_counts[2] = {static_cast<NCDF_SIZE>(endh - starth + 1), static_cast<NCDF_SIZE>(maxEdgesPerCell)}; // Do a partial read in each subrange #ifdef PNETCDF_FILE success = NCFUNCREQG(_vara_int)(_fileId, verticesOnCellVarId, read_starts, read_counts, &(vertices_on_local_cells[indexInArray]), &requests[idxReq++]); #else success = NCFUNCAG(_vara_int)(_fileId, verticesOnCellVarId, read_starts, read_counts, &(vertices_on_local_cells[indexInArray])); #endif ERRORS(success, "Failed to read verticesOnCell data in a loop"); // Increment the index for next subrange indexInArray += (endh - starth + 1) * maxEdgesPerCell; } #ifdef PNETCDF_FILE // Wait outside the loop success = NCFUNC(wait_all)(_fileId, requests.size(), &requests[0], &statuss[0]); ERRORS(success, "Failed on wait_all."); #endif // Correct local cell vertices array, replace the padded vertices with the last vertices // in the corresponding cells; sometimes the padded vertices are 0, sometimes a large // vertex id. Make sure they are consistent to our padded option for (int local_cell_idx = 0; local_cell_idx < nLocalCells; local_cell_idx++) { int num_edges = num_edges_on_local_cells[local_cell_idx]; int idx_in_local_vert_arr = local_cell_idx * maxEdgesPerCell; int last_vert_idx = vertices_on_local_cells[idx_in_local_vert_arr + num_edges - 1]; for (int i = num_edges; i < maxEdgesPerCell; i++) vertices_on_local_cells[idx_in_local_vert_arr + i] = last_vert_idx; } // Create local vertices EntityHandle start_vertex; ErrorCode rval = create_local_vertices(vertices_on_local_cells, start_vertex); ERRORR(rval, "Failed to create local vertices for MPAS mesh."); // Create local edges (unless NO_EDGES read option is set) if (!noEdges) { rval = create_local_edges(start_vertex, num_edges_on_local_cells); ERRORR(rval, "Failed to create local edges for MPAS mesh."); } // Create local cells, either unpadded or padded if (noMixedElements) { rval = create_padded_local_cells(vertices_on_local_cells, start_vertex, faces); ERRORR(rval, "Failed to create padded local cells for MPAS mesh."); } else { rval = create_local_cells(vertices_on_local_cells, num_edges_on_local_cells, start_vertex, faces); ERRORR(rval, "Failed to create local cells for MPAS mesh."); } // Set tag for numCellGroups Tag numCellGroupsTag = 0; rval = mbImpl->tag_get_handle("__NUM_CELL_GROUPS", 1, MB_TYPE_INTEGER, numCellGroupsTag, MB_TAG_SPARSE | MB_TAG_CREAT); ERRORR(rval, "Failed to get __NUM_CELL_GROUPS tag."); rval = mbImpl->tag_set_data(numCellGroupsTag, &_fileSet, 1, &numCellGroups); ERRORR(rval, "Failed to set data for __NUM_CELL_GROUPS tag."); if (createGatherSet) { EntityHandle gather_set; rval = _readNC->readMeshIface->create_gather_set(gather_set); ERRORR(rval, "Failed to create gather set."); // Create gather set vertices EntityHandle start_gather_set_vertex; rval = create_gather_set_vertices(gather_set, start_gather_set_vertex); ERRORR(rval, "Failed to create gather set vertices for MPAS mesh"); // Create gather set edges (unless NO_EDGES read option is set) if (!noEdges) { rval = create_gather_set_edges(gather_set, start_gather_set_vertex); ERRORR(rval, "Failed to create gather set edges for MPAS mesh."); } // Create gather set cells, either unpadded or padded if (noMixedElements) { rval = create_padded_gather_set_cells(gather_set, start_gather_set_vertex); ERRORR(rval, "Failed to create padded gather set cells for MPAS mesh."); } else { rval = create_gather_set_cells(gather_set, start_gather_set_vertex); ERRORR(rval, "Failed to create gather set cells for MPAS mesh."); } } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_padded_gather_set_cells | ( | EntityHandle | gather_set, |
EntityHandle | gather_set_start_vertex | ||
) | [private] |
Create gather set cells with padding (padded cells will have the same number of edges)
Definition at line 1668 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; // Read number of edges on each gather set cell int nEdgesOnCellVarId; int success = NCFUNC(inq_varid)(_fileId, "nEdgesOnCell", &nEdgesOnCellVarId); ERRORS(success, "Failed to get variable id of nEdgesOnCell."); std::vector<int> num_edges_on_gather_set_cells(nCells); NCDF_SIZE read_start = 0; NCDF_SIZE read_count = static_cast<NCDF_SIZE>(nCells); #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_int)(_fileId, nEdgesOnCellVarId, &read_start, &read_count, &num_edges_on_gather_set_cells[0]); ERRORS(success, "Failed to read nEdgesOnCell data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_int)(_fileId, nEdgesOnCellVarId, &read_start, &read_count, &num_edges_on_gather_set_cells[0]); ERRORS(success, "Failed to read nEdgesOnCell data."); #endif // Create gather set cells EntityHandle start_element; EntityHandle* conn_arr_gather_set_cells = NULL; // Don't need to specify allocation number here, because we know enough cells were created before ErrorCode rval = _readNC->readMeshIface->get_element_connect(nCells, maxEdgesPerCell, MBPOLYGON, 0, start_element, conn_arr_gather_set_cells); ERRORR(rval, "Failed to create cells."); // Add cells to the gather set Range gather_set_cells_range(start_element, start_element + nCells - 1); rval = mbImpl->add_entities(gather_set, gather_set_cells_range); ERRORR(rval, "Failed to add cells to the gather set."); // Read vertices on each gather set cell (connectivity) int verticesOnCellVarId; success = NCFUNC(inq_varid)(_fileId, "verticesOnCell", &verticesOnCellVarId); ERRORS(success, "Failed to get variable id of verticesOnCell."); // Utilize the memory storage pointed by conn_arr_gather_set_cells int* vertices_on_gather_set_cells = (int*) conn_arr_gather_set_cells; NCDF_SIZE read_starts[2] = {0, 0}; NCDF_SIZE read_counts[2] = {static_cast<NCDF_SIZE>(nCells), static_cast<NCDF_SIZE>(maxEdgesPerCell)}; #ifdef PNETCDF_FILE // Enter independent I/O mode, since this read is only for the gather processor success = NCFUNC(begin_indep_data)(_fileId); ERRORS(success, "Failed to begin independent I/O mode."); success = NCFUNCG(_vara_int)(_fileId, verticesOnCellVarId, read_starts, read_counts, vertices_on_gather_set_cells); ERRORS(success, "Failed to read verticesOnCell data."); success = NCFUNC(end_indep_data)(_fileId); ERRORS(success, "Failed to end independent I/O mode."); #else success = NCFUNCG(_vara_int)(_fileId, verticesOnCellVarId, read_starts, read_counts, vertices_on_gather_set_cells); ERRORS(success, "Failed to read verticesOnCell data."); #endif // Correct gather set cell vertices array in the same way as local cell vertices array, // replace the padded vertices with the last vertices in the corresponding cells for (int gather_set_cell_idx = 0; gather_set_cell_idx < nCells; gather_set_cell_idx++) { int num_edges = num_edges_on_gather_set_cells[gather_set_cell_idx]; int idx_in_gather_set_vert_arr = gather_set_cell_idx * maxEdgesPerCell; int last_vert_idx = vertices_on_gather_set_cells[idx_in_gather_set_vert_arr + num_edges - 1]; for (int i = num_edges; i < maxEdgesPerCell; i++) vertices_on_gather_set_cells[idx_in_gather_set_vert_arr + i] = last_vert_idx; } // Populate connectivity data for gather set cells // Convert in-place from int (stored in the first half) to EntityHandle // Reading backward is the trick for (int cell_vert = nCells * maxEdgesPerCell - 1; cell_vert >= 0; cell_vert--) { int gather_set_vert_idx = vertices_on_gather_set_cells[cell_vert]; // Global vertex index, 1 based gather_set_vert_idx--; // 1 based -> 0 based // Connectivity array is shifted by where the gather set vertices start conn_arr_gather_set_cells[cell_vert] = gather_set_start_vertex + gather_set_vert_idx; } return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::create_padded_local_cells | ( | const std::vector< int > & | vertices_on_local_cells, |
EntityHandle | start_vertex, | ||
Range & | faces | ||
) | [private] |
Create local cells with padding (padded cells will have the same number of edges)
Definition at line 1389 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; Tag& mGlobalIdTag = _readNC->mGlobalIdTag; // Only one group of cells (each cell is represented by a polygon with maxEdgesPerCell edges) numCellGroups = 1; // Create cells for this cell group EntityHandle start_element; EntityHandle* conn_arr_local_cells = NULL; ErrorCode rval = _readNC->readMeshIface->get_element_connect(nLocalCells, maxEdgesPerCell, MBPOLYGON, 0, start_element, conn_arr_local_cells, // Might have to create gather mesh later (createGatherSet ? nLocalCells + nCells : nLocalCells)); ERRORR(rval, "Failed to create cells."); faces.insert(start_element, start_element + nLocalCells - 1); // Add local cells to the file set Range local_cells_range(start_element, start_element + nLocalCells - 1); rval = _readNC->mbImpl->add_entities(_fileSet, local_cells_range); ERRORR(rval, "Failed to add local cells to the file set."); // Get ptr to GID memory for local cells int count = 0; void* data = NULL; rval = mbImpl->tag_iterate(mGlobalIdTag, local_cells_range.begin(), local_cells_range.end(), count, data); ERRORR(rval, "Failed to iterate global id tag on local cells."); assert(count == nLocalCells); int* gid_data = (int*) data; std::copy(localGidCells.begin(), localGidCells.end(), gid_data); // Set connectivity array with proper local vertices handles // vertices_on_local_cells array was already corrected to have the last vertices padded // no need for extra checks considering for (int local_cell_idx = 0; local_cell_idx < nLocalCells; local_cell_idx++) { for (int i = 0; i < maxEdgesPerCell; i++) { EntityHandle global_vert_idx = vertices_on_local_cells[local_cell_idx * maxEdgesPerCell + i]; // Global vertex index, 1 based int local_vert_idx = localGidVerts.index(global_vert_idx); // Local vertex index, 0 based assert(local_vert_idx != -1); conn_arr_local_cells[local_cell_idx * maxEdgesPerCell + i] = start_vertex + local_vert_idx; } } return MB_SUCCESS; }
virtual std::string moab::NCHelperMPAS::get_mesh_type_name | ( | ) | [inline, private, virtual] |
Implementation of NCHelper::get_mesh_type_name()
Implements moab::NCHelper.
Definition at line 31 of file NCHelperMPAS.hpp.
{ return "MPAS"; }
ErrorCode moab::NCHelperMPAS::init_mesh_vals | ( | ) | [private, virtual] |
Implementation of NCHelper::init_mesh_vals()
Implements moab::NCHelper.
Definition at line 58 of file NCHelperMPAS.cpp.
{ std::vector<std::string>& dimNames = _readNC->dimNames; std::vector<int>& dimLens = _readNC->dimLens; std::map<std::string, ReadNC::VarData>& varInfo = _readNC->varInfo; ErrorCode rval; unsigned int idx; std::vector<std::string>::iterator vit; // Get max edges per cell reported in the MPAS file header if ((vit = std::find(dimNames.begin(), dimNames.end(), "maxEdges")) != dimNames.end()) { idx = vit - dimNames.begin(); maxEdgesPerCell = dimLens[idx]; if (maxEdgesPerCell > DEFAULT_MAX_EDGES_PER_CELL) { ERRORR(MB_FAILURE, "maxEdgesPerCell read from the MPAS file header has exceeded DEFAULT_MAX_EDGES_PER_CELL."); } } // Look for time dimension if ((vit = std::find(dimNames.begin(), dimNames.end(), "Time")) != dimNames.end()) idx = vit - dimNames.begin(); else if ((vit = std::find(dimNames.begin(), dimNames.end(), "time")) != dimNames.end()) idx = vit - dimNames.begin(); else { ERRORR(MB_FAILURE, "Couldn't find 'Time' or 'time' dimension."); } tDim = idx; nTimeSteps = dimLens[idx]; // Get number of cells if ((vit = std::find(dimNames.begin(), dimNames.end(), "nCells")) != dimNames.end()) idx = vit - dimNames.begin(); else { ERRORR(MB_FAILURE, "Couldn't find 'nCells' dimension."); } cDim = idx; nCells = dimLens[idx]; // Get number of edges if ((vit = std::find(dimNames.begin(), dimNames.end(), "nEdges")) != dimNames.end()) idx = vit - dimNames.begin(); else { ERRORR(MB_FAILURE, "Couldn't find 'nEdges' dimension."); } eDim = idx; nEdges = dimLens[idx]; // Get number of vertices vDim = -1; if ((vit = std::find(dimNames.begin(), dimNames.end(), "nVertices")) != dimNames.end()) idx = vit - dimNames.begin(); else { ERRORR(MB_FAILURE, "Couldn't find 'nVertices' dimension."); } vDim = idx; nVertices = dimLens[idx]; // Get number of vertex levels if ((vit = std::find(dimNames.begin(), dimNames.end(), "nVertLevels")) != dimNames.end()) idx = vit - dimNames.begin(); else { ERRORR(MB_FAILURE, "Couldn't find 'nVertLevels' dimension."); } levDim = idx; nLevels = dimLens[idx]; // Dimension numbers for other optional levels std::vector<unsigned int> opt_lev_dims; // Get number of vertex levels P1 if ((vit = std::find(dimNames.begin(), dimNames.end(), "nVertLevelsP1")) != dimNames.end()) { idx = vit - dimNames.begin(); opt_lev_dims.push_back(idx); } // Get number of vertex levels P2 if ((vit = std::find(dimNames.begin(), dimNames.end(), "nVertLevelsP2")) != dimNames.end()) { idx = vit - dimNames.begin(); opt_lev_dims.push_back(idx); } // Get number of soil levels if ((vit = std::find(dimNames.begin(), dimNames.end(), "nSoilLevels")) != dimNames.end()) { idx = vit - dimNames.begin(); opt_lev_dims.push_back(idx); } std::map<std::string, ReadNC::VarData>::iterator vmit; // Store time coordinate values in tVals if (nTimeSteps > 0) { // Note, two possible types for xtime variable: double(Time) or char(Time, StrLen) if ((vmit = varInfo.find("xtime")) != varInfo.end() && (*vmit).second.varDims.size() == 1) { // If xtime variable is double type, read time coordinate values to tVals rval = read_coordinate("xtime", 0, nTimeSteps - 1, tVals); ERRORR(rval, "Trouble reading 'xtime' variable."); } else { // If xtime variable does not exist, or it is string type, set dummy values to tVals for (int t = 0; t < nTimeSteps; t++) tVals.push_back((double)t); } } // For each variable, determine the entity location type and number of levels for (vmit = varInfo.begin(); vmit != varInfo.end(); ++vmit) { ReadNC::VarData& vd = (*vmit).second; vd.entLoc = ReadNC::ENTLOCSET; if (std::find(vd.varDims.begin(), vd.varDims.end(), tDim) != vd.varDims.end()) { if (std::find(vd.varDims.begin(), vd.varDims.end(), vDim) != vd.varDims.end()) vd.entLoc = ReadNC::ENTLOCVERT; else if (std::find(vd.varDims.begin(), vd.varDims.end(), eDim) != vd.varDims.end()) vd.entLoc = ReadNC::ENTLOCEDGE; else if (std::find(vd.varDims.begin(), vd.varDims.end(), cDim) != vd.varDims.end()) vd.entLoc = ReadNC::ENTLOCFACE; } vd.numLev = 1; if (std::find(vd.varDims.begin(), vd.varDims.end(), levDim) != vd.varDims.end()) vd.numLev = nLevels; else { // If nVertLevels dimension is not found, try other optional levels such as nVertLevelsP1 for (unsigned int i = 0; i < opt_lev_dims.size(); i++) { if (std::find(vd.varDims.begin(), vd.varDims.end(), opt_lev_dims[i]) != vd.varDims.end()) { vd.numLev = dimLens[opt_lev_dims[i]]; break; } } } // Hack: ignore variables with more than 3 dimensions, e.g. tracers(Time, nCells, nVertLevels, nTracers) if (vd.varDims.size() > 3) ignoredVarNames.insert(vd.varName); } // Hack: create dummy variables for dimensions (like nCells) with no corresponding coordinate variables rval = create_dummy_variables(); ERRORR(rval, "Failed to create dummy variables."); return MB_SUCCESS; }
ErrorCode moab::NCHelperMPAS::read_ucd_variable_to_nonset | ( | std::vector< ReadNC::VarData > & | vdatas, |
std::vector< int > & | tstep_nums | ||
) | [private, virtual] |
Implementation of UcdNCHelper::read_ucd_variable_to_nonset()
Implements moab::UcdNCHelper.
Definition at line 806 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; bool& noEdges = _readNC->noEdges; DebugOutput& dbgOut = _readNC->dbgOut; ErrorCode rval = read_ucd_variable_to_nonset_allocate(vdatas, tstep_nums); ERRORR(rval, "Trouble allocating read variables."); // Finally, read into that space int success; Range* pLocalGid = NULL; for (unsigned int i = 0; i < vdatas.size(); i++) { // Skip edge variables, if specified by the read options if (noEdges && vdatas[i].entLoc == ReadNC::ENTLOCEDGE) continue; switch (vdatas[i].entLoc) { case ReadNC::ENTLOCVERT: pLocalGid = &localGidVerts; break; case ReadNC::ENTLOCFACE: pLocalGid = &localGidCells; break; case ReadNC::ENTLOCEDGE: pLocalGid = &localGidEdges; break; default: ERRORR(MB_FAILURE, "Unrecognized entity location type."); break; } std::size_t sz = vdatas[i].sz; for (unsigned int t = 0; t < tstep_nums.size(); t++) { // Set readStart for each timestep along time dimension vdatas[i].readStarts[0] = tstep_nums[t]; switch (vdatas[i].varDataType) { case NC_BYTE: case NC_CHAR: { ERRORR(MB_FAILURE, "not implemented"); break; } case NC_DOUBLE: { std::vector<double> tmpdoubledata(sz); // In the case of ucd mesh, and on multiple proc, // we need to read as many times as subranges we have in the // localGid range; // basically, we have to give a different point // for data to start, for every subrange :( size_t indexInDoubleArray = 0; size_t ic = 0; for (Range::pair_iterator pair_iter = pLocalGid->pair_begin(); pair_iter != pLocalGid->pair_end(); pair_iter++, ic++) { EntityHandle starth = pair_iter->first; EntityHandle endh = pair_iter->second; // Inclusive vdatas[i].readStarts[1] = (NCDF_SIZE) (starth - 1); vdatas[i].readCounts[1] = (NCDF_SIZE) (endh - starth + 1); success = NCFUNCAG(_vara_double)(_fileId, vdatas[i].varId, &(vdatas[i].readStarts[0]), &(vdatas[i].readCounts[0]), &(tmpdoubledata[indexInDoubleArray])); ERRORS(success, "Failed to read double data in loop"); // We need to increment the index in double array for the // next subrange indexInDoubleArray += (endh - starth + 1) * 1 * vdatas[i].numLev; } assert(ic == pLocalGid->psize()); if (vdatas[i].entLoc == ReadNC::ENTLOCFACE && numCellGroups > 1) { // For a cell variable that is NOT on one contiguous chunk of faces, allocate tag space for // each cell group, and utilize cellHandleToGlobalID map to read tag data Range::iterator iter = facesOwned.begin(); while (iter != facesOwned.end()) { int count; void* ptr; rval = mbImpl->tag_iterate(vdatas[i].varTags[t], iter, facesOwned.end(), count, ptr); ERRORR(rval, "Failed to iterate tag on owned faces."); for (int j = 0; j < count; j++) { int global_cell_idx = cellHandleToGlobalID[*(iter + j)]; // Global cell index, 1 based int local_cell_idx = localGidCells.index(global_cell_idx); // Local cell index, 0 based assert(local_cell_idx != -1); for (int level = 0; level < vdatas[i].numLev; level++) ((double*) ptr)[j * vdatas[i].numLev + level] = tmpdoubledata[local_cell_idx * vdatas[i].numLev + level]; } iter += count; } } else { void* data = vdatas[i].varDatas[t]; for (std::size_t idx = 0; idx != tmpdoubledata.size(); idx++) ((double*) data)[idx] = tmpdoubledata[idx]; } break; } case NC_FLOAT: { ERRORR(MB_FAILURE, "not implemented"); break; } case NC_INT: { ERRORR(MB_FAILURE, "not implemented"); break; } case NC_SHORT: { ERRORR(MB_FAILURE, "not implemented"); break; } default: success = 1; } if (success) ERRORR(MB_FAILURE, "Trouble reading variable."); } } for (unsigned int i = 0; i < vdatas.size(); i++) { if (noEdges && vdatas[i].entLoc == ReadNC::ENTLOCEDGE) continue; for (unsigned int t = 0; t < tstep_nums.size(); t++) { dbgOut.tprintf(2, "Converting variable %s, time step %d\n", vdatas[i].varName.c_str(), tstep_nums[t]); ErrorCode tmp_rval = convert_variable(vdatas[i], t); if (MB_SUCCESS != tmp_rval) rval = tmp_rval; } } // Debug output, if requested if (1 == dbgOut.get_verbosity()) { dbgOut.printf(1, "Read variables: %s", vdatas.begin()->varName.c_str()); for (unsigned int i = 1; i < vdatas.size(); i++) dbgOut.printf(1, ", %s ", vdatas[i].varName.c_str()); dbgOut.tprintf(1, "\n"); } return rval; }
ErrorCode moab::NCHelperMPAS::read_ucd_variable_to_nonset_allocate | ( | std::vector< ReadNC::VarData > & | vdatas, |
std::vector< int > & | tstep_nums | ||
) | [private, virtual] |
Implementation of UcdNCHelper::read_ucd_variable_to_nonset_allocate()
Implements moab::UcdNCHelper.
Definition at line 512 of file NCHelperMPAS.cpp.
{ Interface*& mbImpl = _readNC->mbImpl; std::vector<int>& dimLens = _readNC->dimLens; bool& noEdges = _readNC->noEdges; DebugOutput& dbgOut = _readNC->dbgOut; ErrorCode rval = MB_SUCCESS; Range* range = NULL; // Get vertices in set Range verts; rval = mbImpl->get_entities_by_dimension(_fileSet, 0, verts); ERRORR(rval, "Trouble getting vertices in set."); assert("Should only have a single vertex subrange, since they were read in one shot" && verts.psize() == 1); // Get edges in set Range edges; rval = mbImpl->get_entities_by_dimension(_fileSet, 1, edges); ERRORR(rval, "Trouble getting edges in set."); // Get faces in set Range faces; rval = mbImpl->get_entities_by_dimension(_fileSet, 2, faces); ERRORR(rval, "Trouble getting faces in set."); // Note, for MPAS faces.psize() can be more than 1 #ifdef USE_MPI bool& isParallel = _readNC->isParallel; if (isParallel) { ParallelComm*& myPcomm = _readNC->myPcomm; rval = myPcomm->filter_pstatus(faces, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1, &facesOwned); ERRORR(rval, "Trouble getting owned faces in set."); } else facesOwned = faces; // not running in parallel, but still with MPI #else facesOwned = faces; #endif for (unsigned int i = 0; i < vdatas.size(); i++) { // Skip edge variables, if specified by the read options if (noEdges && vdatas[i].entLoc == ReadNC::ENTLOCEDGE) continue; // Support non-set variables with 3 dimensions like (Time, nCells, nVertLevels), or // 2 dimensions like (Time, nCells) assert(3 == vdatas[i].varDims.size() || 2 == vdatas[i].varDims.size()); // For a non-set variable, time should be the first dimension assert(tDim == vdatas[i].varDims[0]); // Set up readStarts and readCounts vdatas[i].readStarts.resize(3); vdatas[i].readCounts.resize(3); // First: Time vdatas[i].readStarts[0] = 0; // This value is timestep dependent, will be set later vdatas[i].readCounts[0] = 1; // Next: nVertices / nCells / nEdges switch (vdatas[i].entLoc) { case ReadNC::ENTLOCVERT: // Vertices // Start from the first localGidVerts // Actually, this will be reset later on in a loop vdatas[i].readStarts[1] = localGidVerts[0] - 1; vdatas[i].readCounts[1] = nLocalVertices; range = &verts; break; case ReadNC::ENTLOCFACE: // Faces // Start from the first localGidCells // Actually, this will be reset later on in a loop vdatas[i].readStarts[1] = localGidCells[0] - 1; vdatas[i].readCounts[1] = nLocalCells; range = &facesOwned; break; case ReadNC::ENTLOCEDGE: // Edges // Start from the first localGidEdges // Actually, this will be reset later on in a loop vdatas[i].readStarts[1] = localGidEdges[0] - 1; vdatas[i].readCounts[1] = nLocalEdges; range = &edges; break; default: ERRORR(MB_FAILURE, "Unexpected entity location type for MPAS non-set variable."); break; } // Finally: nVertLevels or other optional levels, it is possible that there // is no level dimension for this non-set variable vdatas[i].readStarts[2] = 0; vdatas[i].readCounts[2] = vdatas[i].numLev; // Get variable size vdatas[i].sz = 1; for (std::size_t idx = 0; idx != 3; idx++) vdatas[i].sz *= vdatas[i].readCounts[idx]; for (unsigned int t = 0; t < tstep_nums.size(); t++) { dbgOut.tprintf(2, "Reading variable %s, time step %d\n", vdatas[i].varName.c_str(), tstep_nums[t]); if (tstep_nums[t] >= dimLens[tDim]) { ERRORR(MB_INDEX_OUT_OF_RANGE, "Wrong value for a timestep number."); } // Get the tag to read into if (!vdatas[i].varTags[t]) { rval = get_tag_to_nonset(vdatas[i], tstep_nums[t], vdatas[i].varTags[t], vdatas[i].numLev); ERRORR(rval, "Trouble getting tag."); } // Get ptr to tag space if (vdatas[i].entLoc == ReadNC::ENTLOCFACE && numCellGroups > 1) { // For a cell variable that is NOT on one contiguous chunk of faces, defer its tag space allocation vdatas[i].varDatas[t] = NULL; } else { assert(1 == range->psize()); void* data; int count; rval = mbImpl->tag_iterate(vdatas[i].varTags[t], range->begin(), range->end(), count, data); ERRORR(rval, "Failed to iterate tag."); assert((unsigned)count == range->size()); vdatas[i].varDatas[t] = data; } } } return rval; }
ErrorCode moab::NCHelperMPAS::redistribute_local_cells | ( | int | start_cell_index | ) | [private] |
Redistribute local cells after trivial partition (e.g. Zoltan partition, if applicable)
Definition at line 953 of file NCHelperMPAS.cpp.
{ // If possible, apply Zoltan partition if (_readNC->partMethod == ScdParData::RCBZOLTAN) { #if defined(USE_MPI) && defined(HAVE_ZOLTAN) // Read x coordinates of cell centers int xCellVarId; int success = NCFUNC(inq_varid)(_fileId, "xCell", &xCellVarId); ERRORS(success, "Failed to get variable id of xCell."); std::vector<double> xCell(nLocalCells); NCDF_SIZE read_start = static_cast<NCDF_SIZE>(start_cell_idx - 1); NCDF_SIZE read_count = static_cast<NCDF_SIZE>(nLocalCells); success = NCFUNCAG(_vara_double)(_fileId, xCellVarId, &read_start, &read_count, &xCell[0]); ERRORS(success, "Failed to read xCell data."); // Read y coordinates of cell centers int yCellVarId; success = NCFUNC(inq_varid)(_fileId, "yCell", &yCellVarId); ERRORS(success, "Failed to get variable id of yCell."); std::vector<double> yCell(nLocalCells); success = NCFUNCAG(_vara_double)(_fileId, yCellVarId, &read_start, &read_count, &yCell[0]); ERRORS(success, "Failed to read yCell data."); // Read z coordinates of cell centers int zCellVarId; success = NCFUNC(inq_varid)(_fileId, "zCell", &zCellVarId); ERRORS(success, "Failed to get variable id of zCell."); std::vector<double> zCell(nLocalCells); success = NCFUNCAG(_vara_double)(_fileId, zCellVarId, &read_start, &read_count, &zCell[0]); ERRORS(success, "Failed to read zCell data."); // Zoltan partition using RCB; maybe more studies would be good, as to which partition // is better Interface*& mbImpl = _readNC->mbImpl; DebugOutput& dbgOut = _readNC->dbgOut; MBZoltan* mbZTool = new MBZoltan(mbImpl, false, 0, NULL); ErrorCode rval = mbZTool->repartition(xCell, yCell, zCell, start_cell_idx, "RCB", localGidCells); delete mbZTool; ERRORR(rval, "Error in Zoltan partitioning."); dbgOut.tprintf(1, "After Zoltan partitioning, localGidCells.psize() = %d\n", (int)localGidCells.psize()); dbgOut.tprintf(1, " localGidCells.size() = %d\n", (int)localGidCells.size()); // This is important: local cells are now redistributed, so nLocalCells might be different! nLocalCells = localGidCells.size(); return MB_SUCCESS; #endif } // By default, apply trivial partition localGidCells.insert(start_cell_idx, start_cell_idx + nLocalCells - 1); return MB_SUCCESS; }
std::map<EntityHandle, int> moab::NCHelperMPAS::cellHandleToGlobalID [private] |
Definition at line 80 of file NCHelperMPAS.hpp.
bool moab::NCHelperMPAS::createGatherSet [private] |
Definition at line 79 of file NCHelperMPAS.hpp.
Range moab::NCHelperMPAS::facesOwned [private] |
Definition at line 81 of file NCHelperMPAS.hpp.
int moab::NCHelperMPAS::maxEdgesPerCell [private] |
Definition at line 77 of file NCHelperMPAS.hpp.
int moab::NCHelperMPAS::numCellGroups [private] |
Definition at line 78 of file NCHelperMPAS.hpp.