1
6
7 #include "NCWriteMPAS.hpp"
8 #include "MBTagConventions.hpp"
9
10 namespace moab
11 {
12
13 NCWriteMPAS::~NCWriteMPAS()
14 {
15
16 }
17
18 ErrorCode NCWriteMPAS::collect_mesh_info()
19 {
20 Interface*& mbImpl = _writeNC->mbImpl;
21 std::vector< std::string >& dimNames = _writeNC->dimNames;
22 std::vector< int >& dimLens = _writeNC->dimLens;
23 Tag& mGlobalIdTag = _writeNC->mGlobalIdTag;
24
25 ErrorCode rval;
26
27
28 std::vector< std::string >::iterator vecIt;
29 if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "Time" ) ) != dimNames.end() )
30 tDim = vecIt - dimNames.begin();
31 else if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "time" ) ) != dimNames.end() )
32 tDim = vecIt - dimNames.begin();
33 else
34 {
35 MB_SET_ERR( MB_FAILURE, "Couldn't find 'Time' or 'time' dimension" );
36 }
37 nTimeSteps = dimLens[tDim];
38
39
40 if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "nVertLevels" ) ) != dimNames.end() )
41 levDim = vecIt - dimNames.begin();
42 else
43 {
44 MB_SET_ERR( MB_FAILURE, "Couldn't find 'nVertLevels' dimension" );
45 }
46 nLevels = dimLens[levDim];
47
48
49 rval = mbImpl->get_entities_by_dimension( _fileSet, 0, localVertsOwned );MB_CHK_SET_ERR( rval, "Trouble getting local vertices in current file set" );
50 assert( !localVertsOwned.empty() );
51
52
53 rval = mbImpl->get_entities_by_dimension( _fileSet, 1, localEdgesOwned );MB_CHK_SET_ERR( rval, "Trouble getting local edges in current file set" );
54
55
56
57 rval = mbImpl->get_entities_by_dimension( _fileSet, 2, localCellsOwned );MB_CHK_SET_ERR( rval, "Trouble getting local cells in current file set" );
58 assert( !localCellsOwned.empty() );
59
60 #ifdef MOAB_HAVE_MPI
61 bool& isParallel = _writeNC->isParallel;
62 if( isParallel )
63 {
64 ParallelComm*& myPcomm = _writeNC->myPcomm;
65 int rank = myPcomm->proc_config().proc_rank();
66 int procs = myPcomm->proc_config().proc_size();
67 if( procs > 1 )
68 {
69 #ifndef NDEBUG
70 unsigned int num_local_verts = localVertsOwned.size();
71 #endif
72 rval = myPcomm->filter_pstatus( localVertsOwned, PSTATUS_NOT_OWNED, PSTATUS_NOT );MB_CHK_SET_ERR( rval, "Trouble getting owned vertices in current file set" );
73
74
75
76 if( procs - 1 == rank )
77 assert( "PARALLEL_RESOLVE_SHARED_ENTS option is set" && localVertsOwned.size() < num_local_verts );
78
79 if( !localEdgesOwned.empty() )
80 {
81 rval = myPcomm->filter_pstatus( localEdgesOwned, PSTATUS_NOT_OWNED, PSTATUS_NOT );MB_CHK_SET_ERR( rval, "Trouble getting owned edges in current file set" );
82 }
83
84 rval = myPcomm->filter_pstatus( localCellsOwned, PSTATUS_NOT_OWNED, PSTATUS_NOT );MB_CHK_SET_ERR( rval, "Trouble getting owned cells in current file set" );
85 }
86 }
87 #endif
88
89 std::vector< int > gids( localVertsOwned.size() );
90 rval = mbImpl->tag_get_data( mGlobalIdTag, localVertsOwned, &gids[0] );MB_CHK_SET_ERR( rval, "Trouble getting global IDs on local vertices" );
91
92
93 std::copy( gids.rbegin(), gids.rend(), range_inserter( localGidVertsOwned ) );
94
95 if( !localEdgesOwned.empty() )
96 {
97 gids.resize( localEdgesOwned.size() );
98 rval = mbImpl->tag_get_data( mGlobalIdTag, localEdgesOwned, &gids[0] );MB_CHK_SET_ERR( rval, "Trouble getting global IDs on local edges" );
99
100
101 std::copy( gids.rbegin(), gids.rend(), range_inserter( localGidEdgesOwned ) );
102 }
103
104 gids.resize( localCellsOwned.size() );
105 rval = mbImpl->tag_get_data( mGlobalIdTag, localCellsOwned, &gids[0] );MB_CHK_SET_ERR( rval, "Trouble getting global IDs on local cells" );
106
107
108 std::copy( gids.rbegin(), gids.rend(), range_inserter( localGidCellsOwned ) );
109
110 return MB_SUCCESS;
111 }
112
113 ErrorCode NCWriteMPAS::collect_variable_data( std::vector< std::string >& var_names, std::vector< int >& tstep_nums )
114 {
115 NCWriteHelper::collect_variable_data( var_names, tstep_nums );
116
117 std::vector< std::string >& dimNames = _writeNC->dimNames;
118 std::vector< int >& dimLens = _writeNC->dimLens;
119
120
121 std::vector< unsigned int > opt_lev_dims;
122
123 unsigned int lev_idx;
124 std::vector< std::string >::iterator vecIt;
125
126
127 if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "nVertLevelsP1" ) ) != dimNames.end() )
128 {
129 lev_idx = vecIt - dimNames.begin();
130 opt_lev_dims.push_back( lev_idx );
131 }
132
133
134 if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "nVertLevelsP2" ) ) != dimNames.end() )
135 {
136 lev_idx = vecIt - dimNames.begin();
137 opt_lev_dims.push_back( lev_idx );
138 }
139
140
141 if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "nSoilLevels" ) ) != dimNames.end() )
142 {
143 lev_idx = vecIt - dimNames.begin();
144 opt_lev_dims.push_back( lev_idx );
145 }
146
147 std::map< std::string, WriteNC::VarData >& varInfo = _writeNC->varInfo;
148
149 for( size_t i = 0; i < var_names.size(); i++ )
150 {
151 std::string varname = var_names[i];
152 std::map< std::string, WriteNC::VarData >::iterator vit = varInfo.find( varname );
153 if( vit == varInfo.end() ) MB_SET_ERR( MB_FAILURE, "Can't find variable " << varname );
154
155 WriteNC::VarData& currentVarData = vit->second;
156 std::vector< int >& varDims = currentVarData.varDims;
157
158
159 if( localEdgesOwned.empty() && currentVarData.entLoc == WriteNC::ENTLOCEDGE ) continue;
160
161
162 if( std::find( varDims.begin(), varDims.end(), levDim ) == varDims.end() )
163 {
164 for( unsigned int j = 0; j < opt_lev_dims.size(); j++ )
165 {
166 if( std::find( varDims.begin(), varDims.end(), opt_lev_dims[j] ) != varDims.end() )
167 {
168 currentVarData.numLev = dimLens[opt_lev_dims[j]];
169 break;
170 }
171 }
172 }
173
174
175
176 if( WriteNC::ENTLOCSET == currentVarData.entLoc ) continue;
177
178
179 currentVarData.writeStarts.resize( 3 );
180 currentVarData.writeCounts.resize( 3 );
181 unsigned int dim_idx = 0;
182
183
184 if( currentVarData.has_tsteps )
185 {
186
187
188
189 assert( 3 == varDims.size() || 2 == varDims.size() );
190
191
192 assert( tDim == varDims[0] );
193
194 currentVarData.writeStarts[dim_idx] = 0;
195 currentVarData.writeCounts[dim_idx] = 1;
196 dim_idx++;
197 }
198 else
199 {
200
201
202
203 assert( 2 == varDims.size() || 1 == varDims.size() );
204 }
205
206
207 switch( currentVarData.entLoc )
208 {
209 case WriteNC::ENTLOCVERT:
210
211
212
213 currentVarData.writeStarts[dim_idx] = localGidVertsOwned[0] - 1;
214 currentVarData.writeCounts[dim_idx] = localGidVertsOwned.size();
215 break;
216 case WriteNC::ENTLOCFACE:
217
218
219
220 currentVarData.writeStarts[dim_idx] = localGidCellsOwned[0] - 1;
221 currentVarData.writeCounts[dim_idx] = localGidCellsOwned.size();
222 break;
223 case WriteNC::ENTLOCEDGE:
224
225
226
227 currentVarData.writeStarts[dim_idx] = localGidEdgesOwned[0] - 1;
228 currentVarData.writeCounts[dim_idx] = localGidEdgesOwned.size();
229 break;
230 default:
231 MB_SET_ERR( MB_FAILURE, "Unexpected entity location type for variable " << varname );
232 }
233 dim_idx++;
234
235
236
237 if( currentVarData.numLev > 0 )
238 {
239
240
241
242 assert( 3 == varDims.size() || 2 == varDims.size() );
243
244 currentVarData.writeStarts[dim_idx] = 0;
245 currentVarData.writeCounts[dim_idx] = currentVarData.numLev;
246 dim_idx++;
247 }
248 else
249 {
250
251
252
253 assert( 2 == varDims.size() || 1 == varDims.size() );
254 }
255
256
257 currentVarData.sz = 1;
258 for( std::size_t idx = 0; idx < dim_idx; idx++ )
259 currentVarData.sz *= currentVarData.writeCounts[idx];
260 }
261
262 return MB_SUCCESS;
263 }
264
265 ErrorCode NCWriteMPAS::write_nonset_variables( std::vector< WriteNC::VarData >& vdatas, std::vector< int >& tstep_nums )
266 {
267 Interface*& mbImpl = _writeNC->mbImpl;
268
269 int success;
270
271
272 for( unsigned int i = 0; i < vdatas.size(); i++ )
273 {
274 WriteNC::VarData& variableData = vdatas[i];
275
276
277 if( localEdgesOwned.empty() && variableData.entLoc == WriteNC::ENTLOCEDGE ) continue;
278
279
280 Range* pLocalEntsOwned = NULL;
281 Range* pLocalGidEntsOwned = NULL;
282 switch( variableData.entLoc )
283 {
284 case WriteNC::ENTLOCVERT:
285
286 pLocalEntsOwned = &localVertsOwned;
287 pLocalGidEntsOwned = &localGidVertsOwned;
288 break;
289 case WriteNC::ENTLOCEDGE:
290
291 pLocalEntsOwned = &localEdgesOwned;
292 pLocalGidEntsOwned = &localGidEdgesOwned;
293 break;
294 case WriteNC::ENTLOCFACE:
295
296 pLocalEntsOwned = &localCellsOwned;
297 pLocalGidEntsOwned = &localGidCellsOwned;
298 break;
299 default:
300 MB_SET_ERR( MB_FAILURE, "Unexpected entity location type for variable " << variableData.varName );
301 }
302
303 unsigned int num_timesteps;
304 unsigned int ents_idx = 0;
305 if( variableData.has_tsteps )
306 {
307
308
309
310 num_timesteps = tstep_nums.size();
311 ents_idx++;
312 }
313 else
314 {
315
316
317
318 num_timesteps = 1;
319 }
320
321 unsigned int num_lev;
322 if( variableData.numLev > 0 )
323 {
324
325
326
327 num_lev = variableData.numLev;
328 }
329 else
330 {
331
332
333
334 num_lev = 1;
335 }
336
337 for( unsigned int t = 0; t < num_timesteps; t++ )
338 {
339
340
341
342 if( tDim == variableData.varDims[0] ) variableData.writeStarts[0] = t;
343 std::vector< double > tag_data( pLocalEntsOwned->size() * num_lev );
344 ErrorCode rval = mbImpl->tag_get_data( variableData.varTags[t], *pLocalEntsOwned, &tag_data[0] );MB_CHK_SET_ERR( rval, "Trouble getting tag data on owned entities" );
345
346 #ifdef MOAB_HAVE_PNETCDF
347 size_t nb_writes = pLocalGidEntsOwned->psize();
348 std::vector< int > requests( nb_writes ), statuss( nb_writes );
349 size_t idxReq = 0;
350 #endif
351
352
353
354 switch( variableData.varDataType )
355 {
356 case NC_DOUBLE: {
357 size_t indexInDoubleArray = 0;
358 size_t ic = 0;
359 for( Range::pair_iterator pair_iter = pLocalGidEntsOwned->pair_begin();
360 pair_iter != pLocalGidEntsOwned->pair_end(); ++pair_iter, ic++ )
361 {
362 EntityHandle starth = pair_iter->first;
363 EntityHandle endh = pair_iter->second;
364 variableData.writeStarts[ents_idx] = (NCDF_SIZE)( starth - 1 );
365 variableData.writeCounts[ents_idx] = (NCDF_SIZE)( endh - starth + 1 );
366
367
368 #ifdef MOAB_HAVE_PNETCDF
369
370 success =
371 NCFUNCREQP( _vara_double )( _fileId, variableData.varId, &( variableData.writeStarts[0] ),
372 &( variableData.writeCounts[0] ),
373 &( tag_data[indexInDoubleArray] ), &requests[idxReq++] );
374 #else
375 success = NCFUNCAP(
376 _vara_double )( _fileId, variableData.varId, &( variableData.writeStarts[0] ),
377 &( variableData.writeCounts[0] ), &( tag_data[indexInDoubleArray] ) );
378 #endif
379 if( success )
380 MB_SET_ERR( MB_FAILURE,
381 "Failed to write double data in a loop for variable " << variableData.varName );
382
383
384 indexInDoubleArray += ( endh - starth + 1 ) * num_lev;
385 }
386 assert( ic == pLocalGidEntsOwned->psize() );
387 #ifdef MOAB_HAVE_PNETCDF
388 success = ncmpi_wait_all( _fileId, requests.size(), &requests[0], &statuss[0] );
389 if( success ) MB_SET_ERR( MB_FAILURE, "Failed on wait_all" );
390 #endif
391 break;
392 }
393 default:
394 MB_SET_ERR( MB_NOT_IMPLEMENTED, "Writing non-double data is not implemented yet" );
395 }
396 }
397 }
398
399 return MB_SUCCESS;
400 }
401
402 }