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