/[escript]/branches/arrayview_from_1695_trunk/escript/src/DataTagged.cpp
ViewVC logotype

Contents of /branches/arrayview_from_1695_trunk/escript/src/DataTagged.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1714 - (show annotations)
Thu Aug 21 00:01:55 2008 UTC (12 years, 1 month ago) by jfenwick
File size: 25654 byte(s)
Branch commit

Moved getSliceRegion() and getSliceRange() into DataTypes
Data.cpp - modified not to rely on operator() from DataArrayView
         - Used more const& to avoid copies


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,DataTypes::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 DataTypes::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 DataTypes::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
140 // Create a new object by copying tags
141 DataTagged::DataTagged(const FunctionSpace& what,
142 const DataTypes::ShapeType& shape,
143 const DataTypes::ValueType& defaultvalue,
144 const DataTagged& tagsource)
145 : DataAbstract(what)
146 {
147 // This constructor has not been unit tested yet
148
149 if (defaultvalue.size()!=DataTypes::noValues(shape)) {
150 throw DataException("Programming error - defaultvalue does not match supplied shape.");
151 }
152
153 setShape(shape);
154
155 int numtags=tagsource.getTagLookup().size();
156 // m_offsetLookup.reserve(tagsource.getTagLookup().size());
157 m_data.resize(defaultvalue.size(),0.); // since this is tagged data, we should have blocksize=1
158
159 // need to set the default value ....
160 for (int i=0; i<defaultvalue.size(); i++) {
161 m_data[i]=defaultvalue[i];
162 }
163
164
165 // create the data view
166 DataArrayView temp(m_data,shape);
167 setPointDataView(temp);
168
169 DataTagged::DataMapType::const_iterator i;
170 for (i=tagsource.getTagLookup().begin();i!=tagsource.getTagLookup().end();i++) {
171 addTag(i->first);
172 }
173 }
174
175
176 DataAbstract*
177 DataTagged::getSlice(const DataTypes::RegionType& region) const
178 {
179 return new DataTagged(*this, region);
180 }
181
182 DataTagged::DataTagged(const DataTagged& other,
183 const DataTypes::RegionType& region)
184 : DataAbstract(other.getFunctionSpace())
185 {
186 // slice constructor
187
188 // get the shape of the slice to copy from other
189 DataTypes::ShapeType regionShape(DataTypes::getResultSliceShape(region));
190 DataTypes::RegionLoopRangeType regionLoopRange=getSliceRegionLoopRange(region);
191
192 // allocate enough space in this for all values
193 // (need to add one to allow for the default value)
194 int len = DataTypes::noValues(regionShape)*(other.m_offsetLookup.size()+1);
195 m_data.resize(len,0.0,len);
196
197 // create the data view
198 DataArrayView temp(m_data,regionShape);
199 setPointDataView(temp);
200
201 // copy the default value from other to this
202 getDefaultValue().copySlice(other.getDefaultValue(), regionLoopRange);
203
204 // loop through the tag values copying these
205 DataMapType::const_iterator pos;
206 DataTypes::ValueType::size_type tagOffset=getPointDataView().noValues();
207 for (pos=other.m_offsetLookup.begin();pos!=other.m_offsetLookup.end();pos++){
208 getPointDataView().copySlice(tagOffset,other.getPointDataView(),pos->second,regionLoopRange);
209 m_offsetLookup.insert(DataMapType::value_type(pos->first,tagOffset));
210 tagOffset+=getPointDataView().noValues();
211 }
212 }
213
214 void
215 DataTagged::setSlice(const DataAbstract* other,
216 const DataTypes::RegionType& region)
217 {
218
219 // other must be another DataTagged object
220 // Data:setSlice implementation should ensure this
221 const DataTagged* otherTemp=dynamic_cast<const DataTagged*>(other);
222 if (otherTemp==0) {
223 throw DataException("Programming error - casting to DataTagged.");
224 }
225
226 // determine shape of the specified region
227 DataTypes::ShapeType regionShape(DataTypes::getResultSliceShape(region));
228
229 // modify region specification as needed to match rank of this object
230 DataTypes::RegionLoopRangeType regionLoopRange=getSliceRegionLoopRange(region);
231
232 // ensure rank/shape of this object is compatible with specified region
233 if (getPointDataView().getRank()!=region.size()) {
234 throw DataException("Error - Invalid slice region.");
235 }
236 if (otherTemp->getPointDataView().getRank()>0 && !other->getPointDataView().checkShape(regionShape)) {
237 throw DataException (other->getPointDataView().createShapeErrorMessage(
238 "Error - Couldn't copy slice due to shape mismatch.",regionShape));
239 }
240
241 // copy slice from other default value to this default value
242 getDefaultValue().copySliceFrom(otherTemp->getDefaultValue(), regionLoopRange);
243
244 // loop through tag values in other, adding any which aren't in this, using default value
245 DataMapType::const_iterator pos;
246 for (pos=otherTemp->m_offsetLookup.begin();pos!=otherTemp->m_offsetLookup.end();pos++) {
247 if (!isCurrentTag(pos->first)) {
248 addTag(pos->first);
249 }
250 }
251
252 // loop through the tag values copying slices from other to this
253 for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++) {
254 getDataPointByTag(pos->first).copySliceFrom(otherTemp->getDataPointByTag(pos->first), regionLoopRange);
255 }
256
257 }
258
259 int
260 DataTagged::getTagNumber(int dpno)
261 {
262 //
263 // Get the number of samples and data-points per sample
264 int numSamples = getNumSamples();
265 int numDataPointsPerSample = getNumDPPSample();
266 int numDataPoints = numSamples * numDataPointsPerSample;
267
268 if (numDataPointsPerSample==0) {
269 throw DataException("DataTagged::getTagNumber error: no data-points associated with this object.");
270 }
271
272 if (dpno<0 || dpno>numDataPoints-1) {
273 throw DataException("DataTagged::getTagNumber error: invalid data-point number supplied.");
274 }
275
276 //
277 // Determine the sample number which corresponds to this data-point number
278 int sampleNo = dpno / numDataPointsPerSample;
279
280 //
281 // Determine the tag number which corresponds to this sample number
282 int tagNo = getFunctionSpace().getTagFromSampleNo(sampleNo);
283
284 //
285 // return the tag number
286 return(tagNo);
287 }
288
289 void
290 DataTagged::setTaggedValues(const TagListType& tagKeys,
291 const ValueListType& values)
292 {
293 addTaggedValues(tagKeys,values);
294 }
295
296 void
297 DataTagged::setTaggedValue(int tagKey,
298 const DataArrayView& value)
299 {
300 if (!getPointDataView().checkShape(value.getShape())) {
301 throw DataException(getPointDataView().createShapeErrorMessage(
302 "Error - Cannot setTaggedValue due to shape mismatch.", value.getShape()));
303 }
304 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
305 if (pos==m_offsetLookup.end()) {
306 // tag couldn't be found so use addTaggedValue
307 addTaggedValue(tagKey,value);
308 } else {
309 // copy the values into the data array at the offset determined by m_offsetLookup
310 int offset=pos->second;
311 for (int i=0; i<getPointDataView().noValues(); i++) {
312 m_data[offset+i]=value.getData(i);
313 }
314 }
315 }
316
317 void
318 DataTagged::addTaggedValues(const TagListType& tagKeys,
319 const ValueListType& values)
320 {
321 if (values.size()==0) {
322 // copy the current default value for each of the tags
323 TagListType::const_iterator iT;
324 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
325 // the point data view for DataTagged points at the default value
326 addTaggedValue(*iT,getPointDataView());
327 }
328 } else if (values.size()==1 && tagKeys.size()>1) {
329 // assume the one given value will be used for all tag values
330 TagListType::const_iterator iT;
331 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
332 addTaggedValue(*iT,values[0]);
333 }
334 } else {
335 if (tagKeys.size()!=values.size()) {
336 stringstream temp;
337 temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
338 << " doesn't match number of values: " << values.size();
339 throw DataException(temp.str());
340 } else {
341 unsigned int i;
342 for (i=0;i<tagKeys.size();i++) {
343 addTaggedValue(tagKeys[i],values[i]);
344 }
345 }
346 }
347 }
348
349 void
350 DataTagged::addTaggedValue(int tagKey,
351 const DataArrayView& value)
352 {
353 if (!getPointDataView().checkShape(value.getShape())) {
354 throw DataException(getPointDataView().createShapeErrorMessage(
355 "Error - Cannot addTaggedValue due to shape mismatch.", value.getShape()));
356 }
357 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
358 if (pos!=m_offsetLookup.end()) {
359 // tag already exists so use setTaggedValue
360 setTaggedValue(tagKey,value);
361 } else {
362 // save the key and the location of its data in the lookup tab
363 m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
364 // add the data given in "value" at the end of m_data
365 // need to make a temp copy of m_data, resize m_data, then copy
366 // all the old values plus the value to be added back into m_data
367 ValueType m_data_temp(m_data);
368 int oldSize=m_data.size();
369 int newSize=m_data.size()+value.noValues();
370 m_data.resize(newSize,0.,newSize);
371 for (int i=0;i<oldSize;i++) {
372 m_data[i]=m_data_temp[i];
373 }
374 for (int i=0;i<value.noValues();i++) {
375 m_data[oldSize+i]=value.getData(i);
376 }
377 }
378 }
379
380
381 void
382 DataTagged::addTag(int tagKey)
383 {
384 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
385 if (pos!=m_offsetLookup.end()) {
386 // tag already exists so use setTaggedValue
387 // setTaggedValue(tagKey,value);
388 } else {
389 // save the key and the location of its data in the lookup tab
390 m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
391 // add the data given in "value" at the end of m_data
392 // need to make a temp copy of m_data, resize m_data, then copy
393 // all the old values plus the value to be added back into m_data
394 ValueType m_data_temp(m_data);
395 int oldSize=m_data.size();
396 int newSize=m_data.size()+getNoValues();
397 m_data.resize(newSize,0.,newSize);
398 for (int i=0;i<oldSize;i++) {
399 m_data[i]=m_data_temp[i];
400 }
401 for (int i=0;i<getNoValues();i++) {
402 m_data[oldSize+i]=m_data[m_defaultValueOffset+i];
403 }
404 }
405 }
406
407
408 double*
409 DataTagged::getSampleDataByTag(int tag)
410 {
411 DataMapType::iterator pos(m_offsetLookup.find(tag));
412 if (pos==m_offsetLookup.end()) {
413 // tag couldn't be found so return the default value
414 return &(m_data[0]);
415 } else {
416 // return the data-point corresponding to the given tag
417 return &(m_data[pos->second]);
418 }
419 }
420
421 string
422 DataTagged::toString() const
423 {
424 stringstream temp;
425 DataMapType::const_iterator i;
426 temp << "Tag(Default)" << endl;
427 temp << getDefaultValue().toString() << endl;
428 // create a temporary view as the offset will be changed
429 DataArrayView tempView(getPointDataView().getData(), getPointDataView().getShape());
430 for (i=m_offsetLookup.begin();i!=m_offsetLookup.end();++i) {
431 temp << "Tag(" << i->first << ")" << endl;
432 tempView.setOffset(i->second);
433 temp << tempView.toString() << endl;
434 }
435 return temp.str();
436 }
437
438 DataTypes::ValueType::size_type
439 DataTagged::getPointOffset(int sampleNo,
440 int dataPointNo) const
441 {
442 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
443 DataMapType::const_iterator pos(m_offsetLookup.find(tagKey));
444 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
445 if (pos!=m_offsetLookup.end()) {
446 offset=pos->second;
447 }
448 return offset;
449 }
450
451 DataArrayView
452 DataTagged::getDataPointByTag(int tag) const
453 {
454 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
455 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
456 if (pos!=m_offsetLookup.end()) {
457 offset=pos->second;
458 }
459 DataArrayView temp(getPointDataView());
460 temp.setOffset(offset);
461 return temp;
462 }
463
464
465 DataTypes::ValueType::const_reference
466 DataTagged::getDataByTag(int tag, DataTypes::ValueType::size_type i) const
467 {
468 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
469 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
470 if (pos!=m_offsetLookup.end()) {
471 offset=pos->second;
472 }
473 DataArrayView temp(getPointDataView());
474 temp.setOffset(offset);
475 return temp.getData()[offset+i];
476 }
477
478
479 DataTypes::ValueType::reference
480 DataTagged::getDataByTag(int tag, DataTypes::ValueType::size_type i)
481 {
482 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
483 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
484 if (pos!=m_offsetLookup.end()) {
485 offset=pos->second;
486 }
487 DataArrayView temp(getPointDataView());
488 temp.setOffset(offset);
489 return temp.getData()[offset+i];
490 }
491
492
493
494
495
496
497 DataArrayView
498 DataTagged::getDataPoint(int sampleNo,
499 int dataPointNo)
500 {
501 EsysAssert(validSampleNo(sampleNo),"(getDataPoint) Invalid sampleNo: " << sampleNo);
502 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
503 return getDataPointByTag(tagKey);
504 }
505
506 int
507 DataTagged::archiveData(ofstream& archiveFile,
508 const DataTypes::ValueType::size_type noValues) const
509 {
510 return(m_data.archiveData(archiveFile, noValues));
511 }
512
513 int
514 DataTagged::extractData(ifstream& archiveFile,
515 const DataTypes::ValueType::size_type noValues)
516 {
517 return(m_data.extractData(archiveFile, noValues));
518 }
519 void
520 DataTagged::symmetric(DataAbstract* ev)
521 {
522 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
523 if (temp_ev==0) {
524 throw DataException("Error - DataTagged::symmetric casting to DataTagged failed (probably a programming error).");
525 }
526 const DataTagged::DataMapType& thisLookup=getTagLookup();
527 DataTagged::DataMapType::const_iterator i;
528 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
529 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
530 temp_ev->addTag(i->first);
531 DataArrayView thisView=getDataPointByTag(i->first);
532 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
533 DataArrayView::symmetric(thisView,0,evView,0);
534 }
535 DataArrayView::symmetric(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
536 }
537 void
538 DataTagged::nonsymmetric(DataAbstract* ev)
539 {
540 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
541 if (temp_ev==0) {
542 throw DataException("Error - DataTagged::nonsymmetric casting to DataTagged failed (probably a programming error).");
543 }
544 const DataTagged::DataMapType& thisLookup=getTagLookup();
545 DataTagged::DataMapType::const_iterator i;
546 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
547 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
548 temp_ev->addTag(i->first);
549 DataArrayView thisView=getDataPointByTag(i->first);
550 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
551 DataArrayView::nonsymmetric(thisView,0,evView,0);
552 }
553 DataArrayView::nonsymmetric(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
554 }
555 void
556 DataTagged::trace(DataAbstract* ev, int axis_offset)
557 {
558 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
559 if (temp_ev==0) {
560 throw DataException("Error - DataTagged::trace casting to DataTagged failed (probably a programming error).");
561 }
562 const DataTagged::DataMapType& thisLookup=getTagLookup();
563 DataTagged::DataMapType::const_iterator i;
564 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
565 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
566 temp_ev->addTag(i->first);
567 DataArrayView thisView=getDataPointByTag(i->first);
568 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
569 DataArrayView::trace(thisView,0,evView,0, axis_offset);
570 }
571 DataArrayView::trace(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis_offset);
572 }
573
574 void
575 DataTagged::transpose(DataAbstract* ev, int axis_offset)
576 {
577 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
578 if (temp_ev==0) {
579 throw DataException("Error - DataTagged::transpose casting to DataTagged failed (probably a programming error).");
580 }
581 const DataTagged::DataMapType& thisLookup=getTagLookup();
582 DataTagged::DataMapType::const_iterator i;
583 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
584 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
585 temp_ev->addTag(i->first);
586 DataArrayView thisView=getDataPointByTag(i->first);
587 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
588 DataArrayView::transpose(thisView,0,evView,0, axis_offset);
589 }
590 DataArrayView::transpose(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis_offset);
591 }
592
593 void
594 DataTagged::swapaxes(DataAbstract* ev, int axis0, int axis1)
595 {
596 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
597 if (temp_ev==0) {
598 throw DataException("Error - DataTagged::swapaxes casting to DataTagged failed (probably a programming error).");
599 }
600 const DataTagged::DataMapType& thisLookup=getTagLookup();
601 DataTagged::DataMapType::const_iterator i;
602 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
603 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
604 temp_ev->addTag(i->first);
605 DataArrayView thisView=getDataPointByTag(i->first);
606 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
607 DataArrayView::swapaxes(thisView,0,evView,0,axis0,axis1);
608 }
609 DataArrayView::swapaxes(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis0,axis1);
610 }
611
612 void
613 DataTagged::eigenvalues(DataAbstract* ev)
614 {
615 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
616 if (temp_ev==0) {
617 throw DataException("Error - DataTagged::eigenvalues casting to DataTagged failed (propably a programming error).");
618 }
619 const DataTagged::DataMapType& thisLookup=getTagLookup();
620 DataTagged::DataMapType::const_iterator i;
621 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
622 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
623 temp_ev->addTag(i->first);
624 DataArrayView thisView=getDataPointByTag(i->first);
625 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
626 DataArrayView::eigenvalues(thisView,0,evView,0);
627 }
628 DataArrayView::eigenvalues(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
629 }
630 void
631 DataTagged::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
632 {
633 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
634 if (temp_ev==0) {
635 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
636 }
637 DataTagged* temp_V=dynamic_cast<DataTagged*>(V);
638 if (temp_V==0) {
639 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
640 }
641 const DataTagged::DataMapType& thisLookup=getTagLookup();
642 DataTagged::DataMapType::const_iterator i;
643 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
644 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
645 temp_ev->addTag(i->first);
646 temp_V->addTag(i->first);
647 DataArrayView thisView=getDataPointByTag(i->first);
648 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
649 DataArrayView VView=temp_V->getDataPointByTag(i->first);
650 DataArrayView::eigenvalues_and_eigenvectors(thisView,0,evView,0,VView,0,tol);
651 }
652 DataArrayView::eigenvalues_and_eigenvectors(getDefaultValue(),0,
653 temp_ev->getDefaultValue(),0,
654 temp_V->getDefaultValue(),0,
655 tol);
656
657
658 }
659
660 void
661 DataTagged::setToZero(){
662 DataTypes::ValueType::size_type n=m_data.size();
663 for (int i=0; i<n ;++i) m_data[i]=0.;
664 }
665
666 void
667 DataTagged::dump(const std::string fileName) const
668 {
669 #ifdef PASO_MPI
670 throw DataException("Error - DataTagged:: dump is not implemented for MPI yet.");
671 #endif
672 #ifdef USE_NETCDF
673 const int ldims=DataTypes::maxRank+1;
674 const NcDim* ncdims[ldims];
675 NcVar *var, *tags_var;
676 int rank = getPointDataView().getRank();
677 int type= getFunctionSpace().getTypeCode();
678 int ndims =0;
679 long dims[ldims];
680 const double* d_ptr=&(m_data[0]);
681 DataTypes::ShapeType shape = getPointDataView().getShape();
682
683 // netCDF error handler
684 NcError err(NcError::verbose_nonfatal);
685 // Create the file.
686 NcFile dataFile(fileName.c_str(), NcFile::Replace);
687 // check if writing was successful
688 if (!dataFile.is_valid())
689 throw DataException("Error - DataTagged:: opening of netCDF file for output failed.");
690 if (!dataFile.add_att("type_id",1) )
691 throw DataException("Error - DataTagged:: appending data type to netCDF file failed.");
692 if (!dataFile.add_att("rank",rank) )
693 throw DataException("Error - DataTagged:: appending rank attribute to netCDF file failed.");
694 if (!dataFile.add_att("function_space_type",type))
695 throw DataException("Error - DataTagged:: appending function space attribute to netCDF file failed.");
696 ndims=rank+1;
697 if ( rank >0 ) {
698 dims[0]=shape[0];
699 if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
700 throw DataException("Error - DataTagged:: appending ncdimsion 0 to netCDF file failed.");
701 }
702 if ( rank >1 ) {
703 dims[1]=shape[1];
704 if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
705 throw DataException("Error - DataTagged:: appending ncdimsion 1 to netCDF file failed.");
706 }
707 if ( rank >2 ) {
708 dims[2]=shape[2];
709 if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
710 throw DataException("Error - DataTagged:: appending ncdimsion 2 to netCDF file failed.");
711 }
712 if ( rank >3 ) {
713 dims[3]=shape[3];
714 if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
715 throw DataException("Error - DataTagged:: appending ncdimsion 3 to netCDF file failed.");
716 }
717 const DataTagged::DataMapType& thisLookup=getTagLookup();
718 DataTagged::DataMapType::const_iterator i;
719 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
720 int ntags=1;
721 for (i=thisLookup.begin();i!=thisLookupEnd;i++) ntags++;
722 int* tags =(int*) esysUtils::malloc(ntags*sizeof(int));
723 int c=1;
724 tags[0]=-1;
725 for (i=thisLookup.begin();i!=thisLookupEnd;i++) tags[c++]=i->first;
726 dims[rank]=ntags;
727 if (! (ncdims[rank] = dataFile.add_dim("num_tags", dims[rank])) )
728 {
729 esysUtils::free(tags);
730 throw DataException("Error - DataTagged:: appending num_tags to netCDF file failed.");
731 }
732 if (! ( tags_var = dataFile.add_var("tags", ncInt, ncdims[rank])) )
733 {
734 esysUtils::free(tags);
735 throw DataException("Error - DataTagged:: appending tags to netCDF file failed.");
736 }
737 if (! (tags_var->put(tags,dims[rank])) )
738 {
739 esysUtils::free(tags);
740 throw DataException("Error - DataTagged:: copy tags to netCDF buffer failed.");
741 }
742 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
743 {
744 esysUtils::free(tags);
745 throw DataException("Error - DataTagged:: appending variable to netCDF file failed.");
746 }
747 if (! (var->put(d_ptr,dims)) )
748 {
749 esysUtils::free(tags);
750 throw DataException("Error - DataTagged:: copy data to netCDF buffer failed.");
751 }
752 #else
753 throw DataException("Error - DataTagged:: dump is not configured with netCDF. Please contact your installation manager.");
754 #endif
755 }
756 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26