 |
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
|
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
1 /*
2 * addncdata.cpp
3 * this tool will take an existing h5m file and add data from some chunk description files
4 * generated from e3sm runs;
5 * will support mainly showing the chunks in ViSit
6 *
7 * example of usage:
8 * ./mbaddchunk -i penta3d.h5m -n chunks_on_proc.txt -o penta3d_ch.h5m
9 *
10 *
11 * Basically, will output a new h5m file (penta3d_ch.h5m), which has 2 extra tags, corresponding to
12 * the chunks number and processors it sits on
13 *
14 *
15 * file penta3d.h5m is obtained from the pentagons file, with a python script
16 *
17 */
18
19 #include "moab/ProgOptions.hpp"
20 #include "moab/Core.hpp"
21
22 #include <cmath>
23 #include <sstream>
24 #include <iostream>
25 #include <iomanip>
26 #include <fstream>
27
28 using namespace moab;
29 using namespace std;
30
31 int main( int argc, char* argv[] )
32 {
33
34 ProgOptions opts;
35
36 std::string inputfile( "penta3d.h5m" ), outfile( "penta3d_ch.h5m" ), chunkfile_name, gsmapfile;
37
38 opts.addOpt< std::string >( "input,i", "input mesh filename", &inputfile );
39 opts.addOpt< std::string >( "chunkFile,n", "chunk file from cam run", &chunkfile_name );
40 opts.addOpt< std::string >( "gsMAPfile,g", "gsmap file", &gsmapfile );
41
42 opts.addOpt< std::string >( "output,o", "output mesh filename", &outfile );
43
44 opts.parseCommandLine( argc, argv );
45
46 ErrorCode rval;
47 Core* mb = new Core();
48
49 rval = mb->load_file( inputfile.c_str() );MB_CHK_SET_ERR( rval, "can't load input file" );
50
51 std::cout << " opened " << inputfile << " with initial h5m data.\n";
52 // open the netcdf file, and see if it has that variable we are looking for
53
54 Range nodes;
55 rval = mb->get_entities_by_dimension( 0, 0, nodes );MB_CHK_SET_ERR( rval, "can't get nodes" );
56
57 Range edges;
58 rval = mb->get_entities_by_dimension( 0, 1, edges );MB_CHK_SET_ERR( rval, "can't get edges" );
59
60 Range cells;
61 rval = mb->get_entities_by_dimension( 0, 2, cells );MB_CHK_SET_ERR( rval, "can't get cells" );
62
63 std::cout << " it has " << nodes.size() << " vertices " << edges.size() << " edges " << cells.size() << " cells\n";
64
65 // construct maps between global id and handles
66 std::map< int, EntityHandle > vGidHandle;
67 std::map< int, EntityHandle > eGidHandle;
68 std::map< int, EntityHandle > cGidHandle;
69 std::vector< int > gids;
70 Tag gid;
71 rval = mb->tag_get_handle( "GLOBAL_ID", gid );MB_CHK_SET_ERR( rval, "can't get global id tag" );
72 gids.resize( nodes.size() );
73 rval = mb->tag_get_data( gid, nodes, &gids[0] );MB_CHK_SET_ERR( rval, "can't get global id on vertices" );
74 int i = 0;
75 for( Range::iterator vit = nodes.begin(); vit != nodes.end(); vit++ )
76 {
77 vGidHandle[gids[i++]] = *vit;
78 }
79
80 gids.resize( edges.size() );
81 rval = mb->tag_get_data( gid, edges, &gids[0] );MB_CHK_SET_ERR( rval, "can't get global id on edges" );
82 i = 0;
83 for( Range::iterator vit = edges.begin(); vit != edges.end(); vit++ )
84 {
85 eGidHandle[gids[i++]] = *vit;
86 }
87
88 gids.resize( cells.size() );
89 rval = mb->tag_get_data( gid, cells, &gids[0] );MB_CHK_SET_ERR( rval, "can't get global id on cells" );
90 i = 0;
91 for( Range::iterator vit = cells.begin(); vit != cells.end(); vit++ )
92 {
93 cGidHandle[gids[i++]] = *vit;
94 }
95
96 if( chunkfile_name.length() > 0 )
97 {
98
99 // Open chunk file
100 ifstream inFile;
101
102 inFile.open( chunkfile_name.c_str() );
103 if( !inFile )
104 {
105 cout << "Unable to open chunk file";
106 exit( 1 ); // terminate with error
107 }
108 Tag pTag, cTag;
109 int def_val = -1;
110 rval = mb->tag_get_handle( "ProcID", 1, MB_TYPE_INTEGER, pTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define processor tag" );
111 rval = mb->tag_get_handle( "ChunkID", 1, MB_TYPE_INTEGER, cTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define chunk tag" );
112
113 int proc, lcid, ncols;
114 while( inFile >> proc )
115 {
116 inFile >> lcid >> ncols;
117 int Gid;
118 for( i = 0; i < ncols; i++ )
119 {
120 inFile >> Gid;
121 EntityHandle cell = cGidHandle[Gid];
122 rval = mb->tag_set_data( pTag, &cell, 1, &proc );MB_CHK_SET_ERR( rval, "can't set proc tag" );
123 rval = mb->tag_set_data( cTag, &cell, 1, &lcid );MB_CHK_SET_ERR( rval, "can't set chunk tag" );
124 }
125 }
126
127 inFile.close();
128 }
129
130 if( gsmapfile.length() > 0 )
131 {
132
133 // Open chunk file
134 ifstream inFile;
135
136 inFile.open( gsmapfile.c_str() );
137 if( !inFile )
138 {
139 cout << "Unable to open gsmap file";
140 exit( 1 ); // terminate with error
141 }
142 Tag pTag, cTag;
143 int def_val = -1;
144 std::string procTagName = gsmapfile + "_proc";
145 rval =
146 mb->tag_get_handle( procTagName.c_str(), 1, MB_TYPE_INTEGER, pTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define processor tag" );
147 std::string segTagName = gsmapfile + "_seg";
148 rval =
149 mb->tag_get_handle( segTagName.c_str(), 1, MB_TYPE_INTEGER, cTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define segment tag" );
150
151 int compid, ngseg, gsize;
152 inFile >> compid >> ngseg >> gsize;
153 for( i = 1; i <= ngseg; i++ )
154 {
155 int start, len, pe;
156 inFile >> start >> len >> pe;
157 int Gid;
158 for( int j = 0; j < len; j++ )
159 {
160 Gid = start + j;
161 EntityHandle cell = cGidHandle[Gid];
162 rval = mb->tag_set_data( pTag, &cell, 1, &pe );MB_CHK_SET_ERR( rval, "can't set proc tag" );
163 rval = mb->tag_set_data( cTag, &cell, 1, &i );MB_CHK_SET_ERR( rval, "can't set segment tag" );
164 }
165 }
166
167 inFile.close();
168 }
169
170 rval = mb->write_file( outfile.c_str() );MB_CHK_SET_ERR( rval, "can't write file" );
171 std::cout << " wrote file " << outfile << "\n";
172 return 0;
173 }
1.9.1.