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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 126 - (show annotations)
Fri Jul 22 03:53:08 2005 UTC (14 years ago) by jgs
Original Path: trunk/esys2/escript/src/Data/DataExpanded.cpp
File size: 15188 byte(s)
Merge of development branch back to main trunk on 2005-07-22

1 // $Id$
2 /*
3 ******************************************************************************
4 * *
5 * COPYRIGHT ACcESS 2004 - All Rights Reserved *
6 * *
7 * This software is the property of ACcESS. No part of this code *
8 * may be copied in any form or by any means without the expressed written *
9 * consent of ACcESS. Copying, use or modification of this software *
10 * by any unauthorised person is illegal unless that person has a software *
11 * license agreement with ACcESS. *
12 * *
13 ******************************************************************************
14 */
15
16 #include "escript/Data/DataException.h"
17 #include "escript/Data/DataExpanded.h"
18 #include "escript/Data/DataConstant.h"
19 #include "escript/Data/DataTagged.h"
20 #include "escript/Data/DataArrayView.h"
21
22 #include <boost/python/extract.hpp>
23
24 #include <iostream>
25
26 using namespace std;
27 using namespace boost::python;
28 using namespace boost;
29
30 namespace escript {
31
32 DataExpanded::DataExpanded(const boost::python::numeric::array& value,
33 const FunctionSpace& what)
34 : DataAbstract(what)
35 {
36 DataArrayView::ShapeType tempShape;
37 //
38 // extract the shape of the python numarray
39 for (int i=0; i<value.getrank(); i++) {
40 tempShape.push_back(extract<int>(value.getshape()[i]));
41 }
42 //
43 // initialise the data array for this object
44 initialise(tempShape,what.getNumSamples(),what.getNumDPPSample());
45 //
46 // copy the given value to every data point
47 copy(value);
48 }
49
50 DataExpanded::DataExpanded(const DataExpanded& other)
51 : DataAbstract(other.getFunctionSpace()),
52 m_data(other.m_data)
53 {
54 //
55 // create the view for the data
56 DataArrayView temp(m_data.getData(),other.getPointDataView().getShape());
57 setPointDataView(temp);
58 }
59
60 DataExpanded::DataExpanded(const DataConstant& other)
61 : DataAbstract(other.getFunctionSpace())
62 {
63 //
64 // initialise the data array for this object
65 initialise(other.getPointDataView().getShape(),other.getNumSamples(),other.getNumDPPSample());
66 //
67 // DataConstant only has one value, copy this to every data point
68 copy(other.getPointDataView());
69 }
70
71 DataExpanded::DataExpanded(const DataTagged& other)
72 : DataAbstract(other.getFunctionSpace())
73 {
74 //
75 // initialise the data array for this object
76 initialise(other.getPointDataView().getShape(),other.getNumSamples(),other.getNumDPPSample());
77 //
78 // for each data point in this object, extract and copy the corresponding data
79 // value from the given DataTagged object
80 int i,j;
81 DataArrayView::ValueType::size_type numRows=m_data.getNumRows();
82 DataArrayView::ValueType::size_type numCols=m_data.getNumCols();
83 #pragma omp parallel for private(i,j) schedule(static)
84 for (i=0;i<numRows;i++) {
85 for (j=0;j<numCols;j++) {
86 try {
87 getPointDataView().copy(getPointOffset(i,j),
88 other.getPointDataView(),
89 other.getPointOffset(i,j));
90 }
91 catch (std::exception& e) {
92 cout << e.what() << endl;
93 }
94 }
95 }
96 }
97
98 DataExpanded::DataExpanded(const DataExpanded& other,
99 const DataArrayView::RegionType& region)
100 : DataAbstract(other.getFunctionSpace())
101 {
102 //
103 // get the shape of the slice
104 DataArrayView::ShapeType shape(DataArrayView::getResultSliceShape(region));
105 //
106 // initialise this Data object to the shape of the slice
107 initialise(shape,other.getNumSamples(),other.getNumDPPSample());
108 //
109 // copy the data
110 DataArrayView::RegionLoopRangeType region_loop_range=getSliceRegionLoopRange(region);
111 DataArrayView::ValueType::size_type numRows=m_data.getNumRows();
112 DataArrayView::ValueType::size_type numCols=m_data.getNumCols();
113 int i,j;
114 #pragma omp parallel for private(i,j) schedule(static)
115 for (i=0;i<numRows;i++) {
116 for (j=0;j<numCols;j++) {
117 try {
118 getPointDataView().copySlice(getPointOffset(i,j),
119 other.getPointDataView(),
120 other.getPointOffset(i,j),
121 region_loop_range);
122 }
123 catch (std::exception& e) {
124 cout << e.what() << endl;
125 }
126 }
127 }
128 }
129
130 DataExpanded::DataExpanded(const DataArrayView& value,
131 const FunctionSpace& what)
132 : DataAbstract(what)
133 {
134 //
135 // get the shape of the given data value
136 DataArrayView::ShapeType tempShape=value.getShape();
137 //
138 // initialise this Data object to the shape of the given data value
139 initialise(tempShape,what.getNumSamples(),what.getNumDPPSample());
140 //
141 // copy the given value to every data point
142 copy(value);
143 }
144
145 DataExpanded::DataExpanded(const FunctionSpace& what,
146 const DataArrayView::ShapeType &shape,
147 const DataArrayView::ValueType &data)
148 : DataAbstract(what)
149 {
150 //
151 // copy the data in the correct format
152 m_data.getData()=data;
153 //
154 // create the view of the data
155 DataArrayView tempView(m_data.getData(),shape);
156 setPointDataView(tempView);
157 }
158
159 DataExpanded::~DataExpanded()
160 {
161 }
162
163 void
164 DataExpanded::reshapeDataPoint(const DataArrayView::ShapeType& shape)
165 {
166 if (getPointDataView().getRank()!=0) {
167 stringstream temp;
168 temp << "Error - Can only reshape Data with data points of rank 0. "
169 << "This Data has data points with rank: "
170 << getPointDataView().getRank();
171 throw DataException(temp.str());
172 }
173 //
174 // create the new DataBlocks2D data array, and a corresponding DataArrayView
175 DataBlocks2D newData(getNumSamples(),getNumDPPSample(),DataArrayView::noValues(shape));
176 DataArrayView newView(newData.getData(),shape);
177 //
178 // Copy the original data to every value for the new shape
179 int i,j;
180 int nRows=m_data.getNumRows();
181 int nCols=m_data.getNumCols();
182 #pragma omp parallel for private(i,j) schedule(static)
183 for (i=0;i<nRows;i++) {
184 for (j=0;j<nCols;j++) {
185 // NOTE: An exception may be thown from this call if
186 // DOASSERT is on which of course will play
187 // havoc with the omp threads. Run single threaded
188 // if using DOASSERT.
189 newView.copy(newData.index(i,j),m_data(i,j));
190 }
191 }
192 // swap the new data array for the original
193 m_data.Swap(newData);
194 // set the corresponding DataArrayView
195 DataArrayView temp(m_data.getData(),shape);
196 setPointDataView(temp);
197 }
198
199 DataAbstract*
200 DataExpanded::getSlice(const DataArrayView::RegionType& region) const
201 {
202 return new DataExpanded(*this,region);
203 }
204
205 void
206 DataExpanded::setSlice(const DataAbstract* value,
207 const DataArrayView::RegionType& region)
208 {
209 const DataExpanded* tempDataExp=dynamic_cast<const DataExpanded*>(value);
210 if (tempDataExp==0) {
211 throw DataException("Programming error - casting to DataExpanded.");
212 }
213 //
214 // get shape of slice
215 DataArrayView::ShapeType shape(DataArrayView::getResultSliceShape(region));
216 DataArrayView::RegionLoopRangeType region_loop_range=getSliceRegionLoopRange(region);
217 //
218 // check shape
219 if (getPointDataView().getRank()!=region.size()) {
220 throw DataException("Error - Invalid slice region.");
221 }
222 if (tempDataExp->getPointDataView().getRank()>0 and !value->getPointDataView().checkShape(shape)) {
223 throw DataException (value->getPointDataView().createShapeErrorMessage(
224 "Error - Couldn't copy slice due to shape mismatch.",shape));
225 }
226 //
227 // copy the data from the slice into this object
228 DataArrayView::ValueType::size_type numRows=m_data.getNumRows();
229 DataArrayView::ValueType::size_type numCols=m_data.getNumCols();
230 int i, j;
231 #pragma omp parallel for private(i,j) schedule(static)
232 for (i=0;i<numRows;i++) {
233 for (j=0;j<numCols;j++) {
234 getPointDataView().copySliceFrom(getPointOffset(i,j),
235 tempDataExp->getPointDataView(),
236 tempDataExp->getPointOffset(i,j),
237 region_loop_range);
238 }
239 }
240 }
241
242 void
243 DataExpanded::copy(const DataArrayView& value)
244 {
245 //
246 // copy a single value to every data point in this object
247 int nRows=m_data.getNumRows();
248 int nCols=m_data.getNumCols();
249 int i,j;
250 #pragma omp parallel for private(i,j) schedule(static)
251 for (i=0;i<nRows;i++) {
252 for (j=0;j<nCols;j++) {
253 // NOTE: An exception may be thown from this call if
254 // DOASSERT is on which of course will play
255 // havoc with the omp threads. Run single threaded
256 // if using DOASSERT.
257 getPointDataView().copy(m_data.index(i,j),value);
258 }
259 }
260 }
261
262 void
263 DataExpanded::copy(const boost::python::numeric::array& value)
264 {
265 //
266 // first convert the numarray into a DataArray object
267 DataArray temp(value);
268 //
269 // check the input shape matches this shape
270 if (!getPointDataView().checkShape(temp.getView().getShape())) {
271 throw DataException(getPointDataView().createShapeErrorMessage(
272 "Error - (DataExpanded) Cannot copy due to shape mismatch.",
273 temp.getView().getShape()));
274 }
275 //
276 // now copy over the data
277 copy(temp.getView());
278 }
279
280 void
281 DataExpanded::initialise(const DataArrayView::ShapeType& shape,
282 int noSamples,
283 int noDataPointsPerSample)
284 {
285 //
286 // resize data array to the required size
287 m_data.resize(noSamples,noDataPointsPerSample,DataArrayView::noValues(shape));
288 //
289 // create the data view of the data array
290 DataArrayView temp(m_data.getData(),shape);
291 setPointDataView(temp);
292 }
293
294 string
295 DataExpanded::toString() const
296 {
297 stringstream temp;
298 //
299 // create a temporary view as the offset will be changed
300 DataArrayView tempView(getPointDataView().getData(),getPointDataView().getShape(),getPointDataView().getOffset());
301 for (int i=0;i<m_data.getNumRows();i++) {
302 for (int j=0;j<m_data.getNumCols();j++) {
303 tempView.setOffset(m_data.index(i,j));
304 stringstream suffix;
305 suffix << "(" << i << "," << j << ")";
306 temp << tempView.toString(suffix.str());
307 if (!(i==(m_data.getNumRows()-1) && j==(m_data.getNumCols()-1))) {
308 temp << endl;
309 }
310 }
311 }
312 return temp.str();
313 }
314
315 DataArrayView::ValueType::size_type
316 DataExpanded::getPointOffset(int sampleNo,
317 int dataPointNo) const
318 {
319 return m_data.index(sampleNo,dataPointNo);
320 }
321
322 DataArrayView
323 DataExpanded::getDataPoint(int sampleNo,
324 int dataPointNo)
325 {
326 DataArrayView temp(m_data.getData(),getPointDataView().getShape(),m_data.index(sampleNo,dataPointNo));
327 return temp;
328 }
329
330 DataArrayView::ValueType::size_type
331 DataExpanded::getLength() const
332 {
333 return m_data.size();
334 }
335
336 void
337 DataExpanded::setRefValue(int ref,
338 const DataArray& value)
339 {
340 //
341 // Get the number of samples and data-points per sample.
342 int numSamples = getNumSamples();
343 int numDPPSample = getNumDPPSample();
344
345 //
346 // Determine the sample number which corresponds to this reference number.
347 int sampleNo = -1;
348 int tempRef = -1;
349 for (int n=0; n<numSamples; n++) {
350 tempRef = getFunctionSpace().getReferenceNoFromSampleNo(n);
351 if (tempRef == ref) {
352 sampleNo = n;
353 break;
354 }
355 }
356 if (sampleNo == -1) {
357 throw DataException("DataExpanded::setRefValue error: invalid ref number supplied.");
358 }
359
360 for (int n=0; n<numDPPSample; n++) {
361 //
362 // Get *each* data-point in the sample in turn.
363 DataArrayView pointView = getDataPoint(sampleNo, n);
364 //
365 // Assign the values in the DataArray to this data-point.
366 pointView.copy(value.getView());
367 }
368 }
369
370 void
371 DataExpanded::getRefValue(int ref,
372 DataArray& value)
373 {
374 //
375 // Get the number of samples and data-points per sample.
376 int numSamples = getNumSamples();
377 int numDPPSample = getNumDPPSample();
378
379 //
380 // Determine the sample number which corresponds to this reference number
381 int sampleNo = -1;
382 int tempRef = -1;
383 for (int n=0; n<numSamples; n++) {
384 tempRef = getFunctionSpace().getReferenceNoFromSampleNo(n);
385 if (tempRef == ref) {
386 sampleNo = n;
387 break;
388 }
389 }
390 if (sampleNo == -1) {
391 throw DataException("DataExpanded::getRefValue error: invalid ref number supplied.");
392 }
393
394 //
395 // Get the *first* data-point associated with this sample number.
396 DataArrayView pointView = getDataPoint(sampleNo, 0);
397
398 //
399 // Load the values from this data-point into the DataArray
400 value.getView().copy(pointView);
401 }
402
403 int
404 DataExpanded::archiveData(ofstream& archiveFile,
405 const DataArrayView::ValueType::size_type noValues) const
406 {
407 return(m_data.archiveData(archiveFile, noValues));
408 }
409
410 int
411 DataExpanded::extractData(ifstream& archiveFile,
412 const DataArrayView::ValueType::size_type noValues)
413 {
414 return(m_data.extractData(archiveFile, noValues));
415 }
416
417 void
418 DataExpanded::copyAll(const boost::python::numeric::array& value) {
419 //
420 // Get the number of samples and data-points per sample.
421 int numSamples = getNumSamples();
422 int numDataPointsPerSample = getNumDPPSample();
423 int dataPointRank = getPointDataView().getRank();
424 ShapeType dataPointShape = getPointDataView().getShape();
425 //
426 // check rank:
427 if (value.getrank()!=dataPointRank+1)
428 throw DataException("Rank of numarray does not match Data object rank");
429 if (value.getshape()[0]!=numSamples*numDataPointsPerSample)
430 throw DataException("leading dimension of numarray is too small");
431 //
432 int dataPoint = 0;
433 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
434 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
435 DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);
436 if (dataPointRank==0) {
437 dataPointView()=extract<double>(value[dataPoint]);
438 } else if (dataPointRank==1) {
439 for (int i=0; i<dataPointShape[0]; i++) {
440 dataPointView(i)=extract<double>(value[dataPoint][i]);
441 }
442 } else if (dataPointRank==2) {
443 for (int i=0; i<dataPointShape[0]; i++) {
444 for (int j=0; j<dataPointShape[1]; j++) {
445 dataPointView(i,j)=extract<double>(value[dataPoint][i][j]);
446 }
447 }
448 } else if (dataPointRank==3) {
449 for (int i=0; i<dataPointShape[0]; i++) {
450 for (int j=0; j<dataPointShape[1]; j++) {
451 for (int k=0; k<dataPointShape[2]; k++) {
452 dataPointView(i,j,k)=extract<double>(value[dataPoint][i][j][k]);
453 }
454 }
455 }
456 } else if (dataPointRank==4) {
457 for (int i=0; i<dataPointShape[0]; i++) {
458 for (int j=0; j<dataPointShape[1]; j++) {
459 for (int k=0; k<dataPointShape[2]; k++) {
460 for (int l=0; l<dataPointShape[3]; l++) {
461 dataPointView(i,j,k,l)=extract<double>(value[dataPoint][i][j][k][l]);
462 }
463 }
464 }
465 }
466 }
467 dataPoint++;
468 }
469 }
470 }
471 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26