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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1946 - (show annotations)
Wed Oct 29 05:48:53 2008 UTC (10 years, 6 months ago) by jfenwick
Original Path: trunk/escript/src/DataTagged.cpp
File size: 35005 byte(s)
A cleanup of some of the problems I found doing a Wall compile.

Removed some commented out lines.
Swapped some member initialisers.
Removed virtual qualifiers from some methods in FunctionSpace.
Fixed some unused or (possibly) uninitialised variables.


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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26