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 |