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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/esys2/escript/src/Data/DataExpanded.cpp revision 97 by jgs, Tue Dec 14 05:39:33 2004 UTC trunk/escript/src/DataExpanded.cpp revision 474 by jgs, Mon Jan 30 04:23:44 2006 UTC
# Line 13  Line 13 
13   ******************************************************************************   ******************************************************************************
14  */  */
15    
16  #include "escript/Data/DataException.h"  #include "DataException.h"
17  #include "escript/Data/DataExpanded.h"  #include "DataExpanded.h"
18  #include "escript/Data/DataConstant.h"  #include "DataConstant.h"
19  #include "escript/Data/DataTagged.h"  #include "DataTagged.h"
20  #include "escript/Data/DataArrayView.h"  #include "DataArrayView.h"
21    
22  #include <boost/python/extract.hpp>  #include <boost/python/extract.hpp>
23    
# Line 35  DataExpanded::DataExpanded(const boost:: Line 35  DataExpanded::DataExpanded(const boost::
35  {  {
36    DataArrayView::ShapeType tempShape;    DataArrayView::ShapeType tempShape;
37    //    //
38    // extract the shape from the python array    // extract the shape of the python numarray
39    for (int i=0; i<value.getrank(); ++i) {    for (int i=0; i<value.getrank(); i++) {
     //cout << extract<int>(value.getshape()[i]) << endl;  
40      tempShape.push_back(extract<int>(value.getshape()[i]));      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());    initialise(tempShape,what.getNumSamples(),what.getNumDPPSample());
45    //    //
46    // copy the given value to every data point    // copy the given value to every data point
# Line 49  DataExpanded::DataExpanded(const boost:: Line 50  DataExpanded::DataExpanded(const boost::
50  DataExpanded::DataExpanded(const DataExpanded& other)  DataExpanded::DataExpanded(const DataExpanded& other)
51    : DataAbstract(other.getFunctionSpace()),    : DataAbstract(other.getFunctionSpace()),
52    m_data(other.m_data)    m_data(other.m_data)
53  {      {
54    //    //
55    // create the view for the data    // create the view for the data
56    DataArrayView temp(m_data.getData(),other.getPointDataView().getShape());    DataArrayView temp(m_data.getData(),other.getPointDataView().getShape());
# Line 58  DataExpanded::DataExpanded(const DataExp Line 59  DataExpanded::DataExpanded(const DataExp
59    
60  DataExpanded::DataExpanded(const DataConstant& other)  DataExpanded::DataExpanded(const DataConstant& other)
61    : DataAbstract(other.getFunctionSpace())    : DataAbstract(other.getFunctionSpace())
62  {      {
63      //
64      // initialise the data array for this object
65    initialise(other.getPointDataView().getShape(),other.getNumSamples(),other.getNumDPPSample());    initialise(other.getPointDataView().getShape(),other.getNumSamples(),other.getNumDPPSample());
66    //    //
67    // DataConstant only has one value    // DataConstant only has one value, copy this to every data point
68    copy(other.getPointDataView());    copy(other.getPointDataView());
69  }  }
70    
71  DataExpanded::DataExpanded(const DataTagged& other)  DataExpanded::DataExpanded(const DataTagged& other)
72    : DataAbstract(other.getFunctionSpace())    : DataAbstract(other.getFunctionSpace())
73  {      {
74      //
75      // initialise the data array for this object
76    initialise(other.getPointDataView().getShape(),other.getNumSamples(),other.getNumDPPSample());    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;    int i,j;
81    DataArrayView::ValueType::size_type numRows=m_data.getNumRows();    DataArrayView::ValueType::size_type numRows=m_data.getNumRows();
82    DataArrayView::ValueType::size_type numCols=m_data.getNumCols();    DataArrayView::ValueType::size_type numCols=m_data.getNumCols();
83  #pragma omp parallel for private(i,j) schedule(static)    #pragma omp parallel for private(i,j) schedule(static)
84    for (i=0;i<numRows;++i) {    for (i=0;i<numRows;i++) {
85      for (j=0;j<numCols;++j) {      for (j=0;j<numCols;j++) {
86        try {        try {
87          getPointDataView().copy(getPointOffset(i,j), other.getPointDataView(), other.getPointOffset(i,j));          getPointDataView().copy(getPointOffset(i,j),
88                                    other.getPointDataView(),
89                                    other.getPointOffset(i,j));
90        }        }
91        catch (std::exception& e) {        catch (std::exception& e) {
92          cout << e.what() << endl;          cout << e.what() << endl;
# Line 97  DataExpanded::DataExpanded(const DataExp Line 107  DataExpanded::DataExpanded(const DataExp
107    initialise(shape,other.getNumSamples(),other.getNumDPPSample());    initialise(shape,other.getNumSamples(),other.getNumDPPSample());
108    //    //
109    // copy the data    // copy the data
110      DataArrayView::RegionLoopRangeType region_loop_range=getSliceRegionLoopRange(region);
111    DataArrayView::ValueType::size_type numRows=m_data.getNumRows();    DataArrayView::ValueType::size_type numRows=m_data.getNumRows();
112    DataArrayView::ValueType::size_type numCols=m_data.getNumCols();    DataArrayView::ValueType::size_type numCols=m_data.getNumCols();
113    int i,j;    int i,j;
114  #pragma omp parallel for private(i,j) schedule(static)    #pragma omp parallel for private(i,j) schedule(static)
115    for (i=0;i<numRows;++i) {    for (i=0;i<numRows;i++) {
116      for (j=0;j<numCols;++j) {      for (j=0;j<numCols;j++) {
117        try {        try {
118          getPointDataView().copySlice(getPointOffset(i,j),other.getPointDataView(),other.getPointOffset(i,j),region);          getPointDataView().copySlice(getPointOffset(i,j),
119                                         other.getPointDataView(),
120                                         other.getPointOffset(i,j),
121                                         region_loop_range);
122        }        }
123        catch (std::exception& e) {        catch (std::exception& e) {
124          cout << e.what() << endl;          cout << e.what() << endl;
# Line 117  DataExpanded::DataExpanded(const DataArr Line 131  DataExpanded::DataExpanded(const DataArr
131                             const FunctionSpace& what)                             const FunctionSpace& what)
132    : DataAbstract(what)    : DataAbstract(what)
133  {  {
134      //
135      // get the shape of the given data value
136    DataArrayView::ShapeType tempShape=value.getShape();    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());    initialise(tempShape,what.getNumSamples(),what.getNumDPPSample());
140      //
141      // copy the given value to every data point
142    copy(value);    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      // create the view of the data
152      initialise(shape,what.getNumSamples(),what.getNumDPPSample());
153      //
154      // copy the data in the correct format
155      m_data.getData()=data;
156    }
157    
158  DataExpanded::~DataExpanded()  DataExpanded::~DataExpanded()
159  {  {
160  }  }
# Line 129  DataExpanded::~DataExpanded() Line 162  DataExpanded::~DataExpanded()
162  void  void
163  DataExpanded::reshapeDataPoint(const DataArrayView::ShapeType& shape)  DataExpanded::reshapeDataPoint(const DataArrayView::ShapeType& shape)
164  {  {
   //  
   // reshape a rank zero data point  
165    if (getPointDataView().getRank()!=0) {    if (getPointDataView().getRank()!=0) {
166      stringstream temp;      stringstream temp;
167      temp << "Error - Can only reshape Data with data points of rank 0. "      temp << "Error - Can only reshape Data with data points of rank 0. "
# Line 138  DataExpanded::reshapeDataPoint(const Dat Line 169  DataExpanded::reshapeDataPoint(const Dat
169           << getPointDataView().getRank();           << getPointDataView().getRank();
170      throw DataException(temp.str());      throw DataException(temp.str());
171    }    }
172    DataBlocks2D newData(getNumSamples(),getNumDPPSample(), DataArrayView::noValues(shape));    //
173      // create the new DataBlocks2D data array, and a corresponding DataArrayView
174      DataBlocks2D newData(getNumSamples(),getNumDPPSample(),DataArrayView::noValues(shape));
175    DataArrayView newView(newData.getData(),shape);    DataArrayView newView(newData.getData(),shape);
176    //    //
177    // Copy the original data to every value for the new shape    // Copy the original data to every value for the new shape
178    int i,j;    int i,j;
179    int nRows=m_data.getNumRows();    int nRows=m_data.getNumRows();
180    int nCols=m_data.getNumCols();    int nCols=m_data.getNumCols();
181  #pragma omp parallel for private(i,j) schedule(static)    #pragma omp parallel for private(i,j) schedule(static)
182    for (i=0;i<nRows;++i) {    for (i=0;i<nRows;i++) {
183      for (j=0;j<nCols;++j) {      for (j=0;j<nCols;j++) {
       //  
       // Copy the data into the specified offset  
184        // NOTE: An exception may be thown from this call if        // NOTE: An exception may be thown from this call if
185        // DOASSERT is on which of course will play        // DOASSERT is on which of course will play
186        // havoc with the omp threads. Run single threaded        // havoc with the omp threads. Run single threaded
187        // if using DOASSERT.        // if using DOASSERT.
188        newView.copy(newData.index(i,j),m_data.getData()[m_data.index(i,j)]);        newView.copy(newData.index(i,j),m_data(i,j));
189      }      }
190    }    }
191      // swap the new data array for the original
192    m_data.Swap(newData);    m_data.Swap(newData);
193      // set the corresponding DataArrayView
194    DataArrayView temp(m_data.getData(),shape);    DataArrayView temp(m_data.getData(),shape);
195    setPointDataView(temp);    setPointDataView(temp);
196  }  }
# Line 176  DataExpanded::setSlice(const DataAbstrac Line 209  DataExpanded::setSlice(const DataAbstrac
209    if (tempDataExp==0) {    if (tempDataExp==0) {
210      throw DataException("Programming error - casting to DataExpanded.");      throw DataException("Programming error - casting to DataExpanded.");
211    }    }
212      //
213      // get shape of slice
214      DataArrayView::ShapeType shape(DataArrayView::getResultSliceShape(region));
215      DataArrayView::RegionLoopRangeType region_loop_range=getSliceRegionLoopRange(region);
216      //
217      // check shape
218    if (getPointDataView().getRank()!=region.size()) {    if (getPointDataView().getRank()!=region.size()) {
219      throw DataException("Error - Invalid slice region.");      throw DataException("Error - Invalid slice region.");
220    }    }
221    if (!value->getPointDataView().checkShape(DataArrayView::getResultSliceShape(region))) {    if (tempDataExp->getPointDataView().getRank()>0 and !value->getPointDataView().checkShape(shape)) {
222      throw DataException (value->getPointDataView().createShapeErrorMessage(      throw DataException (value->getPointDataView().createShapeErrorMessage(
223          "Error - Couldn't copy slice due to shape mismatch.",DataArrayView::getResultSliceShape(region)));          "Error - Couldn't copy slice due to shape mismatch.",shape));
224    }    }
225    //    //
226    // copy the data    // copy the data from the slice into this object
227    DataArrayView::ValueType::size_type numRows=m_data.getNumRows();    DataArrayView::ValueType::size_type numRows=m_data.getNumRows();
228    DataArrayView::ValueType::size_type numCols=m_data.getNumCols();    DataArrayView::ValueType::size_type numCols=m_data.getNumCols();
229    int i, j;    int i, j;
230  #pragma omp parallel for private(i,j) schedule(static)    #pragma omp parallel for private(i,j) schedule(static)
231    for (i=0;i<numRows;i++) {    for (i=0;i<numRows;i++) {
232      for (j=0;j<numCols;j++) {      for (j=0;j<numCols;j++) {
233        getPointDataView().copySliceFrom(getPointOffset(i,j),tempDataExp->getPointDataView(),tempDataExp->getPointOffset(i,j),region);        getPointDataView().copySliceFrom(getPointOffset(i,j),
234                                           tempDataExp->getPointDataView(),
235                                           tempDataExp->getPointOffset(i,j),
236                                           region_loop_range);
237      }      }
238    }    }
239  }  }
# Line 200  void Line 242  void
242  DataExpanded::copy(const DataArrayView& value)  DataExpanded::copy(const DataArrayView& value)
243  {  {
244    //    //
245    // Copy a single value to every data point    // copy a single value to every data point in this object
246    int nRows=m_data.getNumRows();    int nRows=m_data.getNumRows();
247    int nCols=m_data.getNumCols();    int nCols=m_data.getNumCols();
248    int i, j;    int i,j;
249  #pragma omp parallel for private(i,j) schedule(static)    #pragma omp parallel for private(i,j) schedule(static)
250    for (i=0;i<nRows;++i) {    for (i=0;i<nRows;i++) {
251      for (j=0;j<nCols;++j) {      for (j=0;j<nCols;j++) {
       //  
       // Copy the data into the specified offset  
252        // NOTE: An exception may be thown from this call if        // NOTE: An exception may be thown from this call if
253        // DOASSERT is on which of course will play        // DOASSERT is on which of course will play
254        // havoc with the omp threads. Run single threaded        // havoc with the omp threads. Run single threaded
# Line 222  void Line 262  void
262  DataExpanded::copy(const boost::python::numeric::array& value)  DataExpanded::copy(const boost::python::numeric::array& value)
263  {  {
264    //    //
265    // first convert the numarray into a DataArrayView format    // first convert the numarray into a DataArray object
266    DataArray temp(value);    DataArray temp(value);
267    //    //
268    // check the input shape matches this shape, this will throw an exception    // check the input shape matches this shape
269    if (!getPointDataView().checkShape(temp.getView().getShape())) {    if (!getPointDataView().checkShape(temp.getView().getShape())) {
270      throw DataException(getPointDataView().createShapeErrorMessage(      throw DataException(getPointDataView().createShapeErrorMessage(
271                          "Error - (DataExpanded) Cannot copy due to shape mismatch.",                          "Error - (DataExpanded) Cannot copy due to shape mismatch.",
272                          temp.getView().getShape()));                          temp.getView().getShape()));
273    }    }
274    //    //
275    // now copy over the entire data structure    // now copy over the data
276    copy(temp.getView());    copy(temp.getView());
277  }  }
278    
# Line 242  DataExpanded::initialise(const DataArray Line 282  DataExpanded::initialise(const DataArray
282                           int noDataPointsPerSample)                           int noDataPointsPerSample)
283  {  {
284    //    //
285    // resize data to the required size    // resize data array to the required size
286    m_data.resize(noSamples,noDataPointsPerSample,DataArrayView::noValues(shape));    m_data.resize(noSamples,noDataPointsPerSample,DataArrayView::noValues(shape));
287    //    //
288    // create a point data viewer of the data    // create the data view of the data array
289    DataArrayView temp(m_data.getData(),shape);    DataArrayView temp(m_data.getData(),shape);
290    setPointDataView(temp);    setPointDataView(temp);
291  }  }
# Line 257  DataExpanded::toString() const Line 297  DataExpanded::toString() const
297    //    //
298    // create a temporary view as the offset will be changed    // create a temporary view as the offset will be changed
299    DataArrayView tempView(getPointDataView().getData(),getPointDataView().getShape(),getPointDataView().getOffset());    DataArrayView tempView(getPointDataView().getData(),getPointDataView().getShape(),getPointDataView().getOffset());
300    for (int i=0;i<m_data.getNumRows();++i) {    for (int i=0;i<m_data.getNumRows();i++) {
301      for (int j=0;j<m_data.getNumCols();++j) {      for (int j=0;j<m_data.getNumCols();j++) {
302        tempView.setOffset(m_data.index(i,j));        tempView.setOffset(m_data.index(i,j));
303        stringstream suffix;        stringstream suffix;
304        suffix << "(" << i << "," << j << ")";        suffix << "(" << i << "," << j << ")";
# Line 289  DataExpanded::getDataPoint(int sampleNo, Line 329  DataExpanded::getDataPoint(int sampleNo,
329  DataArrayView::ValueType::size_type  DataArrayView::ValueType::size_type
330  DataExpanded::getLength() const  DataExpanded::getLength() const
331  {  {
332    return m_data.getData().size();    return m_data.size();
333  }  }
334    
335    void
336    DataExpanded::setRefValue(int ref,
337                              const DataArray& value)
338    {
339      //
340      // Get the number of samples and data-points per sample.
341      int numSamples = getNumSamples();
342      int numDPPSample = getNumDPPSample();
343    
344      //
345      // Determine the sample number which corresponds to this reference number.
346      int sampleNo = -1;
347      int tempRef = -1;
348      for (int n=0; n<numSamples; n++) {
349        tempRef = getFunctionSpace().getReferenceNoFromSampleNo(n);
350        if (tempRef == ref) {
351          sampleNo = n;
352          break;
353        }
354      }
355      if (sampleNo == -1) {
356        throw DataException("DataExpanded::setRefValue error: invalid ref number supplied.");
357      }
358    
359      for (int n=0; n<numDPPSample; n++) {
360        //
361        // Get *each* data-point in the sample in turn.
362        DataArrayView pointView = getDataPoint(sampleNo, n);
363        //
364        // Assign the values in the DataArray to this data-point.
365        pointView.copy(value.getView());
366      }
367    }
368    
369    void
370    DataExpanded::getRefValue(int ref,
371                              DataArray& value)
372    {
373      //
374      // Get the number of samples and data-points per sample.
375      int numSamples = getNumSamples();
376      int numDPPSample = getNumDPPSample();
377    
378      //
379      // Determine the sample number which corresponds to this reference number
380      int sampleNo = -1;
381      int tempRef = -1;
382      for (int n=0; n<numSamples; n++) {
383        tempRef = getFunctionSpace().getReferenceNoFromSampleNo(n);
384        if (tempRef == ref) {
385          sampleNo = n;
386          break;
387        }
388      }
389      if (sampleNo == -1) {
390        throw DataException("DataExpanded::getRefValue error: invalid ref number supplied.");
391      }
392    
393      //
394      // Get the *first* data-point associated with this sample number.
395      DataArrayView pointView = getDataPoint(sampleNo, 0);
396    
397      //
398      // Load the values from this data-point into the DataArray
399      value.getView().copy(pointView);
400    }
401    
402    int
403    DataExpanded::archiveData(ofstream& archiveFile,
404                              const DataArrayView::ValueType::size_type noValues) const
405    {
406      return(m_data.archiveData(archiveFile, noValues));
407    }
408    
409    int
410    DataExpanded::extractData(ifstream& archiveFile,
411                              const DataArrayView::ValueType::size_type noValues)
412    {
413      return(m_data.extractData(archiveFile, noValues));
414    }
415    
416    void
417    DataExpanded::copyAll(const boost::python::numeric::array& value) {
418      //
419      // Get the number of samples and data-points per sample.
420      int numSamples = getNumSamples();
421      int numDataPointsPerSample = getNumDPPSample();
422      int dataPointRank = getPointDataView().getRank();
423      ShapeType dataPointShape = getPointDataView().getShape();
424      //
425      // check rank:
426      if (value.getrank()!=dataPointRank+1)
427           throw DataException("Rank of numarray does not match Data object rank");
428      if (value.getshape()[0]!=numSamples*numDataPointsPerSample)
429           throw DataException("leading dimension of numarray is too small");
430      //
431      int dataPoint = 0;
432      for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
433        for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
434          DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);
435          if (dataPointRank==0) {
436             dataPointView()=extract<double>(value[dataPoint]);
437          } else if (dataPointRank==1) {
438             for (int i=0; i<dataPointShape[0]; i++) {
439                dataPointView(i)=extract<double>(value[dataPoint][i]);
440             }
441          } else if (dataPointRank==2) {
442             for (int i=0; i<dataPointShape[0]; i++) {
443               for (int j=0; j<dataPointShape[1]; j++) {
444                 dataPointView(i,j)=extract<double>(value[dataPoint][i][j]);
445               }
446             }
447           } else if (dataPointRank==3) {
448             for (int i=0; i<dataPointShape[0]; i++) {
449               for (int j=0; j<dataPointShape[1]; j++) {
450                 for (int k=0; k<dataPointShape[2]; k++) {
451                     dataPointView(i,j,k)=extract<double>(value[dataPoint][i][j][k]);
452                 }
453               }
454             }
455           } else if (dataPointRank==4) {
456             for (int i=0; i<dataPointShape[0]; i++) {
457               for (int j=0; j<dataPointShape[1]; j++) {
458                 for (int k=0; k<dataPointShape[2]; k++) {
459                   for (int l=0; l<dataPointShape[3]; l++) {
460                     dataPointView(i,j,k,l)=extract<double>(value[dataPoint][i][j][k][l]);
461                   }
462                 }
463               }
464             }
465          }
466          dataPoint++;
467        }
468      }
469    }
470  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.97  
changed lines
  Added in v.474

  ViewVC Help
Powered by ViewVC 1.1.26