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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1802 - (show annotations)
Tue Sep 23 01:03:29 2008 UTC (11 years, 5 months ago) by jfenwick
File size: 34685 byte(s)
Added canTag methods to FunctionSpace and AbstractDomain (and its 
offspring).
This checks to see if the domain supports tags for the given type of 
function space.

Constructors for DataTagged now throw exceptions if you attempt to make 
a DataTagged with a FunctionSpace which does not support tags.
To allow the default constructor to work, NullDomain has a single 
functioncode which "supports" tagging.

Fixed a bug in DataTagged::toString and DataTypes::pointToString.

Added FunctionSpace::getListOfTagsSTL.

algorithm(DataTagged, BinaryFunction) in DataAlgorithm now only 
processes tags known to be in use.
This fixes mantis issue #0000186.

Added comment to Data.h intro warning about holding references if the 
underlying DataAbstract changes.

_python_ unit tests have been updated to test TaggedData with invalid 
FunctionSpaces and to give the correct answers to Lsup etc.


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 #include "DataMaths.h"
25
26 using namespace std;
27
28 namespace escript {
29
30 DataTagged::DataTagged()
31 : DataAbstract(FunctionSpace(),DataTypes::scalarShape)
32 {
33 // default constructor
34
35 // create a scalar default value
36 m_data.resize(1,0.,1);
37 /* DataArrayView temp(m_data,DataTypes::ShapeType());
38 setPointDataView(temp);*/
39 }
40
41 // DataTagged::DataTagged(const TagListType& tagKeys,
42 // const ValueListType& values,
43 // const DataArrayView& defaultValue,
44 // const FunctionSpace& what)
45 // : DataAbstract(what)
46 // {
47 // // constructor
48 //
49 // // initialise the array of data values
50 // // the default value is always the first item in the values list
51 // int len = defaultValue.noValues();
52 // m_data.resize(len,0.,len);
53 // for (int i=0; i<defaultValue.noValues(); i++) {
54 // m_data[i]=defaultValue.getData(i);
55 // }
56 //
57 // // create the data view
58 // DataArrayView temp(m_data,defaultValue.getShape());
59 // setPointDataView(temp);
60 //
61 // // add remaining tags and values
62 // addTaggedValues(tagKeys,values);
63 // }
64
65 DataTagged::DataTagged(const FunctionSpace& what,
66 const DataTypes::ShapeType &shape,
67 const int tags[],
68 const ValueType& data)
69 : DataAbstract(what,shape)
70 {
71 // alternative constructor
72 // not unit_tested tested yet
73 // It is not explicitly unit tested yet, but it is called from DataFactory
74
75 if (!what.canTag())
76 {
77 throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
78 }
79 // copy the data
80 m_data=data;
81
82 // create the view of the data
83 // DataArrayView tempView(m_data,shape);
84 // setPointDataView(tempView);
85
86 // we can't rely on the tag array to give us the number of tags so
87 // use the data we have been passed
88 int valsize=DataTypes::noValues(shape);
89 int ntags=data.size()/valsize;
90
91 // create the tag lookup map
92 // we assume that the first value and first tag are the default value so we skip
93 for (int i=1;i<ntags;++i)
94 {
95 m_offsetLookup.insert(DataMapType::value_type(tags[i],i*valsize));
96 }
97 }
98
99 DataTagged::DataTagged(const FunctionSpace& what,
100 const DataTypes::ShapeType &shape,
101 const TagListType& tags,
102 const ValueType& data)
103 : DataAbstract(what,shape)
104 {
105 // alternative constructor
106
107 if (!what.canTag())
108 {
109 throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
110 }
111
112 // copy the data
113 m_data=data;
114
115 // create the view of the data
116 // DataArrayView tempView(m_data,shape);
117 // setPointDataView(tempView);
118
119 // create the tag lookup map
120
121 // for (int sampleNo=0; sampleNo<getNumSamples(); sampleNo++) {
122 // m_offsetLookup.insert(DataMapType::value_type(sampleNo,tags[sampleNo]));
123 // }
124
125 // The above code looks like it will create a map the wrong way around
126
127 int valsize=DataTypes::noValues(shape);
128 int npoints=(data.size()/valsize)-1;
129 int ntags=tags.size();
130 if (ntags>npoints)
131 { // This throw is not unit tested yet
132 throw DataException("Programming error - Too many tags for the supplied values.");
133 }
134
135 // create the tag lookup map
136 // we assume that the first value is the default value so we skip it (hence the i+1 below)
137 for (int i=0;i<ntags;++i)
138 {
139 m_offsetLookup.insert(DataMapType::value_type(tags[i],(i+1)*valsize));
140 }
141 }
142
143
144 DataTagged::DataTagged(const DataTagged& other)
145 : DataAbstract(other.getFunctionSpace(),other.getShape()),
146 m_data(other.m_data),
147 m_offsetLookup(other.m_offsetLookup)
148 {
149 // copy constructor
150
151 // create the data view
152 // DataArrayView temp(m_data,other.getPointDataView().getShape());
153 // setPointDataView(temp);
154 }
155
156 DataTagged::DataTagged(const DataConstant& other)
157 : DataAbstract(other.getFunctionSpace(),other.getShape())
158 {
159 // copy constructor
160
161 if (!other.getFunctionSpace().canTag())
162 {
163 throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
164 }
165
166 // fill the default value with the constant value item from "other"
167 // const DataArrayView& value=other.getPointDataView();
168 int len = other.getNoValues();
169 m_data.resize(len,0.,len);
170 for (int i=0; i<len; i++) {
171 m_data[i]=other.getVector()[i];
172 }
173
174 // // create the data view
175 // DataArrayView temp(m_data,value.getShape());
176 // setPointDataView(temp);
177 }
178
179
180 // Create a new object by copying tags
181 DataTagged::DataTagged(const FunctionSpace& what,
182 const DataTypes::ShapeType& shape,
183 const DataTypes::ValueType& defaultvalue,
184 const DataTagged* tagsource)
185 : DataAbstract(what,shape)
186 {
187 // This constructor has not been unit tested yet
188
189 if (defaultvalue.size()!=DataTypes::noValues(shape)) {
190 throw DataException("Programming error - defaultvalue does not match supplied shape.");
191 }
192
193
194 if (!what.canTag())
195 {
196 throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
197 }
198
199 if (tagsource!=0)
200 {
201 int numtags=tagsource->getTagLookup().size();
202 // m_offsetLookup.reserve(tagsource.getTagLookup().size());
203 m_data.resize(defaultvalue.size(),0.); // since this is tagged data, we should have blocksize=1
204
205 DataTagged::DataMapType::const_iterator i;
206 for (i=tagsource->getTagLookup().begin();i!=tagsource->getTagLookup().end();i++) {
207 addTag(i->first);
208 }
209 }
210 else
211 {
212 m_data.resize(defaultvalue.size());
213 }
214
215
216
217 // need to set the default value ....
218 for (int i=0; i<defaultvalue.size(); i++) {
219 m_data[i]=defaultvalue[i];
220 }
221 }
222
223 DataAbstract*
224 DataTagged::deepCopy()
225 {
226 return new DataTagged(*this);
227 }
228
229 DataAbstract*
230 DataTagged::getSlice(const DataTypes::RegionType& region) const
231 {
232 return new DataTagged(*this, region);
233 }
234
235 DataTagged::DataTagged(const DataTagged& other,
236 const DataTypes::RegionType& region)
237 : DataAbstract(other.getFunctionSpace(),DataTypes::getResultSliceShape(region))
238 {
239 // slice constructor
240
241 // get the shape of the slice to copy from other
242 DataTypes::ShapeType regionShape(DataTypes::getResultSliceShape(region));
243 DataTypes::RegionLoopRangeType regionLoopRange=DataTypes::getSliceRegionLoopRange(region);
244
245 // allocate enough space in this for all values
246 // (need to add one to allow for the default value)
247 int len = DataTypes::noValues(regionShape)*(other.m_offsetLookup.size()+1);
248 m_data.resize(len,0.0,len);
249
250 // create the data view
251 // DataArrayView temp(m_data,regionShape);
252 // setPointDataView(temp);
253
254 // copy the default value from other to this
255 // getDefaultValue().copySlice(other.getDefaultValue(), regionLoopRange);
256 const DataTypes::ShapeType& otherShape=other.getShape();
257 const DataTypes::ValueType& otherData=other.getVector();
258 DataTypes::copySlice(getVector(),getShape(),getDefaultOffset(),otherData,otherShape,other.getDefaultOffset(), regionLoopRange);
259
260 // loop through the tag values copying these
261 DataMapType::const_iterator pos;
262 DataTypes::ValueType::size_type tagOffset=getNoValues();
263 for (pos=other.m_offsetLookup.begin();pos!=other.m_offsetLookup.end();pos++){
264 // getPointDataView().copySlice(tagOffset,other.getPointDataView(),pos->second,regionLoopRange);
265 DataTypes::copySlice(m_data,getShape(),tagOffset,otherData, otherShape, pos->second, regionLoopRange);
266 m_offsetLookup.insert(DataMapType::value_type(pos->first,tagOffset));
267 tagOffset+=getNoValues();
268 }
269 }
270
271 void
272 DataTagged::setSlice(const DataAbstract* other,
273 const DataTypes::RegionType& region)
274 {
275
276 // other must be another DataTagged object
277 // Data:setSlice implementation should ensure this
278 const DataTagged* otherTemp=dynamic_cast<const DataTagged*>(other);
279 if (otherTemp==0) {
280 throw DataException("Programming error - casting to DataTagged.");
281 }
282
283 // determine shape of the specified region
284 DataTypes::ShapeType regionShape(DataTypes::getResultSliceShape(region));
285
286 // modify region specification as needed to match rank of this object
287 DataTypes::RegionLoopRangeType regionLoopRange=DataTypes::getSliceRegionLoopRange(region);
288
289 // ensure rank/shape of this object is compatible with specified region
290 if (getRank()!=region.size()) {
291 throw DataException("Error - Invalid slice region.");
292 }
293 if (otherTemp->getRank()>0 && !DataTypes::checkShape(other->getShape(),regionShape)) {
294 throw DataException (DataTypes::createShapeErrorMessage(
295 "Error - Couldn't copy slice due to shape mismatch.",regionShape,other->getShape()));
296 }
297
298 const DataTypes::ValueType& otherData=otherTemp->getVector();
299 const DataTypes::ShapeType& otherShape=otherTemp->getShape();
300 // copy slice from other default value to this default value
301 // getDefaultValue().copySliceFrom(otherTemp->getDefaultValue(), regionLoopRange);
302 DataTypes::copySliceFrom(m_data,getShape(),getDefaultOffset(),otherData,otherShape,otherTemp->getDefaultOffset(),regionLoopRange);
303
304 // loop through tag values in other, adding any which aren't in this, using default value
305 DataMapType::const_iterator pos;
306 for (pos=otherTemp->m_offsetLookup.begin();pos!=otherTemp->m_offsetLookup.end();pos++) {
307 if (!isCurrentTag(pos->first)) {
308 addTag(pos->first);
309 }
310 }
311
312 // loop through the tag values copying slices from other to this
313 for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++) {
314 // getDataPointByTag(pos->first).copySliceFrom(otherTemp->getDataPointByTag(pos->first), regionLoopRange);
315 DataTypes::copySliceFrom(m_data,getShape(),getOffsetForTag(pos->first),otherData, otherShape, otherTemp->getOffsetForTag(pos->first), regionLoopRange);
316
317 }
318
319 }
320
321 int
322 DataTagged::getTagNumber(int dpno)
323 {
324 //
325 // Get the number of samples and data-points per sample
326 int numSamples = getNumSamples();
327 int numDataPointsPerSample = getNumDPPSample();
328 int numDataPoints = numSamples * numDataPointsPerSample;
329
330 if (numDataPointsPerSample==0) {
331 throw DataException("DataTagged::getTagNumber error: no data-points associated with this object.");
332 }
333
334 if (dpno<0 || dpno>numDataPoints-1) {
335 throw DataException("DataTagged::getTagNumber error: invalid data-point number supplied.");
336 }
337
338 //
339 // Determine the sample number which corresponds to this data-point number
340 int sampleNo = dpno / numDataPointsPerSample;
341
342 //
343 // Determine the tag number which corresponds to this sample number
344 int tagNo = getFunctionSpace().getTagFromSampleNo(sampleNo);
345
346 //
347 // return the tag number
348 return(tagNo);
349 }
350
351 // void
352 // DataTagged::setTaggedValues(const TagListType& tagKeys,
353 // const ValueListType& values)
354 // {
355 // addTaggedValues(tagKeys,values);
356 // }
357
358 void
359 DataTagged::setTaggedValue(int tagKey,
360 const DataTypes::ShapeType& pointshape,
361 const ValueType& value,
362 int dataOffset)
363 {
364 if (!DataTypes::checkShape(getShape(), pointshape)) {
365 throw DataException(DataTypes::createShapeErrorMessage(
366 "Error - Cannot setTaggedValue due to shape mismatch.", pointshape,getShape()));
367 }
368 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
369 if (pos==m_offsetLookup.end()) {
370 // tag couldn't be found so use addTaggedValue
371 addTaggedValue(tagKey,pointshape, value, dataOffset);
372 } else {
373 // copy the values into the data array at the offset determined by m_offsetLookup
374 int offset=pos->second;
375 for (int i=0; i<getNoValues(); i++) {
376 m_data[offset+i]=value[i+dataOffset];
377 }
378 }
379 }
380
381
382 /*
383 void
384 DataTagged::setTaggedValue(int tagKey,
385 const DataArrayView& value)
386 {
387 if (!getPointDataView().checkShape(value.getShape())) {
388 throw DataException(DataTypes::createShapeErrorMessage(
389 "Error - Cannot setTaggedValue due to shape mismatch.", value.getShape(),getShape()));
390 }
391 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
392 if (pos==m_offsetLookup.end()) {
393 // tag couldn't be found so use addTaggedValue
394 addTaggedValue(tagKey,value);
395 } else {
396 // copy the values into the data array at the offset determined by m_offsetLookup
397 int offset=pos->second;
398 for (int i=0; i<getPointDataView().noValues(); i++) {
399 m_data[offset+i]=value.getData(i);
400 }
401 }
402 }*/
403
404 // void
405 // DataTagged::addTaggedValues(const TagListType& tagKeys,
406 // const ValueListType& values)
407 // {
408 // if (values.size()==0) {
409 // // copy the current default value for each of the tags
410 // TagListType::const_iterator iT;
411 // for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
412 // // the point data view for DataTagged points at the default value
413 // addTaggedValue(*iT,getPointDataView());
414 // }
415 // } else if (values.size()==1 && tagKeys.size()>1) {
416 // // assume the one given value will be used for all tag values
417 // TagListType::const_iterator iT;
418 // for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
419 // addTaggedValue(*iT,values[0]);
420 // }
421 // } else {
422 // if (tagKeys.size()!=values.size()) {
423 // stringstream temp;
424 // temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
425 // << " doesn't match number of values: " << values.size();
426 // throw DataException(temp.str());
427 // } else {
428 // unsigned int i;
429 // for (i=0;i<tagKeys.size();i++) {
430 // addTaggedValue(tagKeys[i],values[i]);
431 // }
432 // }
433 // }
434 // }
435
436
437 void
438 DataTagged::addTaggedValues(const TagListType& tagKeys,
439 const ValueBatchType& values,
440 const ShapeType& vShape)
441 {
442 DataTypes::ValueType t(values.size(),0);
443 for (int i=0;i<values.size();++i)
444 {
445 t[i]=values[i];
446 }
447 addTaggedValues(tagKeys,t,vShape);
448 }
449
450
451 // Note: The check to see if vShape==our shape is done in the addTaggedValue method
452 void
453 DataTagged::addTaggedValues(const TagListType& tagKeys,
454 const ValueType& values,
455 const ShapeType& vShape)
456 {
457 int n=getNoValues();
458 int numVals=values.size()/n;
459 if (values.size()==0) {
460 // copy the current default value for each of the tags
461 TagListType::const_iterator iT;
462 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
463 // the point data view for DataTagged points at the default value
464 addTag(*iT);
465 }
466 } else if (numVals==1 && tagKeys.size()>1) {
467 // assume the one given value will be used for all tag values
468 TagListType::const_iterator iT;
469 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
470 addTaggedValue(*iT, vShape, values,0);
471 }
472 } else {
473 if (tagKeys.size()!=numVals) {
474 stringstream temp;
475 temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
476 << " doesn't match number of values: " << values.size();
477 throw DataException(temp.str());
478 } else {
479 unsigned int i;
480 int offset=0;
481 for (i=0;i<tagKeys.size();i++ ,offset+=n) {
482 addTaggedValue(tagKeys[i],vShape,values,offset);
483 }
484 }
485 }
486 }
487
488
489
490
491 void
492 DataTagged::addTaggedValue(int tagKey,
493 const DataTypes::ShapeType& pointshape,
494 const ValueType& value,
495 int dataOffset)
496 {
497 if (!DataTypes::checkShape(getShape(), pointshape)) {
498 throw DataException(DataTypes::createShapeErrorMessage(
499 "Error - Cannot addTaggedValue due to shape mismatch.", pointshape,getShape()));
500 }
501 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
502 if (pos!=m_offsetLookup.end()) {
503 // tag already exists so use setTaggedValue
504 setTaggedValue(tagKey,pointshape, value, dataOffset);
505 } else {
506 // save the key and the location of its data in the lookup tab
507 m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
508 // add the data given in "value" at the end of m_data
509 // need to make a temp copy of m_data, resize m_data, then copy
510 // all the old values plus the value to be added back into m_data
511 ValueType m_data_temp(m_data);
512 int oldSize=m_data.size();
513 int newSize=m_data.size()+getNoValues();
514 m_data.resize(newSize,0.,newSize);
515 for (int i=0;i<oldSize;i++) {
516 m_data[i]=m_data_temp[i];
517 }
518 for (int i=0;i<getNoValues();i++) {
519 m_data[oldSize+i]=value[i+dataOffset];
520 }
521 }
522 }
523
524
525
526
527 // void
528 // DataTagged::addTaggedValue(int tagKey,
529 // const DataArrayView& value)
530 // {
531 // if (!getPointDataView().checkShape(value.getShape())) {
532 // throw DataException(DataTypes::createShapeErrorMessage(
533 // "Error - Cannot addTaggedValue due to shape mismatch.", value.getShape(),getShape()));
534 // }
535 // DataMapType::iterator pos(m_offsetLookup.find(tagKey));
536 // if (pos!=m_offsetLookup.end()) {
537 // // tag already exists so use setTaggedValue
538 // setTaggedValue(tagKey,value);
539 // } else {
540 // // save the key and the location of its data in the lookup tab
541 // m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
542 // // add the data given in "value" at the end of m_data
543 // // need to make a temp copy of m_data, resize m_data, then copy
544 // // all the old values plus the value to be added back into m_data
545 // ValueType m_data_temp(m_data);
546 // int oldSize=m_data.size();
547 // int newSize=m_data.size()+value.noValues();
548 // m_data.resize(newSize,0.,newSize);
549 // for (int i=0;i<oldSize;i++) {
550 // m_data[i]=m_data_temp[i];
551 // }
552 // for (int i=0;i<value.noValues();i++) {
553 // m_data[oldSize+i]=value.getData(i);
554 // }
555 // }
556 // }
557
558
559 void
560 DataTagged::addTag(int tagKey)
561 {
562 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
563 if (pos!=m_offsetLookup.end()) {
564 // tag already exists so use setTaggedValue
565 // setTaggedValue(tagKey,value);
566 } else {
567 // save the key and the location of its data in the lookup tab
568 m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
569 // add the data given in "value" at the end of m_data
570 // need to make a temp copy of m_data, resize m_data, then copy
571 // all the old values plus the value to be added back into m_data
572 ValueType m_data_temp(m_data);
573 int oldSize=m_data.size();
574 int newSize=m_data.size()+getNoValues();
575 m_data.resize(newSize,0.,newSize);
576 for (int i=0;i<oldSize;i++) {
577 m_data[i]=m_data_temp[i];
578 }
579 for (int i=0;i<getNoValues();i++) {
580 m_data[oldSize+i]=m_data[m_defaultValueOffset+i];
581 }
582 }
583 }
584
585
586 double*
587 DataTagged::getSampleDataByTag(int tag)
588 {
589 DataMapType::iterator pos(m_offsetLookup.find(tag));
590 if (pos==m_offsetLookup.end()) {
591 // tag couldn't be found so return the default value
592 return &(m_data[0]);
593 } else {
594 // return the data-point corresponding to the given tag
595 return &(m_data[pos->second]);
596 }
597 }
598
599 string
600 DataTagged::toString() const
601 {
602 using namespace escript::DataTypes;
603 string empty="";
604 stringstream temp;
605 DataMapType::const_iterator i;
606 temp << "Tag(Default)" << endl;
607 temp << pointToString(m_data,getShape(),getDefaultOffset(),empty) << endl;
608 // create a temporary view as the offset will be changed
609 // DataArrayView tempView(getPointDataView().getData(), getPointDataView().getShape());
610 for (i=m_offsetLookup.begin();i!=m_offsetLookup.end();++i) {
611 temp << "Tag(" << i->first << ")" << endl;
612 temp << pointToString(m_data,getShape(),i->second,empty) << endl;
613 // tempView.setOffset(i->second);
614 // temp << tempView.toString() << endl;
615 }
616 return temp.str();
617 }
618
619 DataTypes::ValueType::size_type
620 DataTagged::getPointOffset(int sampleNo,
621 int dataPointNo) const
622 {
623 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
624 DataMapType::const_iterator pos(m_offsetLookup.find(tagKey));
625 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
626 if (pos!=m_offsetLookup.end()) {
627 offset=pos->second;
628 }
629 return offset;
630 }
631
632 // DataArrayView
633 // DataTagged::getDataPointByTag(int tag) const
634 // {
635 // DataMapType::const_iterator pos(m_offsetLookup.find(tag));
636 // DataTypes::ValueType::size_type offset=m_defaultValueOffset;
637 // if (pos!=m_offsetLookup.end()) {
638 // offset=pos->second;
639 // }
640 // DataArrayView temp(getPointDataView());
641 // temp.setOffset(offset);
642 // return temp;
643 // }
644 //
645
646
647 DataTypes::ValueType::size_type
648 DataTagged::getOffsetForTag(int tag) const
649 {
650 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
651 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
652 if (pos!=m_offsetLookup.end()) {
653 offset=pos->second;
654 }
655 return offset;
656 }
657
658 DataTypes::ValueType::const_reference
659 DataTagged::getDataByTag(int tag, DataTypes::ValueType::size_type i) const
660 {
661 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
662 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
663 if (pos!=m_offsetLookup.end()) {
664 offset=pos->second;
665 }
666 return m_data[offset+i];
667 /* DataArrayView temp(getPointDataView());
668 temp.setOffset(offset);
669 return temp.getData()[offset+i];*/
670 }
671
672
673 DataTypes::ValueType::reference
674 DataTagged::getDataByTag(int tag, DataTypes::ValueType::size_type i)
675 {
676 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
677 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
678 if (pos!=m_offsetLookup.end()) {
679 offset=pos->second;
680 }
681 return m_data[offset+i];
682 /* DataArrayView temp(getPointDataView());
683 temp.setOffset(offset);
684 return temp.getData()[offset+i];*/
685 }
686
687
688
689
690
691
692 // DataArrayView
693 // DataTagged::getDataPoint(int sampleNo,
694 // int dataPointNo)
695 // {
696 // EsysAssert(validSampleNo(sampleNo),"(getDataPoint) Invalid sampleNo: " << sampleNo);
697 // int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
698 // return getDataPointByTag(tagKey);
699 // }
700
701
702 void
703 DataTagged::symmetric(DataAbstract* ev)
704 {
705 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
706 if (temp_ev==0) {
707 throw DataException("Error - DataTagged::symmetric casting to DataTagged failed (probably a programming error).");
708 }
709 const DataTagged::DataMapType& thisLookup=getTagLookup();
710 DataTagged::DataMapType::const_iterator i;
711 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
712 ValueType& evVec=temp_ev->getVector();
713 const ShapeType& evShape=temp_ev->getShape();
714 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
715 temp_ev->addTag(i->first);
716 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
717 // DataArrayView thisView=getDataPointByTag(i->first);
718 // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
719 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
720
721 // DataArrayView::symmetric(thisView,0,evView,0);
722 DataMaths::symmetric(m_data,getShape(),offset,evVec, evShape, evoffset);
723 }
724 // symmetric(m_data,getShape(),getDefaultOffset(),
725 DataMaths::symmetric(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset());
726 }
727
728
729 void
730 DataTagged::nonsymmetric(DataAbstract* ev)
731 {
732 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
733 if (temp_ev==0) {
734 throw DataException("Error - DataTagged::nonsymmetric casting to DataTagged failed (probably a programming error).");
735 }
736 const DataTagged::DataMapType& thisLookup=getTagLookup();
737 DataTagged::DataMapType::const_iterator i;
738 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
739 ValueType& evVec=temp_ev->getVector();
740 const ShapeType& evShape=temp_ev->getShape();
741 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
742 temp_ev->addTag(i->first);
743 /* DataArrayView thisView=getDataPointByTag(i->first);
744 DataArrayView evView=temp_ev->getDataPointByTag(i->first);*/
745 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
746 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
747 DataMaths::nonsymmetric(m_data,getShape(),offset,evVec, evShape, evoffset);
748 }
749 DataMaths::nonsymmetric(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset());
750 }
751
752
753 void
754 DataTagged::trace(DataAbstract* ev, int axis_offset)
755 {
756 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
757 if (temp_ev==0) {
758 throw DataException("Error - DataTagged::trace casting to DataTagged failed (probably a programming error).");
759 }
760 const DataTagged::DataMapType& thisLookup=getTagLookup();
761 DataTagged::DataMapType::const_iterator i;
762 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
763 ValueType& evVec=temp_ev->getVector();
764 const ShapeType& evShape=temp_ev->getShape();
765 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
766 temp_ev->addTag(i->first);
767 // DataArrayView thisView=getDataPointByTag(i->first);
768 // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
769 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
770 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
771 DataMaths::trace(m_data,getShape(),offset,evVec, evShape, evoffset, axis_offset);
772 }
773 DataMaths::trace(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset(),axis_offset);
774 }
775
776 void
777 DataTagged::transpose(DataAbstract* ev, int axis_offset)
778 {
779 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
780 if (temp_ev==0) {
781 throw DataException("Error - DataTagged::transpose casting to DataTagged failed (probably a programming error).");
782 }
783 const DataTagged::DataMapType& thisLookup=getTagLookup();
784 DataTagged::DataMapType::const_iterator i;
785 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
786 ValueType& evVec=temp_ev->getVector();
787 const ShapeType& evShape=temp_ev->getShape();
788 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
789 temp_ev->addTag(i->first);
790 // DataArrayView thisView=getDataPointByTag(i->first);
791 // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
792 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
793 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
794 DataMaths::transpose(m_data,getShape(),offset,evVec, evShape, evoffset, axis_offset);
795 }
796 DataMaths::transpose(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset(),axis_offset);
797 }
798
799 void
800 DataTagged::swapaxes(DataAbstract* ev, int axis0, int axis1)
801 {
802 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
803 if (temp_ev==0) {
804 throw DataException("Error - DataTagged::swapaxes casting to DataTagged failed (probably a programming error).");
805 }
806 const DataTagged::DataMapType& thisLookup=getTagLookup();
807 DataTagged::DataMapType::const_iterator i;
808 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
809 ValueType& evVec=temp_ev->getVector();
810 const ShapeType& evShape=temp_ev->getShape();
811 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
812 temp_ev->addTag(i->first);
813 /* DataArrayView thisView=getDataPointByTag(i->first);
814 DataArrayView evView=temp_ev->getDataPointByTag(i->first);*/
815 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
816 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
817 DataMaths::swapaxes(m_data,getShape(),offset,evVec, evShape, evoffset,axis0,axis1);
818 }
819 DataMaths::swapaxes(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset(),axis0,axis1);
820 }
821
822 void
823 DataTagged::eigenvalues(DataAbstract* ev)
824 {
825 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
826 if (temp_ev==0) {
827 throw DataException("Error - DataTagged::eigenvalues casting to DataTagged failed (propably a programming error).");
828 }
829 const DataTagged::DataMapType& thisLookup=getTagLookup();
830 DataTagged::DataMapType::const_iterator i;
831 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
832 ValueType& evVec=temp_ev->getVector();
833 const ShapeType& evShape=temp_ev->getShape();
834 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
835 temp_ev->addTag(i->first);
836 // DataArrayView thisView=getDataPointByTag(i->first);
837 // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
838 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
839 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
840 DataMaths::eigenvalues(m_data,getShape(),offset,evVec, evShape, evoffset);
841 }
842 DataMaths::eigenvalues(m_data,getShape(),getDefaultOffset(),evVec, evShape, temp_ev->getDefaultOffset());
843 }
844 void
845 DataTagged::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
846 {
847 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
848 if (temp_ev==0) {
849 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
850 }
851 DataTagged* temp_V=dynamic_cast<DataTagged*>(V);
852 if (temp_V==0) {
853 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
854 }
855 const DataTagged::DataMapType& thisLookup=getTagLookup();
856 DataTagged::DataMapType::const_iterator i;
857 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
858 ValueType& evVec=temp_ev->getVector();
859 const ShapeType& evShape=temp_ev->getShape();
860 ValueType& VVec=temp_V->getVector();
861 const ShapeType& VShape=temp_V->getShape();
862 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
863 temp_ev->addTag(i->first);
864 temp_V->addTag(i->first);
865 /* DataArrayView thisView=getDataPointByTag(i->first);
866 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
867 DataArrayView VView=temp_V->getDataPointByTag(i->first);*/
868 DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
869 DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
870 DataTypes::ValueType::size_type Voffset=temp_V->getOffsetForTag(i->first);
871 /* DataArrayView::eigenvalues_and_eigenvectors(thisView,0,evView,0,VView,0,tol);*/
872 DataMaths::eigenvalues_and_eigenvectors(m_data,getShape(),offset,evVec, evShape, evoffset,VVec,VShape,Voffset,tol);
873
874 }
875 DataMaths::eigenvalues_and_eigenvectors(m_data,getShape(),getDefaultOffset(),evVec, evShape,
876 temp_ev->getDefaultOffset(),VVec,VShape,
877 temp_V->getDefaultOffset(), tol);
878
879
880 }
881
882 void
883 DataTagged::setToZero(){
884 DataTypes::ValueType::size_type n=m_data.size();
885 for (int i=0; i<n ;++i) m_data[i]=0.;
886 }
887
888 void
889 DataTagged::dump(const std::string fileName) const
890 {
891 #ifdef PASO_MPI
892 throw DataException("Error - DataTagged:: dump is not implemented for MPI yet.");
893 #endif
894 #ifdef USE_NETCDF
895 const int ldims=DataTypes::maxRank+1;
896 const NcDim* ncdims[ldims];
897 NcVar *var, *tags_var;
898 int rank = getRank();
899 int type= getFunctionSpace().getTypeCode();
900 int ndims =0;
901 long dims[ldims];
902 const double* d_ptr=&(m_data[0]);
903 DataTypes::ShapeType shape = getShape();
904
905 // netCDF error handler
906 NcError err(NcError::verbose_nonfatal);
907 // Create the file.
908 NcFile dataFile(fileName.c_str(), NcFile::Replace);
909 // check if writing was successful
910 if (!dataFile.is_valid())
911 throw DataException("Error - DataTagged:: opening of netCDF file for output failed.");
912 if (!dataFile.add_att("type_id",1) )
913 throw DataException("Error - DataTagged:: appending data type to netCDF file failed.");
914 if (!dataFile.add_att("rank",rank) )
915 throw DataException("Error - DataTagged:: appending rank attribute to netCDF file failed.");
916 if (!dataFile.add_att("function_space_type",type))
917 throw DataException("Error - DataTagged:: appending function space attribute to netCDF file failed.");
918 ndims=rank+1;
919 if ( rank >0 ) {
920 dims[0]=shape[0];
921 if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
922 throw DataException("Error - DataTagged:: appending ncdimension 0 to netCDF file failed.");
923 }
924 if ( rank >1 ) {
925 dims[1]=shape[1];
926 if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
927 throw DataException("Error - DataTagged:: appending ncdimension 1 to netCDF file failed.");
928 }
929 if ( rank >2 ) {
930 dims[2]=shape[2];
931 if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
932 throw DataException("Error - DataTagged:: appending ncdimension 2 to netCDF file failed.");
933 }
934 if ( rank >3 ) {
935 dims[3]=shape[3];
936 if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
937 throw DataException("Error - DataTagged:: appending ncdimension 3 to netCDF file failed.");
938 }
939 const DataTagged::DataMapType& thisLookup=getTagLookup();
940 DataTagged::DataMapType::const_iterator i;
941 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
942 int ntags=1;
943 for (i=thisLookup.begin();i!=thisLookupEnd;i++) ntags++;
944 int* tags =(int*) esysUtils::malloc(ntags*sizeof(int));
945 int c=1;
946 tags[0]=-1;
947 for (i=thisLookup.begin();i!=thisLookupEnd;i++) tags[c++]=i->first;
948 dims[rank]=ntags;
949 if (! (ncdims[rank] = dataFile.add_dim("num_tags", dims[rank])) )
950 {
951 esysUtils::free(tags);
952 throw DataException("Error - DataTagged:: appending num_tags to netCDF file failed.");
953 }
954 if (! ( tags_var = dataFile.add_var("tags", ncInt, ncdims[rank])) )
955 {
956 esysUtils::free(tags);
957 throw DataException("Error - DataTagged:: appending tags to netCDF file failed.");
958 }
959 if (! (tags_var->put(tags,dims[rank])) )
960 {
961 esysUtils::free(tags);
962 throw DataException("Error - DataTagged:: copy tags to netCDF buffer failed.");
963 }
964 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
965 {
966 esysUtils::free(tags);
967 throw DataException("Error - DataTagged:: appending variable to netCDF file failed.");
968 }
969 if (! (var->put(d_ptr,dims)) )
970 {
971 esysUtils::free(tags);
972 throw DataException("Error - DataTagged:: copy data to netCDF buffer failed.");
973 }
974 #else
975 throw DataException("Error - DataTagged:: dump is not configured with netCDF. Please contact your installation manager.");
976 #endif
977 }
978
979 DataTypes::ValueType&
980 DataTagged::getVector()
981 {
982 return m_data;
983 }
984
985 const DataTypes::ValueType&
986 DataTagged::getVector() const
987 {
988 return m_data;
989 }
990
991 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26