/[escript]/trunk/escript/src/DataTagged.cpp
ViewVC logotype

Contents of /trunk/escript/src/DataTagged.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1628 - (show annotations)
Fri Jul 11 13:12:46 2008 UTC (11 years, 2 months ago) by phornby
File size: 23289 byte(s)

Merge in /branches/windows_from_1456_trunk_1620_merged_in branch.

You will find a preserved pre-merge trunk in tags under tags/trunk_at_1625.
That will be useful for diffing & checking on my stupidity.

Here is a list of the conflicts and their resolution at this
point in time.


=================================================================================
(LLWS == looks like white space).

finley/src/Assemble_addToSystemMatrix.c - resolve to branch - unused var. may be wrong.....
finley/src/CPPAdapter/SystemMatrixAdapter.cpp - resolve to branch - LLWS
finley/src/CPPAdapter/MeshAdapter.cpp - resolve to branch - LLWS
paso/src/PCG.c - resolve to branch - unused var fixes.
paso/src/SolverFCT.c - resolve to branch - LLWS
paso/src/FGMRES.c - resolve to branch - LLWS
paso/src/Common.h - resolve to trunk version. It's omp.h's include... not sure it's needed,
but for the sake of saftey.....
paso/src/Functions.c - resolve to branch version, indentation/tab removal and return error
on bad unimplemented Paso_FunctionCall.
paso/src/SolverFCT_solve.c - resolve to branch version, unused vars
paso/src/SparseMatrix_MatrixVector.c - resolve to branch version, unused vars.
escript/src/Utils.cpp - resloved to branch, needs WinSock2.h
escript/src/DataExpanded.cpp - resolved to branch version - LLWS
escript/src/DataFactory.cpp - resolve to branch version
=================================================================================

This currently passes tests on linux (debian), but is not checked on windows or Altix yet.

This checkin is to make a trunk I can check out for windows to do tests on it.

Known outstanding problem is in the operator=() method of exceptions
causing warning messages on the intel compilers.

May the God of doughnuts have mercy on my soul.


1
2 /* $Id$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 #include "DataTagged.h"
17 #include "esysUtils/esys_malloc.h"
18
19 #include "DataConstant.h"
20 #include "DataException.h"
21 #ifdef USE_NETCDF
22 #include <netcdfcpp.h>
23 #endif
24
25 using namespace std;
26
27 namespace escript {
28
29 DataTagged::DataTagged()
30 : DataAbstract(FunctionSpace())
31 {
32 // default constructor
33
34 // create a scalar default value
35 m_data.resize(1,0.,1);
36 DataArrayView temp(m_data,DataArrayView::ShapeType());
37 setPointDataView(temp);
38 }
39
40 DataTagged::DataTagged(const TagListType& tagKeys,
41 const ValueListType& values,
42 const DataArrayView& defaultValue,
43 const FunctionSpace& what)
44 : DataAbstract(what)
45 {
46 // constructor
47
48 // initialise the array of data values
49 // the default value is always the first item in the values list
50 int len = defaultValue.noValues();
51 m_data.resize(len,0.,len);
52 for (int i=0; i<defaultValue.noValues(); i++) {
53 m_data[i]=defaultValue.getData(i);
54 }
55
56 // create the data view
57 DataArrayView temp(m_data,defaultValue.getShape());
58 setPointDataView(temp);
59
60 // add remaining tags and values
61 addTaggedValues(tagKeys,values);
62 }
63
64 DataTagged::DataTagged(const FunctionSpace& what,
65 const DataArrayView::ShapeType &shape,
66 const int tags[],
67 const ValueType& data)
68 : DataAbstract(what)
69 {
70 // alternative constructor
71 // not unit_tested tested yet
72
73 // copy the data
74 m_data=data;
75
76 // create the view of the data
77 DataArrayView tempView(m_data,shape);
78 setPointDataView(tempView);
79
80 // create the tag lookup map
81 for (int sampleNo=0; sampleNo<getNumSamples(); sampleNo++) {
82 m_offsetLookup.insert(DataMapType::value_type(sampleNo,tags[sampleNo]));
83 }
84 }
85
86 DataTagged::DataTagged(const FunctionSpace& what,
87 const DataArrayView::ShapeType &shape,
88 const TagListType& tags,
89 const ValueType& data)
90 : DataAbstract(what)
91 {
92 // alternative constructor
93 // not unit_tested tested yet
94
95 // copy the data
96 m_data=data;
97
98 // create the view of the data
99 DataArrayView tempView(m_data,shape);
100 setPointDataView(tempView);
101
102 // create the tag lookup map
103 for (int sampleNo=0; sampleNo<getNumSamples(); sampleNo++) {
104 m_offsetLookup.insert(DataMapType::value_type(sampleNo,tags[sampleNo]));
105 }
106 }
107
108
109 DataTagged::DataTagged(const DataTagged& other)
110 : DataAbstract(other.getFunctionSpace()),
111 m_data(other.m_data),
112 m_offsetLookup(other.m_offsetLookup)
113 {
114 // copy constructor
115
116 // create the data view
117 DataArrayView temp(m_data,other.getPointDataView().getShape());
118 setPointDataView(temp);
119 }
120
121 DataTagged::DataTagged(const DataConstant& other)
122 : DataAbstract(other.getFunctionSpace())
123 {
124 // copy constructor
125
126 // fill the default value with the constant value item from "other"
127 const DataArrayView& value=other.getPointDataView();
128 int len = value.noValues();
129 m_data.resize(len,0.,len);
130 for (int i=0; i<value.noValues(); i++) {
131 m_data[i]=value.getData(i);
132 }
133
134 // create the data view
135 DataArrayView temp(m_data,value.getShape());
136 setPointDataView(temp);
137 }
138
139 DataAbstract*
140 DataTagged::getSlice(const DataArrayView::RegionType& region) const
141 {
142 return new DataTagged(*this, region);
143 }
144
145 DataTagged::DataTagged(const DataTagged& other,
146 const DataArrayView::RegionType& region)
147 : DataAbstract(other.getFunctionSpace())
148 {
149 // slice constructor
150
151 // get the shape of the slice to copy from other
152 DataArrayView::ShapeType regionShape(DataArrayView::getResultSliceShape(region));
153 DataArrayView::RegionLoopRangeType regionLoopRange=getSliceRegionLoopRange(region);
154
155 // allocate enough space in this for all values
156 // (need to add one to allow for the default value)
157 int len = DataArrayView::noValues(regionShape)*(other.m_offsetLookup.size()+1);
158 m_data.resize(len,0.0,len);
159
160 // create the data view
161 DataArrayView temp(m_data,regionShape);
162 setPointDataView(temp);
163
164 // copy the default value from other to this
165 getDefaultValue().copySlice(other.getDefaultValue(), regionLoopRange);
166
167 // loop through the tag values copying these
168 DataMapType::const_iterator pos;
169 DataArrayView::ValueType::size_type tagOffset=getPointDataView().noValues();
170 for (pos=other.m_offsetLookup.begin();pos!=other.m_offsetLookup.end();pos++){
171 getPointDataView().copySlice(tagOffset,other.getPointDataView(),pos->second,regionLoopRange);
172 m_offsetLookup.insert(DataMapType::value_type(pos->first,tagOffset));
173 tagOffset+=getPointDataView().noValues();
174 }
175 }
176
177 void
178 DataTagged::setSlice(const DataAbstract* other,
179 const DataArrayView::RegionType& region)
180 {
181
182 // other must be another DataTagged object
183 // Data:setSlice implementation should ensure this
184 const DataTagged* otherTemp=dynamic_cast<const DataTagged*>(other);
185 if (otherTemp==0) {
186 throw DataException("Programming error - casting to DataTagged.");
187 }
188
189 // determine shape of the specified region
190 DataArrayView::ShapeType regionShape(DataArrayView::getResultSliceShape(region));
191
192 // modify region specification as needed to match rank of this object
193 DataArrayView::RegionLoopRangeType regionLoopRange=getSliceRegionLoopRange(region);
194
195 // ensure rank/shape of this object is compatible with specified region
196 if (getPointDataView().getRank()!=region.size()) {
197 throw DataException("Error - Invalid slice region.");
198 }
199 if (otherTemp->getPointDataView().getRank()>0 && !other->getPointDataView().checkShape(regionShape)) {
200 throw DataException (other->getPointDataView().createShapeErrorMessage(
201 "Error - Couldn't copy slice due to shape mismatch.",regionShape));
202 }
203
204 // copy slice from other default value to this default value
205 getDefaultValue().copySliceFrom(otherTemp->getDefaultValue(), regionLoopRange);
206
207 // loop through tag values in other, adding any which aren't in this, using default value
208 DataMapType::const_iterator pos;
209 for (pos=otherTemp->m_offsetLookup.begin();pos!=otherTemp->m_offsetLookup.end();pos++) {
210 if (!isCurrentTag(pos->first)) {
211 addTaggedValue(pos->first,getDefaultValue());
212 }
213 }
214
215 // loop through the tag values copying slices from other to this
216 for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++) {
217 getDataPointByTag(pos->first).copySliceFrom(otherTemp->getDataPointByTag(pos->first), regionLoopRange);
218 }
219
220 }
221
222 int
223 DataTagged::getTagNumber(int dpno)
224 {
225 //
226 // Get the number of samples and data-points per sample
227 int numSamples = getNumSamples();
228 int numDataPointsPerSample = getNumDPPSample();
229 int numDataPoints = numSamples * numDataPointsPerSample;
230
231 if (numDataPointsPerSample==0) {
232 throw DataException("DataTagged::getTagNumber error: no data-points associated with this object.");
233 }
234
235 if (dpno<0 || dpno>numDataPoints-1) {
236 throw DataException("DataTagged::getTagNumber error: invalid data-point number supplied.");
237 }
238
239 //
240 // Determine the sample number which corresponds to this data-point number
241 int sampleNo = dpno / numDataPointsPerSample;
242
243 //
244 // Determine the tag number which corresponds to this sample number
245 int tagNo = getFunctionSpace().getTagFromSampleNo(sampleNo);
246
247 //
248 // return the tag number
249 return(tagNo);
250 }
251
252 void
253 DataTagged::setTaggedValues(const TagListType& tagKeys,
254 const ValueListType& values)
255 {
256 addTaggedValues(tagKeys,values);
257 }
258
259 void
260 DataTagged::setTaggedValue(int tagKey,
261 const DataArrayView& value)
262 {
263 if (!getPointDataView().checkShape(value.getShape())) {
264 throw DataException(getPointDataView().createShapeErrorMessage(
265 "Error - Cannot setTaggedValue due to shape mismatch.", value.getShape()));
266 }
267 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
268 if (pos==m_offsetLookup.end()) {
269 // tag couldn't be found so use addTaggedValue
270 addTaggedValue(tagKey,value);
271 } else {
272 // copy the values into the data array at the offset determined by m_offsetLookup
273 int offset=pos->second;
274 for (int i=0; i<getPointDataView().noValues(); i++) {
275 m_data[offset+i]=value.getData(i);
276 }
277 }
278 }
279
280 void
281 DataTagged::addTaggedValues(const TagListType& tagKeys,
282 const ValueListType& values)
283 {
284 if (values.size()==0) {
285 // copy the current default value for each of the tags
286 TagListType::const_iterator iT;
287 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
288 // the point data view for DataTagged points at the default value
289 addTaggedValue(*iT,getPointDataView());
290 }
291 } else if (values.size()==1 && tagKeys.size()>1) {
292 // assume the one given value will be used for all tag values
293 TagListType::const_iterator iT;
294 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
295 addTaggedValue(*iT,values[0]);
296 }
297 } else {
298 if (tagKeys.size()!=values.size()) {
299 stringstream temp;
300 temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
301 << " doesn't match number of values: " << values.size();
302 throw DataException(temp.str());
303 } else {
304 unsigned int i;
305 for (i=0;i<tagKeys.size();i++) {
306 addTaggedValue(tagKeys[i],values[i]);
307 }
308 }
309 }
310 }
311
312 void
313 DataTagged::addTaggedValue(int tagKey,
314 const DataArrayView& value)
315 {
316 if (!getPointDataView().checkShape(value.getShape())) {
317 throw DataException(getPointDataView().createShapeErrorMessage(
318 "Error - Cannot addTaggedValue due to shape mismatch.", value.getShape()));
319 }
320 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
321 if (pos!=m_offsetLookup.end()) {
322 // tag already exists so use setTaggedValue
323 setTaggedValue(tagKey,value);
324 } else {
325 // save the key and the location of its data in the lookup tab
326 m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
327 // add the data given in "value" at the end of m_data
328 // need to make a temp copy of m_data, resize m_data, then copy
329 // all the old values plus the value to be added back into m_data
330 ValueType m_data_temp(m_data);
331 int oldSize=m_data.size();
332 int newSize=m_data.size()+value.noValues();
333 m_data.resize(newSize,0.,newSize);
334 for (int i=0;i<oldSize;i++) {
335 m_data[i]=m_data_temp[i];
336 }
337 for (int i=0;i<value.noValues();i++) {
338 m_data[oldSize+i]=value.getData(i);
339 }
340 }
341 }
342
343 double*
344 DataTagged::getSampleDataByTag(int tag)
345 {
346 DataMapType::iterator pos(m_offsetLookup.find(tag));
347 if (pos==m_offsetLookup.end()) {
348 // tag couldn't be found so return the default value
349 return &(m_data[0]);
350 } else {
351 // return the data-point corresponding to the given tag
352 return &(m_data[pos->second]);
353 }
354 }
355
356 string
357 DataTagged::toString() const
358 {
359 stringstream temp;
360 DataMapType::const_iterator i;
361 temp << "Tag(Default)" << endl;
362 temp << getDefaultValue().toString() << endl;
363 // create a temporary view as the offset will be changed
364 DataArrayView tempView(getPointDataView().getData(), getPointDataView().getShape());
365 for (i=m_offsetLookup.begin();i!=m_offsetLookup.end();++i) {
366 temp << "Tag(" << i->first << ")" << endl;
367 tempView.setOffset(i->second);
368 temp << tempView.toString() << endl;
369 }
370 return temp.str();
371 }
372
373 DataArrayView::ValueType::size_type
374 DataTagged::getPointOffset(int sampleNo,
375 int dataPointNo) const
376 {
377 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
378 DataMapType::const_iterator pos(m_offsetLookup.find(tagKey));
379 DataArrayView::ValueType::size_type offset=m_defaultValueOffset;
380 if (pos!=m_offsetLookup.end()) {
381 offset=pos->second;
382 }
383 return offset;
384 }
385
386 DataArrayView
387 DataTagged::getDataPointByTag(int tag) const
388 {
389 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
390 DataArrayView::ValueType::size_type offset=m_defaultValueOffset;
391 if (pos!=m_offsetLookup.end()) {
392 offset=pos->second;
393 }
394 DataArrayView temp(getPointDataView());
395 temp.setOffset(offset);
396 return temp;
397 }
398
399 DataArrayView
400 DataTagged::getDataPoint(int sampleNo,
401 int dataPointNo)
402 {
403 EsysAssert(validSampleNo(sampleNo),"(getDataPoint) Invalid sampleNo: " << sampleNo);
404 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
405 return getDataPointByTag(tagKey);
406 }
407
408 int
409 DataTagged::archiveData(ofstream& archiveFile,
410 const DataArrayView::ValueType::size_type noValues) const
411 {
412 return(m_data.archiveData(archiveFile, noValues));
413 }
414
415 int
416 DataTagged::extractData(ifstream& archiveFile,
417 const DataArrayView::ValueType::size_type noValues)
418 {
419 return(m_data.extractData(archiveFile, noValues));
420 }
421 void
422 DataTagged::symmetric(DataAbstract* ev)
423 {
424 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
425 if (temp_ev==0) {
426 throw DataException("Error - DataTagged::symmetric casting to DataTagged failed (probably a programming error).");
427 }
428 const DataTagged::DataMapType& thisLookup=getTagLookup();
429 DataTagged::DataMapType::const_iterator i;
430 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
431 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
432 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
433 DataArrayView thisView=getDataPointByTag(i->first);
434 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
435 DataArrayView::symmetric(thisView,0,evView,0);
436 }
437 DataArrayView::symmetric(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
438 }
439 void
440 DataTagged::nonsymmetric(DataAbstract* ev)
441 {
442 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
443 if (temp_ev==0) {
444 throw DataException("Error - DataTagged::nonsymmetric casting to DataTagged failed (probably a programming error).");
445 }
446 const DataTagged::DataMapType& thisLookup=getTagLookup();
447 DataTagged::DataMapType::const_iterator i;
448 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
449 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
450 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
451 DataArrayView thisView=getDataPointByTag(i->first);
452 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
453 DataArrayView::nonsymmetric(thisView,0,evView,0);
454 }
455 DataArrayView::nonsymmetric(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
456 }
457 void
458 DataTagged::trace(DataAbstract* ev, int axis_offset)
459 {
460 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
461 if (temp_ev==0) {
462 throw DataException("Error - DataTagged::trace casting to DataTagged failed (probably a programming error).");
463 }
464 const DataTagged::DataMapType& thisLookup=getTagLookup();
465 DataTagged::DataMapType::const_iterator i;
466 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
467 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
468 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
469 DataArrayView thisView=getDataPointByTag(i->first);
470 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
471 DataArrayView::trace(thisView,0,evView,0, axis_offset);
472 }
473 DataArrayView::trace(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis_offset);
474 }
475
476 void
477 DataTagged::transpose(DataAbstract* ev, int axis_offset)
478 {
479 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
480 if (temp_ev==0) {
481 throw DataException("Error - DataTagged::transpose casting to DataTagged failed (probably a programming error).");
482 }
483 const DataTagged::DataMapType& thisLookup=getTagLookup();
484 DataTagged::DataMapType::const_iterator i;
485 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
486 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
487 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
488 DataArrayView thisView=getDataPointByTag(i->first);
489 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
490 DataArrayView::transpose(thisView,0,evView,0, axis_offset);
491 }
492 DataArrayView::transpose(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis_offset);
493 }
494
495 void
496 DataTagged::swapaxes(DataAbstract* ev, int axis0, int axis1)
497 {
498 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
499 if (temp_ev==0) {
500 throw DataException("Error - DataTagged::swapaxes casting to DataTagged failed (probably a programming error).");
501 }
502 const DataTagged::DataMapType& thisLookup=getTagLookup();
503 DataTagged::DataMapType::const_iterator i;
504 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
505 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
506 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
507 DataArrayView thisView=getDataPointByTag(i->first);
508 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
509 DataArrayView::swapaxes(thisView,0,evView,0,axis0,axis1);
510 }
511 DataArrayView::swapaxes(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis0,axis1);
512 }
513
514 void
515 DataTagged::eigenvalues(DataAbstract* ev)
516 {
517 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
518 if (temp_ev==0) {
519 throw DataException("Error - DataTagged::eigenvalues casting to DataTagged failed (propably a programming error).");
520 }
521 const DataTagged::DataMapType& thisLookup=getTagLookup();
522 DataTagged::DataMapType::const_iterator i;
523 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
524 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
525 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
526 DataArrayView thisView=getDataPointByTag(i->first);
527 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
528 DataArrayView::eigenvalues(thisView,0,evView,0);
529 }
530 DataArrayView::eigenvalues(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
531 }
532 void
533 DataTagged::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
534 {
535 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
536 if (temp_ev==0) {
537 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
538 }
539 DataTagged* temp_V=dynamic_cast<DataTagged*>(V);
540 if (temp_V==0) {
541 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
542 }
543 const DataTagged::DataMapType& thisLookup=getTagLookup();
544 DataTagged::DataMapType::const_iterator i;
545 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
546 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
547 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
548 temp_V->addTaggedValue(i->first,temp_V->getDefaultValue());
549 DataArrayView thisView=getDataPointByTag(i->first);
550 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
551 DataArrayView VView=temp_V->getDataPointByTag(i->first);
552 DataArrayView::eigenvalues_and_eigenvectors(thisView,0,evView,0,VView,0,tol);
553 }
554 DataArrayView::eigenvalues_and_eigenvectors(getDefaultValue(),0,
555 temp_ev->getDefaultValue(),0,
556 temp_V->getDefaultValue(),0,
557 tol);
558
559
560 }
561
562 void
563 DataTagged::setToZero(){
564 DataArrayView::ValueType::size_type n=m_data.size();
565 for (int i=0; i<n ;++i) m_data[i]=0.;
566 }
567
568 void
569 DataTagged::dump(const std::string fileName) const
570 {
571 #ifdef PASO_MPI
572 throw DataException("Error - DataTagged:: dump is not implemented for MPI yet.");
573 #endif
574 #ifdef USE_NETCDF
575 const int ldims=DataArrayView::maxRank+1;
576 const NcDim* ncdims[ldims];
577 NcVar *var, *tags_var;
578 int rank = getPointDataView().getRank();
579 int type= getFunctionSpace().getTypeCode();
580 int ndims =0;
581 long dims[ldims];
582 const double* d_ptr=&(m_data[0]);
583 DataArrayView::ShapeType shape = getPointDataView().getShape();
584
585 // netCDF error handler
586 NcError err(NcError::verbose_nonfatal);
587 // Create the file.
588 NcFile dataFile(fileName.c_str(), NcFile::Replace);
589 // check if writing was successful
590 if (!dataFile.is_valid())
591 throw DataException("Error - DataTagged:: opening of netCDF file for output failed.");
592 if (!dataFile.add_att("type_id",1) )
593 throw DataException("Error - DataTagged:: appending data type to netCDF file failed.");
594 if (!dataFile.add_att("rank",rank) )
595 throw DataException("Error - DataTagged:: appending rank attribute to netCDF file failed.");
596 if (!dataFile.add_att("function_space_type",type))
597 throw DataException("Error - DataTagged:: appending function space attribute to netCDF file failed.");
598 ndims=rank+1;
599 if ( rank >0 ) {
600 dims[0]=shape[0];
601 if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
602 throw DataException("Error - DataTagged:: appending ncdimsion 0 to netCDF file failed.");
603 }
604 if ( rank >1 ) {
605 dims[1]=shape[1];
606 if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
607 throw DataException("Error - DataTagged:: appending ncdimsion 1 to netCDF file failed.");
608 }
609 if ( rank >2 ) {
610 dims[2]=shape[2];
611 if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
612 throw DataException("Error - DataTagged:: appending ncdimsion 2 to netCDF file failed.");
613 }
614 if ( rank >3 ) {
615 dims[3]=shape[3];
616 if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
617 throw DataException("Error - DataTagged:: appending ncdimsion 3 to netCDF file failed.");
618 }
619 const DataTagged::DataMapType& thisLookup=getTagLookup();
620 DataTagged::DataMapType::const_iterator i;
621 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
622 int ntags=1;
623 for (i=thisLookup.begin();i!=thisLookupEnd;i++) ntags++;
624 int* tags =(int*) esysUtils::malloc(ntags*sizeof(int));
625 int c=1;
626 tags[0]=-1;
627 for (i=thisLookup.begin();i!=thisLookupEnd;i++) tags[c++]=i->first;
628 dims[rank]=ntags;
629 if (! (ncdims[rank] = dataFile.add_dim("num_tags", dims[rank])) )
630 {
631 esysUtils::free(tags);
632 throw DataException("Error - DataTagged:: appending num_tags to netCDF file failed.");
633 }
634 if (! ( tags_var = dataFile.add_var("tags", ncInt, ncdims[rank])) )
635 {
636 esysUtils::free(tags);
637 throw DataException("Error - DataTagged:: appending tags to netCDF file failed.");
638 }
639 if (! (tags_var->put(tags,dims[rank])) )
640 {
641 esysUtils::free(tags);
642 throw DataException("Error - DataTagged:: copy tags to netCDF buffer failed.");
643 }
644 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
645 {
646 esysUtils::free(tags);
647 throw DataException("Error - DataTagged:: appending variable to netCDF file failed.");
648 }
649 if (! (var->put(d_ptr,dims)) )
650 {
651 esysUtils::free(tags);
652 throw DataException("Error - DataTagged:: copy data to netCDF buffer failed.");
653 }
654 #else
655 throw DataException("Error - DataTagged:: dump is not configured with netCDF. Please contact your installation manager.");
656 #endif
657 }
658 } // end of namespace

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26