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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 580 - (show annotations)
Wed Mar 8 05:45:51 2006 UTC (13 years, 8 months ago) by gross
File size: 16171 byte(s)
faster version of the local eigenvalue calculation
1 // $Id$
2
3 /*
4 ******************************************************************************
5 * *
6 * COPYRIGHT ACcESS 2004 - All Rights Reserved *
7 * *
8 * This software is the property of ACcESS. No part of this code *
9 * may be copied in any form or by any means without the expressed written *
10 * consent of ACcESS. Copying, use or modification of this software *
11 * by any unauthorised person is illegal unless that person has a software *
12 * license agreement with ACcESS. *
13 * *
14 ******************************************************************************
15 */
16
17 #include "DataTagged.h"
18
19 #include "DataConstant.h"
20 #include "DataException.h"
21
22 using namespace std;
23
24 namespace escript {
25
26 DataTagged::DataTagged()
27 : DataAbstract(FunctionSpace())
28 {
29 // default constructor
30
31 // create a scalar default value
32 m_data.resize(1,0.,1);
33 DataArrayView temp(m_data,DataArrayView::ShapeType());
34 setPointDataView(temp);
35 }
36
37 DataTagged::DataTagged(const TagListType& tagKeys,
38 const ValueListType& values,
39 const DataArrayView& defaultValue,
40 const FunctionSpace& what)
41 : DataAbstract(what)
42 {
43 // constructor
44
45 // initialise the array of data values
46 // the default value is always the first item in the values list
47 int len = defaultValue.noValues();
48 m_data.resize(len,0.,len);
49 for (int i=0; i<defaultValue.noValues(); i++) {
50 m_data[i]=defaultValue.getData(i);
51 }
52
53 // create the data view
54 DataArrayView temp(m_data,defaultValue.getShape());
55 setPointDataView(temp);
56
57 // add remaining tags and values
58 addTaggedValues(tagKeys,values);
59 }
60
61 DataTagged::DataTagged(const FunctionSpace& what,
62 const DataArrayView::ShapeType &shape,
63 const int tags[],
64 const ValueType& data)
65 : DataAbstract(what)
66 {
67 // alternative constructor
68 // not unit_tested tested yet
69
70 // copy the data
71 m_data=data;
72
73 // create the view of the data
74 DataArrayView tempView(m_data,shape);
75 setPointDataView(tempView);
76
77 // create the tag lookup map
78 for (int sampleNo=0; sampleNo<getNumSamples(); sampleNo++) {
79 m_offsetLookup.insert(DataMapType::value_type(sampleNo,tags[sampleNo]));
80 }
81 }
82
83 DataTagged::DataTagged(const DataTagged& other)
84 : DataAbstract(other.getFunctionSpace()),
85 m_data(other.m_data),
86 m_offsetLookup(other.m_offsetLookup)
87 {
88 // copy constructor
89
90 // create the data view
91 DataArrayView temp(m_data,other.getPointDataView().getShape());
92 setPointDataView(temp);
93 }
94
95 DataTagged::DataTagged(const DataConstant& other)
96 : DataAbstract(other.getFunctionSpace())
97 {
98 // copy constructor
99
100 // fill the default value with the constant value item from "other"
101 const DataArrayView& value=other.getPointDataView();
102 int len = value.noValues();
103 m_data.resize(len,0.,len);
104 for (int i=0; i<value.noValues(); i++) {
105 m_data[i]=value.getData(i);
106 }
107
108 // create the data view
109 DataArrayView temp(m_data,value.getShape());
110 setPointDataView(temp);
111 }
112
113 void
114 DataTagged::reshapeDataPoint(const DataArrayView::ShapeType& shape)
115 {
116 // can only reshape a rank zero data point
117 if (getPointDataView().getRank()!=0) {
118 stringstream temp;
119 temp << "Error - Can only reshape Data with data points of rank 0. "
120 << "This Data has data points with rank: "
121 << getPointDataView().getRank();
122 throw DataException(temp.str());
123 }
124
125 // allocate enough space for all values
126 DataArrayView::ValueType newData(DataArrayView::noValues(shape)*(m_offsetLookup.size()+1));
127 DataArrayView newView(newData,shape);
128 newView.copy(0,getDefaultValue()());
129
130 // loop through the tag values
131 DataMapType::iterator pos;
132 DataArrayView::ValueType::size_type tagOffset=DataArrayView::noValues(shape);
133 for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++){
134 newView.copy(tagOffset,m_data[pos->second]);
135 pos->second=tagOffset;
136 tagOffset+=DataArrayView::noValues(shape);
137 }
138 m_data=newData;
139 DataArrayView temp(m_data,shape);
140 setPointDataView(temp);
141 }
142
143 DataAbstract*
144 DataTagged::getSlice(const DataArrayView::RegionType& region) const
145 {
146 return new DataTagged(*this, region);
147 }
148
149 DataTagged::DataTagged(const DataTagged& other,
150 const DataArrayView::RegionType& region)
151 : DataAbstract(other.getFunctionSpace())
152 {
153 // slice constructor
154
155 // get the shape of the slice to copy from other
156 DataArrayView::ShapeType regionShape(DataArrayView::getResultSliceShape(region));
157 DataArrayView::RegionLoopRangeType regionLoopRange=getSliceRegionLoopRange(region);
158
159 // allocate enough space in this for all values
160 // (need to add one to allow for the default value)
161 int len = DataArrayView::noValues(regionShape)*(other.m_offsetLookup.size()+1);
162 m_data.resize(len,0.0,len);
163
164 // create the data view
165 DataArrayView temp(m_data,regionShape);
166 setPointDataView(temp);
167
168 // copy the default value from other to this
169 getDefaultValue().copySlice(other.getDefaultValue(), regionLoopRange);
170
171 // loop through the tag values copying these
172 DataMapType::const_iterator pos;
173 DataArrayView::ValueType::size_type tagOffset=getPointDataView().noValues();
174 for (pos=other.m_offsetLookup.begin();pos!=other.m_offsetLookup.end();pos++){
175 getPointDataView().copySlice(tagOffset,other.getPointDataView(),pos->second,regionLoopRange);
176 m_offsetLookup.insert(DataMapType::value_type(pos->first,tagOffset));
177 tagOffset+=getPointDataView().noValues();
178 }
179 }
180
181 void
182 DataTagged::setSlice(const DataAbstract* other,
183 const DataArrayView::RegionType& region)
184 {
185
186 // other must be another DataTagged object
187 // Data:setSlice implementation should ensure this
188 const DataTagged* otherTemp=dynamic_cast<const DataTagged*>(other);
189 if (otherTemp==0) {
190 throw DataException("Programming error - casting to DataTagged.");
191 }
192
193 // determine shape of the specified region
194 DataArrayView::ShapeType regionShape(DataArrayView::getResultSliceShape(region));
195
196 // modify region specification as needed to match rank of this object
197 DataArrayView::RegionLoopRangeType regionLoopRange=getSliceRegionLoopRange(region);
198
199 // ensure rank/shape of this object is compatible with specified region
200 if (getPointDataView().getRank()!=region.size()) {
201 throw DataException("Error - Invalid slice region.");
202 }
203 if (otherTemp->getPointDataView().getRank()>0 and !other->getPointDataView().checkShape(regionShape)) {
204 throw DataException (other->getPointDataView().createShapeErrorMessage(
205 "Error - Couldn't copy slice due to shape mismatch.",regionShape));
206 }
207
208 // copy slice from other default value to this default value
209 getDefaultValue().copySliceFrom(otherTemp->getDefaultValue(), regionLoopRange);
210
211 // loop through tag values in other, adding any which aren't in this, using default value
212 DataMapType::const_iterator pos;
213 for (pos=otherTemp->m_offsetLookup.begin();pos!=otherTemp->m_offsetLookup.end();pos++) {
214 if (!isCurrentTag(pos->first)) {
215 addTaggedValue(pos->first,getDefaultValue());
216 }
217 }
218
219 // loop through the tag values copying slices from other to this
220 for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++) {
221 getDataPointByTag(pos->first).copySliceFrom(otherTemp->getDataPointByTag(pos->first), regionLoopRange);
222 }
223
224 }
225
226 int
227 DataTagged::getTagNumber(int dpno)
228 {
229 //
230 // Get the number of samples and data-points per sample
231 int numSamples = getNumSamples();
232 int numDataPointsPerSample = getNumDPPSample();
233 int numDataPoints = numSamples * numDataPointsPerSample;
234
235 if (numDataPointsPerSample==0) {
236 throw DataException("DataTagged::getTagNumber error: no data-points associated with this object.");
237 }
238
239 if (dpno<0 || dpno>numDataPoints-1) {
240 throw DataException("DataTagged::getTagNumber error: invalid data-point number supplied.");
241 }
242
243 //
244 // Determine the sample number which corresponds to this data-point number
245 int sampleNo = dpno / numDataPointsPerSample;
246
247 //
248 // Determine the tag number which corresponds to this sample number
249 int tagNo = getFunctionSpace().getTagFromSampleNo(sampleNo);
250
251 //
252 // return the tag number
253 return(tagNo);
254 }
255
256 void
257 DataTagged::setTaggedValues(const TagListType& tagKeys,
258 const ValueListType& values)
259 {
260 addTaggedValues(tagKeys,values);
261 }
262
263 void
264 DataTagged::setTaggedValue(int tagKey,
265 const DataArrayView& value)
266 {
267 if (!getPointDataView().checkShape(value.getShape())) {
268 throw DataException(getPointDataView().createShapeErrorMessage(
269 "Error - Cannot setTaggedValue due to shape mismatch.", value.getShape()));
270 }
271 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
272 if (pos==m_offsetLookup.end()) {
273 // tag couldn't be found so use addTaggedValue
274 addTaggedValue(tagKey,value);
275 } else {
276 // copy the values into the data array at the offset determined by m_offsetLookup
277 int offset=pos->second;
278 for (int i=0; i<getPointDataView().noValues(); i++) {
279 m_data[offset+i]=value.getData(i);
280 }
281 }
282 }
283
284 void
285 DataTagged::addTaggedValues(const TagListType& tagKeys,
286 const ValueListType& values)
287 {
288 if (values.size()==0) {
289 // copy the current default value for each of the tags
290 TagListType::const_iterator iT;
291 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
292 // the point data view for DataTagged points at the default value
293 addTaggedValue(*iT,getPointDataView());
294 }
295 } else if (values.size()==1 && tagKeys.size()>1) {
296 // assume the one given value will be used for all tag values
297 TagListType::const_iterator iT;
298 for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
299 addTaggedValue(*iT,values[0]);
300 }
301 } else {
302 if (tagKeys.size()!=values.size()) {
303 stringstream temp;
304 temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
305 << " doesn't match number of values: " << values.size();
306 throw DataException(temp.str());
307 } else {
308 for (int i=0;i<tagKeys.size();i++) {
309 addTaggedValue(tagKeys[i],values[i]);
310 }
311 }
312 }
313 }
314
315 void
316 DataTagged::addTaggedValue(int tagKey,
317 const DataArrayView& value)
318 {
319 if (!getPointDataView().checkShape(value.getShape())) {
320 throw DataException(getPointDataView().createShapeErrorMessage(
321 "Error - Cannot addTaggedValue due to shape mismatch.", value.getShape()));
322 }
323 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
324 if (pos!=m_offsetLookup.end()) {
325 // tag already exists so use setTaggedValue
326 setTaggedValue(tagKey,value);
327 } else {
328 // save the key and the location of its data in the lookup tab
329 m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
330 // add the data given in "value" at the end of m_data
331 // need to make a temp copy of m_data, resize m_data, then copy
332 // all the old values plus the value to be added back into m_data
333 ValueType m_data_temp(m_data);
334 int oldSize=m_data.size();
335 int newSize=m_data.size()+value.noValues();
336 m_data.resize(newSize,0.,newSize);
337 for (int i=0;i<oldSize;i++) {
338 m_data[i]=m_data_temp[i];
339 }
340 for (int i=0;i<value.noValues();i++) {
341 m_data[oldSize+i]=value.getData(i);
342 }
343 }
344 }
345
346 double*
347 DataTagged::getSampleDataByTag(int tag)
348 {
349 DataMapType::iterator pos(m_offsetLookup.find(tag));
350 if (pos==m_offsetLookup.end()) {
351 // tag couldn't be found so return the default value
352 return &(m_data[0]);
353 } else {
354 // return the data-point corresponding to the given tag
355 return &(m_data[pos->second]);
356 }
357 }
358
359 string
360 DataTagged::toString() const
361 {
362 stringstream temp;
363 DataMapType::const_iterator i;
364 temp << "Tag(Default)" << endl;
365 temp << getDefaultValue().toString() << endl;
366 // create a temporary view as the offset will be changed
367 DataArrayView tempView(getPointDataView().getData(), getPointDataView().getShape());
368 for (i=m_offsetLookup.begin();i!=m_offsetLookup.end();++i) {
369 temp << "Tag(" << i->first << ")" << endl;
370 tempView.setOffset(i->second);
371 temp << tempView.toString() << endl;
372 }
373 return temp.str();
374 }
375
376 DataArrayView::ValueType::size_type
377 DataTagged::getPointOffset(int sampleNo,
378 int dataPointNo) const
379 {
380 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
381 DataMapType::const_iterator pos(m_offsetLookup.find(tagKey));
382 DataArrayView::ValueType::size_type offset=m_defaultValueOffset;
383 if (pos!=m_offsetLookup.end()) {
384 offset=pos->second;
385 }
386 return offset;
387 }
388
389 DataArrayView
390 DataTagged::getDataPointByTag(int tag) const
391 {
392 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
393 DataArrayView::ValueType::size_type offset=m_defaultValueOffset;
394 if (pos!=m_offsetLookup.end()) {
395 offset=pos->second;
396 }
397 DataArrayView temp(getPointDataView());
398 temp.setOffset(offset);
399 return temp;
400 }
401
402 DataArrayView
403 DataTagged::getDataPoint(int sampleNo,
404 int dataPointNo)
405 {
406 EsysAssert(validSampleNo(sampleNo),"(getDataPoint) Invalid sampleNo: " << sampleNo);
407 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
408 return getDataPointByTag(tagKey);
409 }
410
411 int
412 DataTagged::archiveData(ofstream& archiveFile,
413 const DataArrayView::ValueType::size_type noValues) const
414 {
415 return(m_data.archiveData(archiveFile, noValues));
416 }
417
418 int
419 DataTagged::extractData(ifstream& archiveFile,
420 const DataArrayView::ValueType::size_type noValues)
421 {
422 return(m_data.extractData(archiveFile, noValues));
423 }
424 void
425 DataTagged::eigenvalues(DataAbstract* ev)
426 {
427 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
428 if (temp_ev==0) {
429 throw DataException("Error - DataTagged::eigenvalues casting to DataTagged failed (propably a programming error).");
430 }
431 const DataTagged::DataMapType& thisLookup=getTagLookup();
432 DataTagged::DataMapType::const_iterator i;
433 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
434 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
435 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
436 DataArrayView thisView=getDataPointByTag(i->first);
437 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
438 cout << i->first << thisView(0,0) << "\n";
439 DataArrayView::eigenvalues(thisView,0,evView,0);
440 }
441 DataArrayView::eigenvalues(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
442 }
443 void
444 DataTagged::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
445 {
446 DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
447 if (temp_ev==0) {
448 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
449 }
450 DataTagged* temp_V=dynamic_cast<DataTagged*>(V);
451 if (temp_V==0) {
452 throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
453 }
454 const DataTagged::DataMapType& thisLookup=getTagLookup();
455 DataTagged::DataMapType::const_iterator i;
456 DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
457 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
458 temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
459 temp_V->addTaggedValue(i->first,temp_V->getDefaultValue());
460 DataArrayView thisView=getDataPointByTag(i->first);
461 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
462 DataArrayView VView=temp_V->getDataPointByTag(i->first);
463 DataArrayView::eigenvalues_and_eigenvectors(thisView,0,evView,0,VView,0,tol);
464 }
465 DataArrayView::eigenvalues_and_eigenvectors(getDefaultValue(),0,
466 temp_ev->getDefaultValue(),0,
467 temp_V->getDefaultValue(),0,
468 tol);
469
470
471 }
472
473
474 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26