Branch data Line data Source code
1 : : /*
2 : : * NCWriteHOMME.cpp
3 : : *
4 : : * Created on: April 9, 2014
5 : : */
6 : :
7 : : #include "NCWriteHOMME.hpp"
8 : : #include "MBTagConventions.hpp"
9 : :
10 : : namespace moab
11 : : {
12 : :
13 : 0 : NCWriteHOMME::~NCWriteHOMME()
14 : : {
15 : : // TODO Auto-generated destructor stub
16 [ # # ]: 0 : }
17 : :
18 : 0 : ErrorCode NCWriteHOMME::collect_mesh_info()
19 : : {
20 : 0 : Interface*& mbImpl = _writeNC->mbImpl;
21 : 0 : std::vector< std::string >& dimNames = _writeNC->dimNames;
22 : 0 : std::vector< int >& dimLens = _writeNC->dimLens;
23 : 0 : Tag& mGlobalIdTag = _writeNC->mGlobalIdTag;
24 : :
25 : : ErrorCode rval;
26 : :
27 : : // Look for time dimension
28 : 0 : std::vector< std::string >::iterator vecIt;
29 [ # # ][ # # ]: 0 : if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "time" ) ) != dimNames.end() )
[ # # ]
30 [ # # ]: 0 : tDim = vecIt - dimNames.begin();
31 : : else
32 : : {
33 [ # # ][ # # ]: 0 : MB_SET_ERR( MB_FAILURE, "Couldn't find 'time' dimension" );
[ # # ][ # # ]
[ # # ]
34 : : }
35 [ # # ]: 0 : nTimeSteps = dimLens[tDim];
36 : :
37 : : // Get number of levels
38 [ # # ][ # # ]: 0 : if( ( vecIt = std::find( dimNames.begin(), dimNames.end(), "lev" ) ) != dimNames.end() )
[ # # ]
39 [ # # ]: 0 : levDim = vecIt - dimNames.begin();
40 : : else
41 : : {
42 [ # # ][ # # ]: 0 : MB_SET_ERR( MB_FAILURE, "Couldn't find 'lev' dimension" );
[ # # ][ # # ]
[ # # ]
43 : : }
44 [ # # ]: 0 : nLevels = dimLens[levDim];
45 : :
46 : : // Get local vertices
47 [ # # ][ # # ]: 0 : rval = mbImpl->get_entities_by_dimension( _fileSet, 0, localVertsOwned );MB_CHK_SET_ERR( rval, "Trouble getting local vertices in current file set" );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
48 [ # # ][ # # ]: 0 : assert( !localVertsOwned.empty() );
49 : :
50 : : #ifdef MOAB_HAVE_MPI
51 : 0 : bool& isParallel = _writeNC->isParallel;
52 [ # # ]: 0 : if( isParallel )
53 : : {
54 : 0 : ParallelComm*& myPcomm = _writeNC->myPcomm;
55 [ # # ][ # # ]: 0 : int rank = myPcomm->proc_config().proc_rank();
56 [ # # ][ # # ]: 0 : int procs = myPcomm->proc_config().proc_size();
57 [ # # ]: 0 : if( procs > 1 )
58 : : {
59 : : #ifndef NDEBUG
60 [ # # ]: 0 : unsigned int num_local_verts = localVertsOwned.size();
61 : : #endif
62 [ # # ][ # # ]: 0 : rval = myPcomm->filter_pstatus( localVertsOwned, PSTATUS_NOT_OWNED, PSTATUS_NOT );MB_CHK_SET_ERR( rval, "Trouble getting owned vertices in current set" );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
63 : :
64 : : // Assume that PARALLEL_RESOLVE_SHARED_ENTS option is set
65 : : // Verify that not all local vertices are owned by the last processor
66 [ # # ]: 0 : if( procs - 1 == rank )
67 [ # # ][ # # ]: 0 : assert( "PARALLEL_RESOLVE_SHARED_ENTS option is set" && localVertsOwned.size() < num_local_verts );
68 : : }
69 : : }
70 : : #endif
71 : :
72 [ # # ][ # # ]: 0 : std::vector< int > gids( localVertsOwned.size() );
73 [ # # ][ # # ]: 0 : rval = mbImpl->tag_get_data( mGlobalIdTag, localVertsOwned, &gids[0] );MB_CHK_SET_ERR( rval, "Trouble getting global IDs on local vertices" );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
74 : :
75 : : // Get localGidVertsOwned
76 [ # # ][ # # ]: 0 : std::copy( gids.rbegin(), gids.rend(), range_inserter( localGidVertsOwned ) );
77 : :
78 : 0 : return MB_SUCCESS;
79 : : }
80 : :
81 : 0 : ErrorCode NCWriteHOMME::collect_variable_data( std::vector< std::string >& var_names, std::vector< int >& tstep_nums )
82 : : {
83 : 0 : NCWriteHelper::collect_variable_data( var_names, tstep_nums );
84 : :
85 : 0 : std::map< std::string, WriteNC::VarData >& varInfo = _writeNC->varInfo;
86 : :
87 [ # # ]: 0 : for( size_t i = 0; i < var_names.size(); i++ )
88 : : {
89 [ # # ][ # # ]: 0 : std::string varname = var_names[i];
90 [ # # ]: 0 : std::map< std::string, WriteNC::VarData >::iterator vit = varInfo.find( varname );
91 [ # # ][ # # ]: 0 : if( vit == varInfo.end() ) MB_SET_ERR( MB_FAILURE, "Can't find variable " << varname );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
92 : : ;
93 : :
94 [ # # ]: 0 : WriteNC::VarData& currentVarData = vit->second;
95 : : #ifndef NDEBUG
96 : 0 : std::vector< int >& varDims = currentVarData.varDims;
97 : : #endif
98 : :
99 : : // Skip set variables, which were already processed in
100 : : // NCWriteHelper::collect_variable_data()
101 [ # # ]: 0 : if( WriteNC::ENTLOCSET == currentVarData.entLoc ) continue;
102 : :
103 : : // Set up writeStarts and writeCounts (maximum number of dimensions is 3)
104 [ # # ]: 0 : currentVarData.writeStarts.resize( 3 );
105 [ # # ]: 0 : currentVarData.writeCounts.resize( 3 );
106 : 0 : unsigned int dim_idx = 0;
107 : :
108 : : // First: time
109 [ # # ]: 0 : if( currentVarData.has_tsteps )
110 : : {
111 : : // Non-set variables with timesteps
112 : : // 3 dimensions like (time, lev, ncol)
113 : : // 2 dimensions like (time, ncol)
114 [ # # ][ # # ]: 0 : assert( 3 == varDims.size() || 2 == varDims.size() );
115 : :
116 : : // Time should be the first dimension
117 [ # # ][ # # ]: 0 : assert( tDim == varDims[0] );
118 : :
119 [ # # ]: 0 : currentVarData.writeStarts[dim_idx] = 0; // This value is timestep dependent, will be set later
120 [ # # ]: 0 : currentVarData.writeCounts[dim_idx] = 1;
121 : 0 : dim_idx++;
122 : : }
123 : : else
124 : : {
125 : : // Non-set variables without timesteps
126 : : // 2 dimensions like (lev, ncol)
127 : : // 1 dimension like (ncol)
128 [ # # ][ # # ]: 0 : assert( 2 == varDims.size() || 1 == varDims.size() );
129 : : }
130 : :
131 : : // Next: lev
132 [ # # ]: 0 : if( currentVarData.numLev > 0 )
133 : : {
134 : : // Non-set variables with levels
135 : : // 3 dimensions like (time, lev, ncol)
136 : : // 2 dimensions like (lev, ncol)
137 [ # # ][ # # ]: 0 : assert( 3 == varDims.size() || 2 == varDims.size() );
138 : :
139 [ # # ]: 0 : currentVarData.writeStarts[dim_idx] = 0;
140 [ # # ]: 0 : currentVarData.writeCounts[dim_idx] = currentVarData.numLev;
141 : 0 : dim_idx++;
142 : : }
143 : : else
144 : : {
145 : : // Non-set variables without levels
146 : : // 2 dimensions like (time, ncol)
147 : : // 1 dimension like (ncol)
148 [ # # ][ # # ]: 0 : assert( 2 == varDims.size() || 1 == varDims.size() );
149 : : }
150 : :
151 : : // Finally: ncol
152 [ # # ]: 0 : switch( currentVarData.entLoc )
153 : : {
154 : : case WriteNC::ENTLOCVERT:
155 : : // Vertices
156 : : // Start from the first localGidVerts
157 : : // Actually, this will be reset later for writing
158 [ # # ][ # # ]: 0 : currentVarData.writeStarts[dim_idx] = localGidVertsOwned[0] - 1;
159 [ # # ][ # # ]: 0 : currentVarData.writeCounts[dim_idx] = localGidVertsOwned.size();
160 : 0 : break;
161 : : default:
162 [ # # ][ # # ]: 0 : MB_SET_ERR( MB_FAILURE, "Unexpected entity location type for variable " << varname );
[ # # ][ # # ]
[ # # ][ # # ]
163 : : }
164 : 0 : dim_idx++;
165 : :
166 : : // Get variable size
167 : 0 : currentVarData.sz = 1;
168 [ # # ]: 0 : for( std::size_t idx = 0; idx < dim_idx; idx++ )
[ # # # ]
169 [ # # ]: 0 : currentVarData.sz *= currentVarData.writeCounts[idx];
170 : 0 : }
171 : :
172 : 0 : return MB_SUCCESS;
173 : : }
174 : :
175 : 0 : ErrorCode NCWriteHOMME::write_nonset_variables( std::vector< WriteNC::VarData >& vdatas,
176 : : std::vector< int >& tstep_nums )
177 : : {
178 : 0 : Interface*& mbImpl = _writeNC->mbImpl;
179 : :
180 : : int success;
181 : 0 : int num_local_verts_owned = localVertsOwned.size();
182 : :
183 : : // For each indexed variable tag, write a time step data
184 [ # # ]: 0 : for( unsigned int i = 0; i < vdatas.size(); i++ )
185 : : {
186 : 0 : WriteNC::VarData& variableData = vdatas[i];
187 : :
188 : : // Assume this variable is on vertices for the time being
189 [ # # ]: 0 : switch( variableData.entLoc )
190 : : {
191 : : case WriteNC::ENTLOCVERT:
192 : : // Vertices
193 : 0 : break;
194 : : default:
195 [ # # ][ # # ]: 0 : MB_SET_ERR( MB_FAILURE, "Unexpected entity location type for variable " << variableData.varName );
[ # # ][ # # ]
[ # # ][ # # ]
196 : : }
197 : :
198 : : unsigned int num_timesteps;
199 : 0 : unsigned int ncol_idx = 0;
200 [ # # ]: 0 : if( variableData.has_tsteps )
201 : : {
202 : : // Non-set variables with timesteps
203 : : // 3 dimensions like (time, lev, ncol)
204 : : // 2 dimensions like (time, ncol)
205 : 0 : num_timesteps = tstep_nums.size();
206 : 0 : ncol_idx++;
207 : : }
208 : : else
209 : : {
210 : : // Non-set variables without timesteps
211 : : // 2 dimensions like (lev, ncol)
212 : : // 1 dimension like (ncol)
213 : 0 : num_timesteps = 1;
214 : : }
215 : :
216 : : unsigned int num_lev;
217 [ # # ]: 0 : if( variableData.numLev > 0 )
218 : : {
219 : : // Non-set variables with levels
220 : : // 3 dimensions like (time, lev, ncol)
221 : : // 2 dimensions like (lev, ncol)
222 : 0 : num_lev = variableData.numLev;
223 : 0 : ncol_idx++;
224 : : }
225 : : else
226 : : {
227 : : // Non-set variables without levels
228 : : // 2 dimensions like (time, ncol)
229 : : // 1 dimension like (ncol)
230 : 0 : num_lev = 1;
231 : : }
232 : :
233 : : // At each timestep, we need to transpose tag format (ncol, lev) back
234 : : // to NC format (lev, ncol) for writing
235 [ # # ]: 0 : for( unsigned int t = 0; t < num_timesteps; t++ )
236 : : {
237 : : // We will write one time step, and count will be one; start will be different
238 : : // Use tag_get_data instead of tag_iterate to copy tag data, as localVertsOwned
239 : : // might not be contiguous. We should also transpose for level so that means
240 : : // deep copy for transpose
241 [ # # ][ # # ]: 0 : if( tDim == variableData.varDims[0] ) variableData.writeStarts[0] = t; // This is start for time
[ # # ]
242 [ # # ]: 0 : std::vector< double > tag_data( num_local_verts_owned * num_lev );
243 [ # # ][ # # ]: 0 : ErrorCode rval = mbImpl->tag_get_data( variableData.varTags[t], localVertsOwned, &tag_data[0] );MB_CHK_SET_ERR( rval, "Trouble getting tag data on owned vertices" );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
244 : :
245 : : #ifdef MOAB_HAVE_PNETCDF
246 : : size_t nb_writes = localGidVertsOwned.psize();
247 : : std::vector< int > requests( nb_writes ), statuss( nb_writes );
248 : : size_t idxReq = 0;
249 : : #endif
250 : :
251 : : // Now transpose and write copied tag data
252 : : // Use nonblocking put (request aggregation)
253 [ # # ]: 0 : switch( variableData.varDataType )
254 : : {
255 : : case NC_DOUBLE: {
256 [ # # ]: 0 : std::vector< double > tmpdoubledata( num_local_verts_owned * num_lev );
257 [ # # ]: 0 : if( num_lev > 1 )
258 : : {
259 : : // Transpose (ncol, lev) back to (lev, ncol)
260 : : // Note, num_local_verts_owned is not used by jik_to_kji_stride()
261 [ # # ][ # # ]: 0 : jik_to_kji_stride( num_local_verts_owned, 1, num_lev, &tmpdoubledata[0], &tag_data[0],
262 [ # # ]: 0 : localGidVertsOwned );
263 : : }
264 : :
265 : 0 : size_t indexInDoubleArray = 0;
266 : 0 : size_t ic = 0;
267 [ # # ][ # # ]: 0 : for( Range::pair_iterator pair_iter = localGidVertsOwned.pair_begin();
[ # # ][ # # ]
268 [ # # ]: 0 : pair_iter != localGidVertsOwned.pair_end(); ++pair_iter, ic++ )
269 : : {
270 [ # # ]: 0 : EntityHandle starth = pair_iter->first;
271 [ # # ]: 0 : EntityHandle endh = pair_iter->second;
272 [ # # ]: 0 : variableData.writeStarts[ncol_idx] = ( NCDF_SIZE )( starth - 1 );
273 [ # # ]: 0 : variableData.writeCounts[ncol_idx] = ( NCDF_SIZE )( endh - starth + 1 );
274 : :
275 : : // Do a partial write, in each subrange
276 : : #ifdef MOAB_HAVE_PNETCDF
277 : : // Wait outside this loop
278 : : success =
279 : : NCFUNCREQP( _vara_double )( _fileId, variableData.varId, &( variableData.writeStarts[0] ),
280 : : &( variableData.writeCounts[0] ),
281 : : &( tmpdoubledata[indexInDoubleArray] ), &requests[idxReq++] );
282 : : #else
283 : : success = NCFUNCAP(
284 [ # # ]: 0 : _vara_double )( _fileId, variableData.varId, &( variableData.writeStarts[0] ),
285 [ # # ][ # # ]: 0 : &( variableData.writeCounts[0] ), &( tmpdoubledata[indexInDoubleArray] ) );
[ # # ]
286 : : #endif
287 [ # # ]: 0 : if( success )
288 [ # # ][ # # ]: 0 : MB_SET_ERR( MB_FAILURE,
[ # # ][ # # ]
[ # # ][ # # ]
289 : : "Failed to write double data in a loop for variable " << variableData.varName );
290 : : // We need to increment the index in double array for the
291 : : // next subrange
292 : 0 : indexInDoubleArray += ( endh - starth + 1 ) * num_lev;
293 : : }
294 [ # # ][ # # ]: 0 : assert( ic == localGidVertsOwned.psize() );
295 : : #ifdef MOAB_HAVE_PNETCDF
296 : : success = ncmpi_wait_all( _fileId, requests.size(), &requests[0], &statuss[0] );
297 : : if( success ) MB_SET_ERR( MB_FAILURE, "Failed on wait_all" );
298 : : #endif
299 [ # # ]: 0 : break;
300 : : }
301 : : default:
302 [ # # ][ # # ]: 0 : MB_SET_ERR( MB_NOT_IMPLEMENTED, "Writing non-double data is not implemented yet" );
[ # # ][ # # ]
[ # # ][ # # ]
303 : : }
304 : 0 : }
305 : : }
306 : :
307 : 0 : return MB_SUCCESS;
308 : : }
309 : :
310 [ + - ][ + - ]: 228 : } /* namespace moab */
|