MOAB: Mesh Oriented datABase
(version 5.4.1)
|
#include <iostream>
#include <exception>
#include <cmath>
#include <cassert>
#include <vector>
#include <string>
#include <fstream>
#include <iomanip>
#include "moab/ProgOptions.hpp"
#include "moab/Core.hpp"
#include "OfflineMap.h"
#include "netcdfcpp.h"
#include "NetCDFUtilities.h"
#include "DataArray2D.h"
Go to the source code of this file.
Defines | |
#define | DTYPE(a) |
Functions | |
template<typename T > | |
ErrorCode | get_vartag_data (moab::Interface *mbCore, Tag tag, moab::Range &sets, int &out_data_size, std::vector< T > &data) |
void | ReadFileMetaData (std::string &metaFilename, std::map< std::string, std::string > &metadataVals) |
int | main (int argc, char *argv[]) |
#define DTYPE | ( | a | ) |
{ \ ( ( ( a ) == 0 ) ? "FV" : ( ( ( a ) == 1 ) ? "cGLL" : "dGLL" ) ) \ }
Referenced by main().
ErrorCode get_vartag_data | ( | moab::Interface * | mbCore, |
Tag | tag, | ||
moab::Range & | sets, | ||
int & | out_data_size, | ||
std::vector< T > & | data | ||
) |
Definition at line 33 of file h5mtoscrip.cpp.
References ErrorCode, is, MB_CHK_SET_ERR, MB_SUCCESS, moab::Range::size(), T, and moab::Interface::tag_get_by_ptr().
Referenced by main().
{ int* tag_sizes = new int[sets.size()]; const void** tag_data = (const void**)new void*[sets.size()]; ErrorCode rval = mbCore->tag_get_by_ptr( tag, sets, tag_data, tag_sizes );MB_CHK_SET_ERR( rval, "Getting matrix rows failed" ); out_data_size = 0; for( unsigned is = 0; is < sets.size(); ++is ) out_data_size += tag_sizes[is]; data.resize( out_data_size ); int ioffset = 0; for( unsigned index = 0; index < sets.size(); index++ ) { T* m_vals = (T*)tag_data[index]; for( int k = 0; k < tag_sizes[index]; k++ ) { data[ioffset++] = m_vals[k]; } } return moab::MB_SUCCESS; }
int main | ( | int | argc, |
char * | argv[] | ||
) |
Definition at line 81 of file h5mtoscrip.cpp.
References ProgOptions::addOpt(), DTYPE, ErrorCode, moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type_and_tag(), get_vartag_data(), moab::Interface::globalId_tag(), moab::Interface::load_mesh(), MB_CHK_ERR, MB_CHK_SET_ERR, MB_TAG_SPARSE, MB_TYPE_INTEGER, MBENTITYSET, ProgOptions::parseCommandLine(), ProgOptions::printHelp(), ReadFileMetaData(), moab::Range::size(), moab::Interface::tag_get_data(), moab::Interface::tag_get_handle(), and moab::Interface::UNION.
{ moab::ErrorCode rval; int dimension = 2; NcError error2( NcError::verbose_nonfatal ); std::stringstream sstr; ProgOptions opts; std::string h5mfilename, scripfile; bool noMap = false; bool writeXYCoords = false; #ifdef MOAB_HAVE_MPI MPI_Init( &argc, &argv ); #endif opts.addOpt< std::string >( "weights,w", "h5m remapping weights filename", &h5mfilename ); opts.addOpt< std::string >( "scrip,s", "Output SCRIP map filename", &scripfile ); opts.addOpt< int >( "dim,d", "Dimension of entities to use for partitioning", &dimension ); opts.addOpt< void >( "mesh,m", "Only convert the mesh and exclude the remap weight details", &noMap ); opts.addOpt< void >( "coords,c", "Write the center and vertex coordinates in lat/lon format", &writeXYCoords ); opts.parseCommandLine( argc, argv ); if( h5mfilename.empty() || scripfile.empty() ) { opts.printHelp(); exit( 1 ); } moab::Interface* mbCore = new( std::nothrow ) moab::Core; if( NULL == mbCore ) { return 1; } // Set the read options for parallel file loading const std::string partition_set_name = "PARALLEL_PARTITION"; const std::string global_id_name = "GLOBAL_ID"; // Load file rval = mbCore->load_mesh( h5mfilename.c_str() );MB_CHK_ERR( rval ); try { // Temporarily change rval reporting NcError error_temp( NcError::verbose_fatal ); // Open an output file NcFile ncMap( scripfile.c_str(), NcFile::Replace, NULL, 0, NcFile::Offset64Bits ); if( !ncMap.is_valid() ) { _EXCEPTION1( "Unable to open output map file \"%s\"", scripfile.c_str() ); } { // NetCDF-SCRIP Global Attributes std::map< std::string, std::string > mapAttributes; size_t lastindex = h5mfilename.find_last_of( "." ); std::stringstream sstr; sstr << h5mfilename.substr( 0, lastindex ) << ".meta"; std::string metaFilename = sstr.str(); ReadFileMetaData( metaFilename, mapAttributes ); mapAttributes["Command"] = "Converted with MOAB:h5mtoscrip with --w=" + h5mfilename + " and --s=" + scripfile; // Add global attributes std::map< std::string, std::string >::const_iterator iterAttributes = mapAttributes.begin(); for( ; iterAttributes != mapAttributes.end(); iterAttributes++ ) { std::cout << iterAttributes->first << " -- " << iterAttributes->second << std::endl; ncMap.add_att( iterAttributes->first.c_str(), iterAttributes->second.c_str() ); } std::cout << "\n"; } Tag globalIDTag, materialSetTag; globalIDTag = mbCore->globalId_tag(); // materialSetTag = mbCore->material_tag(); rval = mbCore->tag_get_handle( "MATERIAL_SET", 1, MB_TYPE_INTEGER, materialSetTag, MB_TAG_SPARSE );MB_CHK_ERR( rval ); // Get sets entities, by type moab::Range meshsets; rval = mbCore->get_entities_by_type_and_tag( 0, MBENTITYSET, &globalIDTag, NULL, 1, meshsets, moab::Interface::UNION, true );MB_CHK_ERR( rval ); moab::EntityHandle rootset = 0; /////////////////////////////////////////////////////////////////////////// // The metadata in H5M file contains the following data: // // 1. n_a: Total source entities: (number of elements in source mesh) // 2. n_b: Total target entities: (number of elements in target mesh) // 3. nv_a: Max edge size of elements in source mesh // 4. nv_b: Max edge size of elements in target mesh // 5. maxrows: Number of rows in remap weight matrix // 6. maxcols: Number of cols in remap weight matrix // 7. nnz: Number of total nnz in sparse remap weight matrix // 8. np_a: The order of the field description on the source mesh: >= 1 // 9. np_b: The order of the field description on the target mesh: >= 1 // 10. method_a: The type of discretization for field on source mesh: [0 = FV, 1 = cGLL, 2 // = dGLL] // 11. method_b: The type of discretization for field on target mesh: [0 = FV, 1 = cGLL, 2 // = dGLL] // 12. conserved: Flag to specify whether the remap operator has conservation constraints: // [0, 1] // 13. monotonicity: Flags to specify whether the remap operator has monotonicity // constraints: [0, 1, 2] // /////////////////////////////////////////////////////////////////////////// Tag smatMetadataTag; int smat_metadata_glb[13]; rval = mbCore->tag_get_handle( "SMAT_DATA", 13, MB_TYPE_INTEGER, smatMetadataTag, MB_TAG_SPARSE );MB_CHK_ERR( rval ); rval = mbCore->tag_get_data( smatMetadataTag, &rootset, 1, smat_metadata_glb );MB_CHK_ERR( rval ); // std::cout << "Number of mesh sets is " << meshsets.size() << std::endl; #define DTYPE( a ) \ { \ ( ( ( a ) == 0 ) ? "FV" : ( ( ( a ) == 1 ) ? "cGLL" : "dGLL" ) ) \ } // Map dimensions int nA = smat_metadata_glb[0]; int nB = smat_metadata_glb[1]; int nVA = smat_metadata_glb[2]; int nVB = smat_metadata_glb[3]; int nDofB = smat_metadata_glb[4]; int nDofA = smat_metadata_glb[5]; int NNZ = smat_metadata_glb[6]; int nOrdA = smat_metadata_glb[7]; int nOrdB = smat_metadata_glb[8]; int nBasA = smat_metadata_glb[9]; std::string methodA = DTYPE( nBasA ); int nBasB = smat_metadata_glb[10]; std::string methodB = DTYPE( nBasB ); int bConserved = smat_metadata_glb[11]; int bMonotonicity = smat_metadata_glb[12]; EntityHandle source_mesh = 0, target_mesh = 0, overlap_mesh = 0; for( unsigned im = 0; im < meshsets.size(); ++im ) { moab::Range elems; rval = mbCore->get_entities_by_dimension( meshsets[im], 2, elems );MB_CHK_ERR( rval ); if( elems.size() - nA == 0 && source_mesh == 0 ) source_mesh = meshsets[im]; else if( elems.size() - nB == 0 && target_mesh == 0 ) target_mesh = meshsets[im]; else if( overlap_mesh == 0 ) overlap_mesh = meshsets[im]; else continue; } Tag srcIDTag, srcAreaTag, tgtIDTag, tgtAreaTag; rval = mbCore->tag_get_handle( "SourceGIDS", srcIDTag );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "SourceAreas", srcAreaTag );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "TargetGIDS", tgtIDTag );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "TargetAreas", tgtAreaTag );MB_CHK_ERR( rval ); Tag smatRowdataTag, smatColdataTag, smatValsdataTag; rval = mbCore->tag_get_handle( "SMAT_ROWS", smatRowdataTag );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "SMAT_COLS", smatColdataTag );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "SMAT_VALS", smatValsdataTag );MB_CHK_ERR( rval ); Tag srcCenterLon, srcCenterLat, tgtCenterLon, tgtCenterLat; rval = mbCore->tag_get_handle( "SourceCoordCenterLon", srcCenterLon );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "SourceCoordCenterLat", srcCenterLat );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "TargetCoordCenterLon", tgtCenterLon );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "TargetCoordCenterLat", tgtCenterLat );MB_CHK_ERR( rval ); Tag srcVertexLon, srcVertexLat, tgtVertexLon, tgtVertexLat; rval = mbCore->tag_get_handle( "SourceCoordVertexLon", srcVertexLon );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "SourceCoordVertexLat", srcVertexLat );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "TargetCoordVertexLon", tgtVertexLon );MB_CHK_ERR( rval ); rval = mbCore->tag_get_handle( "TargetCoordVertexLat", tgtVertexLat );MB_CHK_ERR( rval ); // Get sets entities, by type moab::Range sets; // rval = mbCore->get_entities_by_type(0, MBENTITYSET, sets);MB_CHK_ERR(rval); rval = mbCore->get_entities_by_type_and_tag( 0, MBENTITYSET, &smatRowdataTag, NULL, 1, sets, moab::Interface::UNION, true );MB_CHK_ERR( rval ); std::vector< int > src_gids, tgt_gids; std::vector< double > src_areas, tgt_areas; int srcID_size, tgtID_size, srcArea_size, tgtArea_size; rval = get_vartag_data( mbCore, srcIDTag, sets, srcID_size, src_gids );MB_CHK_SET_ERR( rval, "Getting source mesh IDs failed" ); rval = get_vartag_data( mbCore, tgtIDTag, sets, tgtID_size, tgt_gids );MB_CHK_SET_ERR( rval, "Getting target mesh IDs failed" ); rval = get_vartag_data( mbCore, srcAreaTag, sets, srcArea_size, src_areas );MB_CHK_SET_ERR( rval, "Getting source mesh areas failed" ); rval = get_vartag_data( mbCore, tgtAreaTag, sets, tgtArea_size, tgt_areas );MB_CHK_SET_ERR( rval, "Getting target mesh areas failed" ); assert( srcArea_size == srcID_size ); assert( tgtArea_size == tgtID_size ); std::vector< double > src_glob_areas( nDofA, 0.0 ), tgt_glob_areas( nDofB, 0.0 ); for( int i = 0; i < srcArea_size; ++i ) { // printf("%d/%d: %d = Found ID %d and area %5.6e\n", i, srcArea_size, nDofA, // src_gids[i], src_areas[i]); assert( i < srcID_size ); assert( src_gids[i] < nDofA ); if( src_areas[i] > src_glob_areas[src_gids[i]] ) src_glob_areas[src_gids[i]] = src_areas[i]; } for( int i = 0; i < tgtArea_size; ++i ) { // printf("%d/%d: %d = Found ID %d and area %5.6e\n", i, tgtArea_size, nDofB, // tgt_gids[i], tgt_areas[i]); assert( i < tgtID_size ); assert( tgt_gids[i] < nDofB ); if( tgt_areas[i] > tgt_glob_areas[tgt_gids[i]] ) tgt_glob_areas[tgt_gids[i]] = tgt_areas[i]; } // Write output dimensions entries int nSrcGridDims = 1; int nDstGridDims = 1; NcDim* dimSrcGridRank = ncMap.add_dim( "src_grid_rank", nSrcGridDims ); NcDim* dimDstGridRank = ncMap.add_dim( "dst_grid_rank", nDstGridDims ); NcVar* varSrcGridDims = ncMap.add_var( "src_grid_dims", ncInt, dimSrcGridRank ); NcVar* varDstGridDims = ncMap.add_var( "dst_grid_dims", ncInt, dimDstGridRank ); if( nA == nDofA ) { varSrcGridDims->put( &nA, 1 ); varSrcGridDims->add_att( "name0", "num_elem" ); } else { varSrcGridDims->put( &nDofA, 1 ); varSrcGridDims->add_att( "name1", "num_dof" ); } if( nB == nDofB ) { varDstGridDims->put( &nB, 1 ); varDstGridDims->add_att( "name0", "num_elem" ); } else { varDstGridDims->put( &nDofB, 1 ); varDstGridDims->add_att( "name1", "num_dof" ); } // Source and Target mesh resolutions NcDim* dimNA = ncMap.add_dim( "n_a", nDofA ); NcDim* dimNB = ncMap.add_dim( "n_b", nDofB ); // Source and Target verticecs per elements const int nva = ( nA == nDofA ? nVA : 1 ); const int nvb = ( nB == nDofB ? nVB : 1 ); NcDim* dimNVA = ncMap.add_dim( "nv_a", nva ); NcDim* dimNVB = ncMap.add_dim( "nv_b", nvb ); // Source and Target verticecs per elements // NcDim * dimNEA = ncMap.add_dim("ne_a", nA); // NcDim * dimNEB = ncMap.add_dim("ne_b", nB); if( writeXYCoords ) { // Write coordinates NcVar* varYCA = ncMap.add_var( "yc_a", ncDouble, dimNA /*dimNA*/ ); NcVar* varYCB = ncMap.add_var( "yc_b", ncDouble, dimNB /*dimNB*/ ); NcVar* varXCA = ncMap.add_var( "xc_a", ncDouble, dimNA /*dimNA*/ ); NcVar* varXCB = ncMap.add_var( "xc_b", ncDouble, dimNB /*dimNB*/ ); NcVar* varYVA = ncMap.add_var( "yv_a", ncDouble, dimNA /*dimNA*/, dimNVA ); NcVar* varYVB = ncMap.add_var( "yv_b", ncDouble, dimNB /*dimNB*/, dimNVB ); NcVar* varXVA = ncMap.add_var( "xv_a", ncDouble, dimNA /*dimNA*/, dimNVA ); NcVar* varXVB = ncMap.add_var( "xv_b", ncDouble, dimNB /*dimNB*/, dimNVB ); varYCA->add_att( "units", "degrees" ); varYCB->add_att( "units", "degrees" ); varXCA->add_att( "units", "degrees" ); varXCB->add_att( "units", "degrees" ); varYVA->add_att( "units", "degrees" ); varYVB->add_att( "units", "degrees" ); varXVA->add_att( "units", "degrees" ); varXVB->add_att( "units", "degrees" ); std::vector< double > src_centerlat, src_centerlon; int srccenter_size; rval = get_vartag_data( mbCore, srcCenterLat, sets, srccenter_size, src_centerlat );MB_CHK_SET_ERR( rval, "Getting source mesh areas failed" ); rval = get_vartag_data( mbCore, srcCenterLon, sets, srccenter_size, src_centerlon );MB_CHK_SET_ERR( rval, "Getting target mesh areas failed" ); std::vector< double > src_glob_centerlat( nDofA, 0.0 ), src_glob_centerlon( nDofA, 0.0 ); for( int i = 0; i < srccenter_size; ++i ) { assert( i < srcID_size ); assert( src_gids[i] < nDofA ); src_glob_centerlat[src_gids[i]] = src_centerlat[i]; src_glob_centerlon[src_gids[i]] = src_centerlon[i]; } std::vector< double > tgt_centerlat, tgt_centerlon; int tgtcenter_size; rval = get_vartag_data( mbCore, tgtCenterLat, sets, tgtcenter_size, tgt_centerlat );MB_CHK_SET_ERR( rval, "Getting source mesh areas failed" ); rval = get_vartag_data( mbCore, tgtCenterLon, sets, tgtcenter_size, tgt_centerlon );MB_CHK_SET_ERR( rval, "Getting target mesh areas failed" ); std::vector< double > tgt_glob_centerlat( nDofB, 0.0 ), tgt_glob_centerlon( nDofB, 0.0 ); for( int i = 0; i < tgtcenter_size; ++i ) { assert( i < tgtID_size ); assert( tgt_gids[i] < nDofB ); tgt_glob_centerlat[tgt_gids[i]] = tgt_centerlat[i]; tgt_glob_centerlon[tgt_gids[i]] = tgt_centerlon[i]; } varYCA->put( &( src_glob_centerlat[0] ), nDofA ); varYCB->put( &( tgt_glob_centerlat[0] ), nDofB ); varXCA->put( &( src_glob_centerlon[0] ), nDofA ); varXCB->put( &( tgt_glob_centerlon[0] ), nDofB ); src_centerlat.clear(); src_centerlon.clear(); tgt_centerlat.clear(); tgt_centerlon.clear(); DataArray2D< double > src_glob_vertexlat( nDofA, nva ), src_glob_vertexlon( nDofA, nva ); if( nva > 1 ) { std::vector< double > src_vertexlat, src_vertexlon; int srcvertex_size; rval = get_vartag_data( mbCore, srcVertexLat, sets, srcvertex_size, src_vertexlat );MB_CHK_SET_ERR( rval, "Getting source mesh areas failed" ); rval = get_vartag_data( mbCore, srcVertexLon, sets, srcvertex_size, src_vertexlon );MB_CHK_SET_ERR( rval, "Getting target mesh areas failed" ); int offset = 0; for( unsigned vIndex = 0; vIndex < src_gids.size(); ++vIndex ) { for( int vNV = 0; vNV < nva; ++vNV ) { assert( offset < srcvertex_size ); src_glob_vertexlat[src_gids[vIndex]][vNV] = src_vertexlat[offset]; src_glob_vertexlon[src_gids[vIndex]][vNV] = src_vertexlon[offset]; offset++; } } } DataArray2D< double > tgt_glob_vertexlat( nDofB, nvb ), tgt_glob_vertexlon( nDofB, nvb ); if( nvb > 1 ) { std::vector< double > tgt_vertexlat, tgt_vertexlon; int tgtvertex_size; rval = get_vartag_data( mbCore, tgtVertexLat, sets, tgtvertex_size, tgt_vertexlat );MB_CHK_SET_ERR( rval, "Getting source mesh areas failed" ); rval = get_vartag_data( mbCore, tgtVertexLon, sets, tgtvertex_size, tgt_vertexlon );MB_CHK_SET_ERR( rval, "Getting target mesh areas failed" ); int offset = 0; for( unsigned vIndex = 0; vIndex < tgt_gids.size(); ++vIndex ) { for( int vNV = 0; vNV < nvb; ++vNV ) { assert( offset < tgtvertex_size ); tgt_glob_vertexlat[tgt_gids[vIndex]][vNV] = tgt_vertexlat[offset]; tgt_glob_vertexlon[tgt_gids[vIndex]][vNV] = tgt_vertexlon[offset]; offset++; } } } varYVA->put( &( src_glob_vertexlat[0][0] ), nDofA, nva ); varYVB->put( &( tgt_glob_vertexlat[0][0] ), nDofB, nvb ); varXVA->put( &( src_glob_vertexlon[0][0] ), nDofA, nva ); varXVB->put( &( tgt_glob_vertexlon[0][0] ), nDofB, nvb ); } // Write areas NcVar* varAreaA = ncMap.add_var( "area_a", ncDouble, dimNA ); varAreaA->put( &( src_glob_areas[0] ), nDofA ); // varAreaA->add_att("units", "steradians"); NcVar* varAreaB = ncMap.add_var( "area_b", ncDouble, dimNB ); varAreaB->put( &( tgt_glob_areas[0] ), nDofB ); // varAreaB->add_att("units", "steradians"); std::vector< int > mat_rows, mat_cols; std::vector< double > mat_vals; int row_sizes, col_sizes, val_sizes; rval = get_vartag_data( mbCore, smatRowdataTag, sets, row_sizes, mat_rows );MB_CHK_SET_ERR( rval, "Getting matrix row data failed" ); assert( row_sizes == NNZ ); rval = get_vartag_data( mbCore, smatColdataTag, sets, col_sizes, mat_cols );MB_CHK_SET_ERR( rval, "Getting matrix col data failed" ); assert( col_sizes == NNZ ); rval = get_vartag_data( mbCore, smatValsdataTag, sets, val_sizes, mat_vals );MB_CHK_SET_ERR( rval, "Getting matrix values failed" ); assert( val_sizes == NNZ ); // Let us form the matrix in-memory and consolidate shared DoF rows from shared-process // contributions SparseMatrix< double > mapMatrix; for( int innz = 0; innz < NNZ; ++innz ) { #ifdef VERBOSE if( fabs( mapMatrix( mat_rows[innz], mat_cols[innz] ) ) > 1e-12 ) { printf( "Adding to existing loc: (%d, %d) = %12.8f\n", mat_rows[innz], mat_cols[innz], mapMatrix( mat_rows[innz], mat_cols[innz] ) ); } #endif mapMatrix( mat_rows[innz], mat_cols[innz] ) += mat_vals[innz]; } // Write SparseMatrix entries DataArray1D< int > vecRow; DataArray1D< int > vecCol; DataArray1D< double > vecS; mapMatrix.GetEntries( vecRow, vecCol, vecS ); int nS = vecS.GetRows(); // Print more information about what we are converting: // Source elements/vertices/type (Discretization ?) // Target elements/vertices/type (Discretization ?) // Overlap elements/types // Rmeapping weights matrix: rows/cols/NNZ // Output the number of sets printf( "Primary sets: %15zu\n", sets.size() ); printf( "Original NNZ: %18d\n", NNZ ); printf( "Consolidated Total NNZ: %8d\n", nS ); printf( "Conservative weights ? %6d\n", ( bConserved > 0 ) ); printf( "Monotone weights ? %10d\n", ( bMonotonicity > 0 ) ); printf( "\n--------------------------------------------------------------\n" ); printf( "%20s %21s %15s\n", "Description", "Source", "Target" ); printf( "--------------------------------------------------------------\n" ); printf( "%25s %15d %15d\n", "Number of elements:", nA, nB ); printf( "%25s %15d %15d\n", "Number of DoFs:", nDofA, nDofB ); printf( "%25s %15d %15d\n", "Maximum vertex/element:", nVA, nVB ); printf( "%25s %15s %15s\n", "Discretization type:", methodA.c_str(), methodB.c_str() ); printf( "%25s %15d %15d\n", "Discretization order:", nOrdA, nOrdB ); // Calculate and write fractional coverage arrays { DataArray1D< double > dFracA( nDofA ); DataArray1D< double > dFracB( nDofB ); for( int i = 0; i < nS; i++ ) { // std::cout << i << " - mat_vals = " << mat_vals[i] << " dFracA = " << mat_vals[i] // / src_glob_areas[vecCol[i]] * tgt_glob_areas[vecRow[i]] << std::endl; dFracA[vecCol[i]] += vecS[i] / src_glob_areas[vecCol[i]] * tgt_glob_areas[vecRow[i]]; dFracB[vecRow[i]] += vecS[i]; } NcVar* varFracA = ncMap.add_var( "frac_a", ncDouble, dimNA ); varFracA->put( &( dFracA[0] ), nDofA ); varFracA->add_att( "name", "fraction of target coverage of source dof" ); varFracA->add_att( "units", "unitless" ); NcVar* varFracB = ncMap.add_var( "frac_b", ncDouble, dimNB ); varFracB->put( &( dFracB[0] ), nDofB ); varFracB->add_att( "name", "fraction of source coverage of target dof" ); varFracB->add_att( "units", "unitless" ); } // Write out data NcDim* dimNS = ncMap.add_dim( "n_s", nS ); NcVar* varRow = ncMap.add_var( "row", ncInt, dimNS ); varRow->add_att( "name", "sparse matrix target dof index" ); varRow->add_att( "first_index", "1" ); NcVar* varCol = ncMap.add_var( "col", ncInt, dimNS ); varCol->add_att( "name", "sparse matrix source dof index" ); varCol->add_att( "first_index", "1" ); NcVar* varS = ncMap.add_var( "S", ncDouble, dimNS ); varS->add_att( "name", "sparse matrix coefficient" ); // Increment vecRow and vecCol: make it 1-based for( int i = 0; i < nS; i++ ) { vecRow[i]++; vecCol[i]++; } varRow->set_cur( (long)0 ); varRow->put( &( vecRow[0] ), nS ); varCol->set_cur( (long)0 ); varCol->put( &( vecCol[0] ), nS ); varS->set_cur( (long)0 ); varS->put( &( vecS[0] ), nS ); ncMap.close(); // rval = mbCore->write_file(scripfile.c_str());MB_CHK_ERR(rval); } catch( std::exception& e ) { std::cout << " exception caught during tree initialization " << e.what() << std::endl; } delete mbCore; #ifdef MOAB_HAVE_MPI MPI_Finalize(); #endif exit( 0 ); }
void ReadFileMetaData | ( | std::string & | metaFilename, |
std::map< std::string, std::string > & | metadataVals | ||
) |
Definition at line 62 of file h5mtoscrip.cpp.
Referenced by main().
{ std::ifstream metafile; std::string line; metafile.open( metaFilename.c_str() ); metadataVals["Title"] = "MOAB-TempestRemap (MBTR) Offline Regridding Weight Converter (h5mtoscrip)"; std::string key, value; while( std::getline( metafile, line ) ) { size_t lastindex = line.find_last_of( "=" ); key = line.substr( 0, lastindex - 1 ); value = line.substr( lastindex + 2, line.length() ); metadataVals[std::string( key )] = std::string( value ); } metafile.close(); }