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

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

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

trunk/esys2/escript/src/Data/Data.cpp revision 119 by jgs, Tue Apr 12 04:45:05 2005 UTC trunk/escript/src/Data.cpp revision 474 by jgs, Mon Jan 30 04:23:44 2006 UTC
# Line 15  Line 15 
15    
16  ******************************************************************************/  ******************************************************************************/
17    
18  #include "escript/Data/Data.h"  #include "Data.h"
19    
20  #include <iostream>  #include <iostream>
21  #include <fstream>  #include <fstream>
# Line 25  Line 25 
25  #include <functional>  #include <functional>
26  #include <math.h>  #include <math.h>
27    
28    #include <boost/python/dict.hpp>
29  #include <boost/python/str.hpp>  #include <boost/python/str.hpp>
30  #include <boost/python/extract.hpp>  #include <boost/python/extract.hpp>
31  #include <boost/python/long.hpp>  #include <boost/python/long.hpp>
32    #include <boost/python/tuple.hpp>
33    
34  #include "escript/Data/DataException.h"  #include "DataException.h"
35  #include "escript/Data/DataExpanded.h"  #include "DataExpanded.h"
36  #include "escript/Data/DataConstant.h"  #include "DataConstant.h"
37  #include "escript/Data/DataTagged.h"  #include "DataTagged.h"
38  #include "escript/Data/DataEmpty.h"  #include "DataEmpty.h"
39  #include "escript/Data/DataArray.h"  #include "DataArray.h"
40  #include "escript/Data/DataAlgorithm.h"  #include "DataProf.h"
41  #include "escript/Data/FunctionSpaceFactory.h"  #include "FunctionSpaceFactory.h"
42  #include "escript/Data/AbstractContinuousDomain.h"  #include "AbstractContinuousDomain.h"
43  #include "escript/Data/UnaryFuncs.h"  #include "UnaryFuncs.h"
44    
45  using namespace std;  using namespace std;
46  using namespace boost::python;  using namespace boost::python;
47  using namespace boost;  using namespace boost;
48  using namespace escript;  using namespace escript;
49    
50    #if defined DOPROF
51    //
52    // global table of profiling data for all Data objects
53    DataProf dataProfTable;
54    #endif
55    
56  Data::Data()  Data::Data()
57  {  {
58    //    //
# Line 52  Data::Data() Line 60  Data::Data()
60    DataAbstract* temp=new DataEmpty();    DataAbstract* temp=new DataEmpty();
61    shared_ptr<DataAbstract> temp_data(temp);    shared_ptr<DataAbstract> temp_data(temp);
62    m_data=temp_data;    m_data=temp_data;
63    #if defined DOPROF
64      // create entry in global profiling table for this object
65      profData = dataProfTable.newData();
66    #endif
67  }  }
68    
69  Data::Data(double value,  Data::Data(double value,
# Line 65  Data::Data(double value, Line 77  Data::Data(double value,
77    }    }
78    DataArray temp(dataPointShape,value);    DataArray temp(dataPointShape,value);
79    initialise(temp.getView(),what,expanded);    initialise(temp.getView(),what,expanded);
80    #if defined DOPROF
81      // create entry in global profiling table for this object
82      profData = dataProfTable.newData();
83    #endif
84  }  }
85    
86  Data::Data(double value,  Data::Data(double value,
# Line 75  Data::Data(double value, Line 91  Data::Data(double value,
91    DataArray temp(dataPointShape,value);    DataArray temp(dataPointShape,value);
92    pair<int,int> dataShape=what.getDataShape();    pair<int,int> dataShape=what.getDataShape();
93    initialise(temp.getView(),what,expanded);    initialise(temp.getView(),what,expanded);
94    #if defined DOPROF
95      // create entry in global profiling table for this object
96      profData = dataProfTable.newData();
97    #endif
98  }  }
99    
100  Data::Data(const Data& inData)  Data::Data(const Data& inData)
101  {  {
102    m_data=inData.m_data;    m_data=inData.m_data;
103    #if defined DOPROF
104      // create entry in global profiling table for this object
105      profData = dataProfTable.newData();
106    #endif
107  }  }
108    
109  Data::Data(const Data& inData,  Data::Data(const Data& inData,
# Line 90  Data::Data(const Data& inData, Line 114  Data::Data(const Data& inData,
114    DataAbstract* tmp = inData.m_data->getSlice(region);    DataAbstract* tmp = inData.m_data->getSlice(region);
115    shared_ptr<DataAbstract> temp_data(tmp);    shared_ptr<DataAbstract> temp_data(tmp);
116    m_data=temp_data;    m_data=temp_data;
117    #if defined DOPROF
118      // create entry in global profiling table for this object
119      profData = dataProfTable.newData();
120    #endif
121  }  }
122    
123  Data::Data(const Data& inData,  Data::Data(const Data& inData,
# Line 99  Data::Data(const Data& inData, Line 127  Data::Data(const Data& inData,
127      m_data=inData.m_data;      m_data=inData.m_data;
128    } else {    } else {
129      Data tmp(0,inData.getPointDataView().getShape(),functionspace,true);      Data tmp(0,inData.getPointDataView().getShape(),functionspace,true);
130      // Note for Lutz, Must use a reference or pointer to a derived object      // Note: Must use a reference or pointer to a derived object
131      // in order to get polymorphic behaviour. Shouldn't really      // in order to get polymorphic behaviour. Shouldn't really
132      // be able to create an instance of AbstractDomain but that was done      // be able to create an instance of AbstractDomain but that was done
133      // as a boost python work around which may no longer be required.      // as a boost:python work around which may no longer be required.
134      const AbstractDomain& inDataDomain=inData.getDomain();      const AbstractDomain& inDataDomain=inData.getDomain();
135      if  (inDataDomain==functionspace.getDomain()) {      if  (inDataDomain==functionspace.getDomain()) {
136        inDataDomain.interpolateOnDomain(tmp,inData);        inDataDomain.interpolateOnDomain(tmp,inData);
# Line 111  Data::Data(const Data& inData, Line 139  Data::Data(const Data& inData,
139      }      }
140      m_data=tmp.m_data;      m_data=tmp.m_data;
141    }    }
142    #if defined DOPROF
143      // create entry in global profiling table for this object
144      profData = dataProfTable.newData();
145    #endif
146  }  }
147    
148  Data::Data(const DataTagged::TagListType& tagKeys,  Data::Data(const DataTagged::TagListType& tagKeys,
# Line 125  Data::Data(const DataTagged::TagListType Line 157  Data::Data(const DataTagged::TagListType
157    if (expanded) {    if (expanded) {
158      expand();      expand();
159    }    }
160    #if defined DOPROF
161      // create entry in global profiling table for this object
162      profData = dataProfTable.newData();
163    #endif
164  }  }
165    
166  Data::Data(const numeric::array& value,  Data::Data(const numeric::array& value,
# Line 132  Data::Data(const numeric::array& value, Line 168  Data::Data(const numeric::array& value,
168             bool expanded)             bool expanded)
169  {  {
170    initialise(value,what,expanded);    initialise(value,what,expanded);
171    #if defined DOPROF
172      // create entry in global profiling table for this object
173      profData = dataProfTable.newData();
174    #endif
175  }  }
176    
177  Data::Data(const DataArrayView& value,  Data::Data(const DataArrayView& value,
# Line 139  Data::Data(const DataArrayView& value, Line 179  Data::Data(const DataArrayView& value,
179             bool expanded)             bool expanded)
180  {  {
181    initialise(value,what,expanded);    initialise(value,what,expanded);
182    #if defined DOPROF
183      // create entry in global profiling table for this object
184      profData = dataProfTable.newData();
185    #endif
186  }  }
187    
188  Data::Data(const object& value,  Data::Data(const object& value,
# Line 147  Data::Data(const object& value, Line 191  Data::Data(const object& value,
191  {  {
192    numeric::array asNumArray(value);    numeric::array asNumArray(value);
193    initialise(asNumArray,what,expanded);    initialise(asNumArray,what,expanded);
194    #if defined DOPROF
195      // create entry in global profiling table for this object
196      profData = dataProfTable.newData();
197    #endif
198  }  }
199    
200  Data::Data(const object& value,  Data::Data(const object& value,
# Line 167  Data::Data(const object& value, Line 215  Data::Data(const object& value,
215      // Create a DataConstant with the same sample shape as other      // Create a DataConstant with the same sample shape as other
216      initialise(temp.getView(),other.getFunctionSpace(),false);      initialise(temp.getView(),other.getFunctionSpace(),false);
217    }    }
218    #if defined DOPROF
219      // create entry in global profiling table for this object
220      profData = dataProfTable.newData();
221    #endif
222    }
223    
224    Data::~Data()
225    {
226    
227  }  }
228    
229  escriptDataC  escriptDataC
# Line 185  Data::getDataC() const Line 242  Data::getDataC() const
242    return temp;    return temp;
243  }  }
244    
245  tuple  const boost::python::tuple
246  Data::getShapeTuple() const  Data::getShapeTuple() const
247  {  {
248    const DataArrayView::ShapeType& shape=getDataPointShape();    const DataArrayView::ShapeType& shape=getDataPointShape();
# Line 353  Data::reshapeDataPoint(const DataArrayVi Line 410  Data::reshapeDataPoint(const DataArrayVi
410  Data  Data
411  Data::wherePositive() const  Data::wherePositive() const
412  {  {
413    #if defined DOPROF
414      profData->where++;
415    #endif
416    return escript::unaryOp(*this,bind2nd(greater<double>(),0.0));    return escript::unaryOp(*this,bind2nd(greater<double>(),0.0));
417  }  }
418    
419  Data  Data
420  Data::whereNegative() const  Data::whereNegative() const
421  {  {
422    #if defined DOPROF
423      profData->where++;
424    #endif
425    return escript::unaryOp(*this,bind2nd(less<double>(),0.0));    return escript::unaryOp(*this,bind2nd(less<double>(),0.0));
426  }  }
427    
428  Data  Data
429  Data::whereNonNegative() const  Data::whereNonNegative() const
430  {  {
431    #if defined DOPROF
432      profData->where++;
433    #endif
434    return escript::unaryOp(*this,bind2nd(greater_equal<double>(),0.0));    return escript::unaryOp(*this,bind2nd(greater_equal<double>(),0.0));
435  }  }
436    
437  Data  Data
438  Data::whereNonPositive() const  Data::whereNonPositive() const
439  {  {
440    #if defined DOPROF
441      profData->where++;
442    #endif
443    return escript::unaryOp(*this,bind2nd(less_equal<double>(),0.0));    return escript::unaryOp(*this,bind2nd(less_equal<double>(),0.0));
444  }  }
445    
446  Data  Data
447  Data::whereZero() const  Data::whereZero() const
448  {  {
449    #if defined DOPROF
450      profData->where++;
451    #endif
452    return escript::unaryOp(*this,bind2nd(equal_to<double>(),0.0));    return escript::unaryOp(*this,bind2nd(equal_to<double>(),0.0));
453  }  }
454    
455  Data  Data
456  Data::whereNonZero() const  Data::whereNonZero() const
457  {  {
458    #if defined DOPROF
459      profData->where++;
460    #endif
461    return escript::unaryOp(*this,bind2nd(not_equal_to<double>(),0.0));    return escript::unaryOp(*this,bind2nd(not_equal_to<double>(),0.0));
462  }  }
463    
464  Data  Data
465  Data::interpolate(const FunctionSpace& functionspace) const  Data::interpolate(const FunctionSpace& functionspace) const
466  {  {
467    #if defined DOPROF
468      profData->interpolate++;
469    #endif
470    return Data(*this,functionspace);    return Data(*this,functionspace);
471  }  }
472    
# Line 410  Data::probeInterpolation(const FunctionS Line 488  Data::probeInterpolation(const FunctionS
488  Data  Data
489  Data::gradOn(const FunctionSpace& functionspace) const  Data::gradOn(const FunctionSpace& functionspace) const
490  {  {
491    #if defined DOPROF
492      profData->grad++;
493    #endif
494    if (functionspace.getDomain()!=getDomain())    if (functionspace.getDomain()!=getDomain())
495      throw DataException("Error - gradient cannot be calculated on different domains.");      throw DataException("Error - gradient cannot be calculated on different domains.");
496    DataArrayView::ShapeType grad_shape=getPointDataView().getShape();    DataArrayView::ShapeType grad_shape=getPointDataView().getShape();
# Line 443  Data::getDataPointShape() const Line 524  Data::getDataPointShape() const
524    return getPointDataView().getShape();    return getPointDataView().getShape();
525  }  }
526    
527    void
528    Data::fillFromNumArray(const boost::python::numeric::array num_array)
529    {
530      //
531      // check rank
532      if (num_array.getrank()<getDataPointRank())
533          throw DataException("Rank of numarray does not match Data object rank");
534    
535      //
536      // check shape of num_array
537      for (int i=0; i<getDataPointRank(); i++) {
538        if (extract<int>(num_array.getshape()[i+1])!=getDataPointShape()[i])
539           throw DataException("Shape of numarray does not match Data object rank");
540      }
541    
542      //
543      // make sure data is expanded:
544      if (!isExpanded()) {
545        expand();
546      }
547    
548      //
549      // and copy over
550      m_data->copyAll(num_array);
551    }
552    
553    const
554  boost::python::numeric::array  boost::python::numeric::array
555  Data::convertToNumArray()  Data::convertToNumArray()
556  {  {
557    //    //
558    // determine the current number of data points    // determine the total number of data points
559    int numSamples = getNumSamples();    int numSamples = getNumSamples();
560    int numDataPointsPerSample = getNumDataPointsPerSample();    int numDataPointsPerSample = getNumDataPointsPerSample();
561    int numDataPoints = numSamples * numDataPointsPerSample;    int numDataPoints = numSamples * numDataPointsPerSample;
# Line 499  Data::convertToNumArray() Line 607  Data::convertToNumArray()
607    int dataPoint = 0;    int dataPoint = 0;
608    for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {    for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
609      for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {      for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
   
610        DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);        DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);
   
611        if (dataPointRank==0) {        if (dataPointRank==0) {
612          numArray[dataPoint]=dataPointView();          numArray[dataPoint]=dataPointView();
613        }        }
   
614        if (dataPointRank==1) {        if (dataPointRank==1) {
615          for (int i=0; i<dataPointShape[0]; i++) {          for (int i=0; i<dataPointShape[0]; i++) {
616            numArray[dataPoint][i]=dataPointView(i);            numArray[dataPoint][i]=dataPointView(i);
617          }          }
618        }        }
   
619        if (dataPointRank==2) {        if (dataPointRank==2) {
620          for (int i=0; i<dataPointShape[0]; i++) {          for (int i=0; i<dataPointShape[0]; i++) {
621            for (int j=0; j<dataPointShape[1]; j++) {            for (int j=0; j<dataPointShape[1]; j++) {
# Line 519  Data::convertToNumArray() Line 623  Data::convertToNumArray()
623            }            }
624          }          }
625        }        }
   
626        if (dataPointRank==3) {        if (dataPointRank==3) {
627          for (int i=0; i<dataPointShape[0]; i++) {          for (int i=0; i<dataPointShape[0]; i++) {
628            for (int j=0; j<dataPointShape[1]; j++) {            for (int j=0; j<dataPointShape[1]; j++) {
# Line 529  Data::convertToNumArray() Line 632  Data::convertToNumArray()
632            }            }
633          }          }
634        }        }
   
635        if (dataPointRank==4) {        if (dataPointRank==4) {
636          for (int i=0; i<dataPointShape[0]; i++) {          for (int i=0; i<dataPointShape[0]; i++) {
637            for (int j=0; j<dataPointShape[1]; j++) {            for (int j=0; j<dataPointShape[1]; j++) {
# Line 541  Data::convertToNumArray() Line 643  Data::convertToNumArray()
643            }            }
644          }          }
645        }        }
   
646        dataPoint++;        dataPoint++;
647        }
648      }
649    
650      //
651      // return the loaded array
652      return numArray;
653    }
654    
655    const
656    boost::python::numeric::array
657    Data::convertToNumArrayFromSampleNo(int sampleNo)
658    {
659      //
660      // Check a valid sample number has been supplied
661      if (sampleNo >= getNumSamples()) {
662        throw DataException("Error - Data::convertToNumArray: invalid sampleNo.");
663      }
664    
665      //
666      // determine the number of data points per sample
667      int numDataPointsPerSample = getNumDataPointsPerSample();
668    
669      //
670      // determine the rank and shape of each data point
671      int dataPointRank = getDataPointRank();
672      DataArrayView::ShapeType dataPointShape = getDataPointShape();
673    
674      //
675      // create the numeric array to be returned
676      boost::python::numeric::array numArray(0.0);
677    
678      //
679      // the rank of the returned numeric array will be the rank of
680      // the data points, plus one. Where the rank of the array is n,
681      // the last n-1 dimensions will be equal to the shape of the
682      // data points, whilst the first dimension will be equal to the
683      // total number of data points. Thus the array will consist of
684      // a serial vector of the data points.
685      int arrayRank = dataPointRank + 1;
686      DataArrayView::ShapeType arrayShape;
687      arrayShape.push_back(numDataPointsPerSample);
688      for (int d=0; d<dataPointRank; d++) {
689         arrayShape.push_back(dataPointShape[d]);
690      }
691    
692      //
693      // resize the numeric array to the shape just calculated
694      if (arrayRank==1) {
695        numArray.resize(arrayShape[0]);
696      }
697      if (arrayRank==2) {
698        numArray.resize(arrayShape[0],arrayShape[1]);
699      }
700      if (arrayRank==3) {
701        numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2]);
702      }
703      if (arrayRank==4) {
704        numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3]);
705      }
706      if (arrayRank==5) {
707        numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3],arrayShape[4]);
708      }
709    
710      //
711      // loop through each data point in turn, loading the values for that data point
712      // into the numeric array.
713      for (int dataPoint = 0; dataPoint < numDataPointsPerSample; dataPoint++) {
714        DataArrayView dataPointView = getDataPoint(sampleNo, dataPoint);
715        if (dataPointRank==0) {
716          numArray[dataPoint]=dataPointView();
717        }
718        if (dataPointRank==1) {
719          for (int i=0; i<dataPointShape[0]; i++) {
720            numArray[dataPoint][i]=dataPointView(i);
721          }
722        }
723        if (dataPointRank==2) {
724          for (int i=0; i<dataPointShape[0]; i++) {
725            for (int j=0; j<dataPointShape[1]; j++) {
726              numArray[dataPoint][i][j] = dataPointView(i,j);
727            }
728          }
729        }
730        if (dataPointRank==3) {
731          for (int i=0; i<dataPointShape[0]; i++) {
732            for (int j=0; j<dataPointShape[1]; j++) {
733              for (int k=0; k<dataPointShape[2]; k++) {
734                numArray[dataPoint][i][j][k]=dataPointView(i,j,k);
735              }
736            }
737          }
738        }
739        if (dataPointRank==4) {
740          for (int i=0; i<dataPointShape[0]; i++) {
741            for (int j=0; j<dataPointShape[1]; j++) {
742              for (int k=0; k<dataPointShape[2]; k++) {
743                for (int l=0; l<dataPointShape[3]; l++) {
744                  numArray[dataPoint][i][j][k][l]=dataPointView(i,j,k,l);
745                }
746              }
747            }
748          }
749        }
750      }
751    
752      //
753      // return the loaded array
754      return numArray;
755    }
756    
757    const
758    boost::python::numeric::array
759    Data::convertToNumArrayFromDPNo(int sampleNo,
760                                    int dataPointNo)
761    {
762      //
763      // Check a valid sample number has been supplied
764      if (sampleNo >= getNumSamples()) {
765        throw DataException("Error - Data::convertToNumArray: invalid sampleNo.");
766      }
767    
768      //
769      // Check a valid data point number has been supplied
770      if (dataPointNo >= getNumDataPointsPerSample()) {
771        throw DataException("Error - Data::convertToNumArray: invalid dataPointNo.");
772      }
773    
774      //
775      // determine the rank and shape of each data point
776      int dataPointRank = getDataPointRank();
777      DataArrayView::ShapeType dataPointShape = getDataPointShape();
778    
779      //
780      // create the numeric array to be returned
781      boost::python::numeric::array numArray(0.0);
782    
783      //
784      // the shape of the returned numeric array will be the same
785      // as that of the data point
786      int arrayRank = dataPointRank;
787      DataArrayView::ShapeType arrayShape = dataPointShape;
788    
789      //
790      // resize the numeric array to the shape just calculated
791      if (arrayRank==0) {
792        numArray.resize(1);
793      }
794      if (arrayRank==1) {
795        numArray.resize(arrayShape[0]);
796      }
797      if (arrayRank==2) {
798        numArray.resize(arrayShape[0],arrayShape[1]);
799      }
800      if (arrayRank==3) {
801        numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2]);
802      }
803      if (arrayRank==4) {
804        numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3]);
805      }
806    
807      //
808      // load the values for the data point into the numeric array.
809      DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);
810      if (dataPointRank==0) {
811        numArray[0]=dataPointView();
812      }
813      if (dataPointRank==1) {
814        for (int i=0; i<dataPointShape[0]; i++) {
815          numArray[i]=dataPointView(i);
816        }
817      }
818      if (dataPointRank==2) {
819        for (int i=0; i<dataPointShape[0]; i++) {
820          for (int j=0; j<dataPointShape[1]; j++) {
821            numArray[i][j] = dataPointView(i,j);
822          }
823        }
824      }
825      if (dataPointRank==3) {
826        for (int i=0; i<dataPointShape[0]; i++) {
827          for (int j=0; j<dataPointShape[1]; j++) {
828            for (int k=0; k<dataPointShape[2]; k++) {
829              numArray[i][j][k]=dataPointView(i,j,k);
830            }
831          }
832        }
833      }
834      if (dataPointRank==4) {
835        for (int i=0; i<dataPointShape[0]; i++) {
836          for (int j=0; j<dataPointShape[1]; j++) {
837            for (int k=0; k<dataPointShape[2]; k++) {
838              for (int l=0; l<dataPointShape[3]; l++) {
839                numArray[i][j][k][l]=dataPointView(i,j,k,l);
840              }
841            }
842          }
843      }      }
844    }    }
845    
# Line 559  Data::integrate() const Line 855  Data::integrate() const
855    int rank = getDataPointRank();    int rank = getDataPointRank();
856    DataArrayView::ShapeType shape = getDataPointShape();    DataArrayView::ShapeType shape = getDataPointShape();
857    
858    #if defined DOPROF
859      profData->integrate++;
860    #endif
861    
862    //    //
863    // calculate the integral values    // calculate the integral values
864    vector<double> integrals(getDataPointSize());    vector<double> integrals(getDataPointSize());
# Line 581  Data::integrate() const Line 881  Data::integrate() const
881      }      }
882    }    }
883    if (rank==2) {    if (rank==2) {
884      bp_array.resize(shape[0],shape[1]);         bp_array.resize(shape[0],shape[1]);
885      for (int i=0; i<shape[0]; i++) {         for (int i=0; i<shape[0]; i++) {
886        for (int j=0; j<shape[1]; j++) {           for (int j=0; j<shape[1]; j++) {
887          index = i + shape[0] * j;             index = i + shape[0] * j;
888          bp_array[i,j] = integrals[index];             bp_array[make_tuple(i,j)] = integrals[index];
889        }           }
890      }         }
891    }    }
892    if (rank==3) {    if (rank==3) {
893      bp_array.resize(shape[0],shape[1],shape[2]);      bp_array.resize(shape[0],shape[1],shape[2]);
# Line 595  Data::integrate() const Line 895  Data::integrate() const
895        for (int j=0; j<shape[1]; j++) {        for (int j=0; j<shape[1]; j++) {
896          for (int k=0; k<shape[2]; k++) {          for (int k=0; k<shape[2]; k++) {
897            index = i + shape[0] * ( j + shape[1] * k );            index = i + shape[0] * ( j + shape[1] * k );
898            bp_array[i,j,k] = integrals[index];            bp_array[make_tuple(i,j,k)] = integrals[index];
899          }          }
900        }        }
901      }      }
# Line 607  Data::integrate() const Line 907  Data::integrate() const
907          for (int k=0; k<shape[2]; k++) {          for (int k=0; k<shape[2]; k++) {
908            for (int l=0; l<shape[3]; l++) {            for (int l=0; l<shape[3]; l++) {
909              index = i + shape[0] * ( j + shape[1] * ( k + shape[2] * l ) );              index = i + shape[0] * ( j + shape[1] * ( k + shape[2] * l ) );
910              bp_array[i,j,k,l] = integrals[index];              bp_array[make_tuple(i,j,k,l)] = integrals[index];
911            }            }
912          }          }
913        }        }
# Line 622  Data::integrate() const Line 922  Data::integrate() const
922  Data  Data
923  Data::sin() const  Data::sin() const
924  {  {
925    #if defined DOPROF
926      profData->unary++;
927    #endif
928    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sin);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sin);
929  }  }
930    
931  Data  Data
932  Data::cos() const  Data::cos() const
933  {  {
934    #if defined DOPROF
935      profData->unary++;
936    #endif
937    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::cos);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::cos);
938  }  }
939    
940  Data  Data
941  Data::tan() const  Data::tan() const
942  {  {
943    #if defined DOPROF
944      profData->unary++;
945    #endif
946    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::tan);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::tan);
947  }  }
948    
949  Data  Data
950  Data::log() const  Data::asin() const
951    {
952    #if defined DOPROF
953      profData->unary++;
954    #endif
955      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::asin);
956    }
957    
958    Data
959    Data::acos() const
960    {
961    #if defined DOPROF
962      profData->unary++;
963    #endif
964      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::acos);
965    }
966    
967    Data
968    Data::atan() const
969    {
970    #if defined DOPROF
971      profData->unary++;
972    #endif
973      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::atan);
974    }
975    
976    Data
977    Data::sinh() const
978  {  {
979    #if defined DOPROF
980      profData->unary++;
981    #endif
982      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sinh);
983    }
984    
985    Data
986    Data::cosh() const
987    {
988    #if defined DOPROF
989      profData->unary++;
990    #endif
991      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::cosh);
992    }
993    
994    Data
995    Data::tanh() const
996    {
997    #if defined DOPROF
998      profData->unary++;
999    #endif
1000      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::tanh);
1001    }
1002    
1003    Data
1004    Data::asinh() const
1005    {
1006    #if defined DOPROF
1007      profData->unary++;
1008    #endif
1009      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::asinh);
1010    }
1011    
1012    Data
1013    Data::acosh() const
1014    {
1015    #if defined DOPROF
1016      profData->unary++;
1017    #endif
1018      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::acosh);
1019    }
1020    
1021    Data
1022    Data::atanh() const
1023    {
1024    #if defined DOPROF
1025      profData->unary++;
1026    #endif
1027      return escript::unaryOp(*this,(Data::UnaryDFunPtr)::atanh);
1028    }
1029    
1030    Data
1031    Data::log10() const
1032    {
1033    #if defined DOPROF
1034      profData->unary++;
1035    #endif
1036    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log10);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log10);
1037  }  }
1038    
1039  Data  Data
1040  Data::ln() const  Data::log() const
1041  {  {
1042    #if defined DOPROF
1043      profData->unary++;
1044    #endif
1045    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log);
1046  }  }
1047    
1048  Data  Data
1049  Data::sign() const  Data::sign() const
1050  {  {
1051    #if defined DOPROF
1052      profData->unary++;
1053    #endif
1054    return escript::unaryOp(*this,escript::fsign);    return escript::unaryOp(*this,escript::fsign);
1055  }  }
1056    
1057  Data  Data
1058  Data::abs() const  Data::abs() const
1059  {  {
1060    #if defined DOPROF
1061      profData->unary++;
1062    #endif
1063    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::fabs);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::fabs);
1064  }  }
1065    
1066  Data  Data
1067  Data::neg() const  Data::neg() const
1068  {  {
1069    #if defined DOPROF
1070      profData->unary++;
1071    #endif
1072    return escript::unaryOp(*this,negate<double>());    return escript::unaryOp(*this,negate<double>());
1073  }  }
1074    
1075  Data  Data
1076  Data::pos() const  Data::pos() const
1077  {  {
1078    return (*this);  #if defined DOPROF
1079      profData->unary++;
1080    #endif
1081      Data result;
1082      // perform a deep copy
1083      result.copy(*this);
1084      return result;
1085  }  }
1086    
1087  Data  Data
1088  Data::exp() const  Data::exp() const
1089  {  {
1090    #if defined DOPROF
1091      profData->unary++;
1092    #endif
1093    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::exp);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::exp);
1094  }  }
1095    
1096  Data  Data
1097  Data::sqrt() const  Data::sqrt() const
1098  {  {
1099    #if defined DOPROF
1100      profData->unary++;
1101    #endif
1102    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sqrt);    return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sqrt);
1103  }  }
1104    
1105  double  double
1106  Data::Lsup() const  Data::Lsup() const
1107  {  {
1108    #if defined DOPROF
1109      profData->reduction1++;
1110    #endif
1111    //    //
1112    // set the initial absolute maximum value to zero    // set the initial absolute maximum value to zero
1113    return algorithm(DataAlgorithmAdapter<AbsMax>(0));    AbsMax abs_max_func;
1114      return algorithm(abs_max_func,0);
1115  }  }
1116    
1117  double  double
1118  Data::Linf() const  Data::Linf() const
1119  {  {
1120    #if defined DOPROF
1121      profData->reduction1++;
1122    #endif
1123    //    //
1124    // set the initial absolute minimum value to max double    // set the initial absolute minimum value to max double
1125    return algorithm(DataAlgorithmAdapter<AbsMin>(numeric_limits<double>::max()));    AbsMin abs_min_func;
1126      return algorithm(abs_min_func,numeric_limits<double>::max());
1127  }  }
1128    
1129  double  double
1130  Data::sup() const  Data::sup() const
1131  {  {
1132    #if defined DOPROF
1133      profData->reduction1++;
1134    #endif
1135    //    //
1136    // set the initial maximum value to min possible double    // set the initial maximum value to min possible double
1137    return algorithm(DataAlgorithmAdapter<FMax>(numeric_limits<double>::max()*-1));    FMax fmax_func;
1138      return algorithm(fmax_func,numeric_limits<double>::max()*-1);
1139  }  }
1140    
1141  double  double
1142  Data::inf() const  Data::inf() const
1143  {  {
1144    #if defined DOPROF
1145      profData->reduction1++;
1146    #endif
1147    //    //
1148    // set the initial minimum value to max possible double    // set the initial minimum value to max possible double
1149    return algorithm(DataAlgorithmAdapter<FMin>(numeric_limits<double>::max()));    FMin fmin_func;
1150      return algorithm(fmin_func,numeric_limits<double>::max());
1151  }  }
1152    
1153  Data  Data
1154  Data::maxval() const  Data::maxval() const
1155  {  {
1156    #if defined DOPROF
1157      profData->reduction2++;
1158    #endif
1159    //    //
1160    // set the initial maximum value to min possible double    // set the initial maximum value to min possible double
1161    return dp_algorithm(DataAlgorithmAdapter<FMax>(numeric_limits<double>::max()*-1));    FMax fmax_func;
1162      return dp_algorithm(fmax_func,numeric_limits<double>::max()*-1);
1163  }  }
1164    
1165  Data  Data
1166  Data::minval() const  Data::minval() const
1167  {  {
1168    #if defined DOPROF
1169      profData->reduction2++;
1170    #endif
1171    //    //
1172    // set the initial minimum value to max possible double    // set the initial minimum value to max possible double
1173    return dp_algorithm(DataAlgorithmAdapter<FMin>(numeric_limits<double>::max()));    FMin fmin_func;
1174  }    return dp_algorithm(fmin_func,numeric_limits<double>::max());
   
 Data  
 Data::length() const  
 {  
   return dp_algorithm(DataAlgorithmAdapter<Length>(0));  
1175  }  }
1176    
1177  Data  Data
1178  Data::trace() const  Data::trace() const
1179  {  {
1180    return dp_algorithm(DataAlgorithmAdapter<Trace>(0));  #if defined DOPROF
1181      profData->reduction2++;
1182    #endif
1183      Trace trace_func;
1184      return dp_algorithm(trace_func,0);
1185  }  }
1186    
1187  Data  Data
1188  Data::transpose(int axis) const  Data::transpose(int axis) const
1189  {  {
1190    #if defined DOPROF
1191      profData->reduction2++;
1192    #endif
1193    // not implemented    // not implemented
1194    throw DataException("Error - Data::transpose not implemented yet.");    throw DataException("Error - Data::transpose not implemented yet.");
1195    return Data();    return Data();
1196  }  }
1197    
1198    const boost::python::tuple
1199    Data::mindp() const
1200    {
1201      // NB: calc_mindp had to be split off from mindp as boost::make_tuple causes an
1202      // abort (for unknown reasons) if there are openmp directives with it in the
1203      // surrounding function
1204    
1205      int SampleNo;
1206      int DataPointNo;
1207    
1208      calc_mindp(SampleNo,DataPointNo);
1209    
1210      return make_tuple(SampleNo,DataPointNo);
1211    }
1212    
1213    void
1214    Data::calc_mindp(int& SampleNo,
1215                     int& DataPointNo) const
1216    {
1217      int i,j;
1218      int lowi=0,lowj=0;
1219      double min=numeric_limits<double>::max();
1220    
1221      Data temp=minval();
1222    
1223      int numSamples=temp.getNumSamples();
1224      int numDPPSample=temp.getNumDataPointsPerSample();
1225    
1226      double next,local_min;
1227      int local_lowi,local_lowj;
1228    
1229      #pragma omp parallel private(next,local_min,local_lowi,local_lowj)
1230      {
1231        local_min=min;
1232        #pragma omp for private(i,j) schedule(static)
1233        for (i=0; i<numSamples; i++) {
1234          for (j=0; j<numDPPSample; j++) {
1235            next=temp.getDataPoint(i,j)();
1236            if (next<local_min) {
1237              local_min=next;
1238              local_lowi=i;
1239              local_lowj=j;
1240            }
1241          }
1242        }
1243        #pragma omp critical
1244        if (local_min<min) {
1245          min=local_min;
1246          lowi=local_lowi;
1247          lowj=local_lowj;
1248        }
1249      }
1250    
1251      SampleNo = lowi;
1252      DataPointNo = lowj;
1253    }
1254    
1255  void  void
1256  Data::saveDX(std::string fileName) const  Data::saveDX(std::string fileName) const
1257  {  {
1258    getDomain().saveDX(fileName,*this);    boost::python::dict args;
1259      args["data"]=boost::python::object(this);
1260      getDomain().saveDX(fileName,args);
1261    return;    return;
1262  }  }
1263    
1264  void  void
1265  Data::saveVTK(std::string fileName) const  Data::saveVTK(std::string fileName) const
1266  {  {
1267    getDomain().saveVTK(fileName,*this);    boost::python::dict args;
1268      args["data"]=boost::python::object(this);
1269      getDomain().saveVTK(fileName,args);
1270    return;    return;
1271  }  }
1272    
1273  Data&  Data&
1274  Data::operator+=(const Data& right)  Data::operator+=(const Data& right)
1275  {  {
1276    #if defined DOPROF
1277      profData->binary++;
1278    #endif
1279    binaryOp(right,plus<double>());    binaryOp(right,plus<double>());
1280    return (*this);    return (*this);
1281  }  }
# Line 777  Data::operator+=(const Data& right) Line 1283  Data::operator+=(const Data& right)
1283  Data&  Data&
1284  Data::operator+=(const boost::python::object& right)  Data::operator+=(const boost::python::object& right)
1285  {  {
1286    #if defined DOPROF
1287      profData->binary++;
1288    #endif
1289    binaryOp(right,plus<double>());    binaryOp(right,plus<double>());
1290    return (*this);    return (*this);
1291  }  }
# Line 784  Data::operator+=(const boost::python::ob Line 1293  Data::operator+=(const boost::python::ob
1293  Data&  Data&
1294  Data::operator-=(const Data& right)  Data::operator-=(const Data& right)
1295  {  {
1296    #if defined DOPROF
1297      profData->binary++;
1298    #endif
1299    binaryOp(right,minus<double>());    binaryOp(right,minus<double>());
1300    return (*this);    return (*this);
1301  }  }
# Line 791  Data::operator-=(const Data& right) Line 1303  Data::operator-=(const Data& right)
1303  Data&  Data&
1304  Data::operator-=(const boost::python::object& right)  Data::operator-=(const boost::python::object& right)
1305  {  {
1306    #if defined DOPROF
1307      profData->binary++;
1308    #endif
1309    binaryOp(right,minus<double>());    binaryOp(right,minus<double>());
1310    return (*this);    return (*this);
1311  }  }
# Line 798  Data::operator-=(const boost::python::ob Line 1313  Data::operator-=(const boost::python::ob
1313  Data&  Data&
1314  Data::operator*=(const Data& right)  Data::operator*=(const Data& right)
1315  {  {
1316    #if defined DOPROF
1317      profData->binary++;
1318    #endif
1319    binaryOp(right,multiplies<double>());    binaryOp(right,multiplies<double>());
1320    return (*this);    return (*this);
1321  }  }
# Line 805  Data::operator*=(const Data& right) Line 1323  Data::operator*=(const Data& right)
1323  Data&  Data&
1324  Data::operator*=(const boost::python::object& right)  Data::operator*=(const boost::python::object& right)
1325  {  {
1326    #if defined DOPROF
1327      profData->binary++;
1328    #endif
1329    binaryOp(right,multiplies<double>());    binaryOp(right,multiplies<double>());
1330    return (*this);    return (*this);
1331  }  }
# Line 812  Data::operator*=(const boost::python::ob Line 1333  Data::operator*=(const boost::python::ob
1333  Data&  Data&
1334  Data::operator/=(const Data& right)  Data::operator/=(const Data& right)
1335  {  {
1336    #if defined DOPROF
1337      profData->binary++;
1338    #endif
1339    binaryOp(right,divides<double>());    binaryOp(right,divides<double>());
1340    return (*this);    return (*this);
1341  }  }
# Line 819  Data::operator/=(const Data& right) Line 1343  Data::operator/=(const Data& right)
1343  Data&  Data&
1344  Data::operator/=(const boost::python::object& right)  Data::operator/=(const boost::python::object& right)
1345  {  {
1346    #if defined DOPROF
1347      profData->binary++;
1348    #endif
1349    binaryOp(right,divides<double>());    binaryOp(right,divides<double>());
1350    return (*this);    return (*this);
1351  }  }
# Line 826  Data::operator/=(const boost::python::ob Line 1353  Data::operator/=(const boost::python::ob
1353  Data  Data
1354  Data::powO(const boost::python::object& right) const  Data::powO(const boost::python::object& right) const
1355  {  {
1356    #if defined DOPROF
1357      profData->binary++;
1358    #endif
1359    Data result;    Data result;
1360    result.copy(*this);    result.copy(*this);
1361    result.binaryOp(right,(Data::BinaryDFunPtr)::pow);    result.binaryOp(right,(Data::BinaryDFunPtr)::pow);
# Line 835  Data::powO(const boost::python::object& Line 1365  Data::powO(const boost::python::object&
1365  Data  Data
1366  Data::powD(const Data& right) const  Data::powD(const Data& right) const
1367  {  {
1368    #if defined DOPROF
1369      profData->binary++;
1370    #endif
1371    Data result;    Data result;
1372    result.copy(*this);    result.copy(*this);
1373    result.binaryOp(right,(Data::BinaryDFunPtr)::pow);    result.binaryOp(right,(Data::BinaryDFunPtr)::pow);
# Line 842  Data::powD(const Data& right) const Line 1375  Data::powD(const Data& right) const
1375  }  }
1376    
1377  //  //
1378  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1379  Data  Data
1380  escript::operator+(const Data& left, const Data& right)  escript::operator+(const Data& left, const Data& right)
1381  {  {
# Line 855  escript::operator+(const Data& left, con Line 1388  escript::operator+(const Data& left, con
1388  }  }
1389    
1390  //  //
1391  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1392  Data  Data
1393  escript::operator-(const Data& left, const Data& right)  escript::operator-(const Data& left, const Data& right)
1394  {  {
# Line 868  escript::operator-(const Data& left, con Line 1401  escript::operator-(const Data& left, con
1401  }  }
1402    
1403  //  //
1404  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1405  Data  Data
1406  escript::operator*(const Data& left, const Data& right)  escript::operator*(const Data& left, const Data& right)
1407  {  {
# Line 881  escript::operator*(const Data& left, con Line 1414  escript::operator*(const Data& left, con
1414  }  }
1415    
1416  //  //
1417  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1418  Data  Data
1419  escript::operator/(const Data& left, const Data& right)  escript::operator/(const Data& left, const Data& right)
1420  {  {
# Line 894  escript::operator/(const Data& left, con Line 1427  escript::operator/(const Data& left, con
1427  }  }
1428    
1429  //  //
1430  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1431  Data  Data
1432  escript::operator+(const Data& left, const boost::python::object& right)  escript::operator+(const Data& left, const boost::python::object& right)
1433  {  {
# Line 910  escript::operator+(const Data& left, con Line 1443  escript::operator+(const Data& left, con
1443  }  }
1444    
1445  //  //
1446  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1447  Data  Data
1448  escript::operator-(const Data& left, const boost::python::object& right)  escript::operator-(const Data& left, const boost::python::object& right)
1449  {  {
# Line 926  escript::operator-(const Data& left, con Line 1459  escript::operator-(const Data& left, con
1459  }  }
1460    
1461  //  //
1462  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1463  Data  Data
1464  escript::operator*(const Data& left, const boost::python::object& right)  escript::operator*(const Data& left, const boost::python::object& right)
1465  {  {
# Line 942  escript::operator*(const Data& left, con Line 1475  escript::operator*(const Data& left, con
1475  }  }
1476    
1477  //  //
1478  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1479  Data  Data
1480  escript::operator/(const Data& left, const boost::python::object& right)  escript::operator/(const Data& left, const boost::python::object& right)
1481  {  {
# Line 958  escript::operator/(const Data& left, con Line 1491  escript::operator/(const Data& left, con
1491  }  }
1492    
1493  //  //
1494  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1495  Data  Data
1496  escript::operator+(const boost::python::object& left, const Data& right)  escript::operator+(const boost::python::object& left, const Data& right)
1497  {  {
# Line 971  escript::operator+(const boost::python:: Line 1504  escript::operator+(const boost::python::
1504  }  }
1505    
1506  //  //
1507  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1508  Data  Data
1509  escript::operator-(const boost::python::object& left, const Data& right)  escript::operator-(const boost::python::object& left, const Data& right)
1510  {  {
# Line 984  escript::operator-(const boost::python:: Line 1517  escript::operator-(const boost::python::
1517  }  }
1518    
1519  //  //
1520  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1521  Data  Data
1522  escript::operator*(const boost::python::object& left, const Data& right)  escript::operator*(const boost::python::object& left, const Data& right)
1523  {  {
# Line 997  escript::operator*(const boost::python:: Line 1530  escript::operator*(const boost::python::
1530  }  }
1531    
1532  //  //
1533  // NOTE: It is essential to specify the namepsace this operator belongs to  // NOTE: It is essential to specify the namespace this operator belongs to
1534  Data  Data
1535  escript::operator/(const boost::python::object& left, const Data& right)  escript::operator/(const boost::python::object& left, const Data& right)
1536  {  {
# Line 1010  escript::operator/(const boost::python:: Line 1543  escript::operator/(const boost::python::
1543  }  }
1544    
1545  //  //
 // NOTE: It is essential to specify the namepsace this operator belongs to  
1546  //bool escript::operator==(const Data& left, const Data& right)  //bool escript::operator==(const Data& left, const Data& right)
1547  //{  //{
1548  //  /*  //  /*
# Line 1072  Data::getItem(const boost::python::objec Line 1604  Data::getItem(const boost::python::objec
1604  Data  Data
1605  Data::getSlice(const DataArrayView::RegionType& region) const  Data::getSlice(const DataArrayView::RegionType& region) const
1606  {  {
1607    #if defined DOPROF
1608      profData->slicing++;
1609    #endif
1610    return Data(*this,region);    return Data(*this,region);
1611  }  }
1612    
# Line 1088  Data::setItemD(const boost::python::obje Line 1623  Data::setItemD(const boost::python::obje
1623                 const Data& value)                 const Data& value)
1624  {  {
1625    const DataArrayView& view=getPointDataView();    const DataArrayView& view=getPointDataView();
1626    
1627    DataArrayView::RegionType slice_region=view.getSliceRegion(key);    DataArrayView::RegionType slice_region=view.getSliceRegion(key);
1628    if (slice_region.size()!=view.getRank()) {    if (slice_region.size()!=view.getRank()) {
1629      throw DataException("Error - slice size does not match Data rank.");      throw DataException("Error - slice size does not match Data rank.");
# Line 1103  void Line 1639  void
1639  Data::setSlice(const Data& value,  Data::setSlice(const Data& value,
1640                 const DataArrayView::RegionType& region)                 const DataArrayView::RegionType& region)
1641  {  {
1642    #if defined DOPROF
1643      profData->slicing++;
1644    #endif
1645    Data tempValue(value);    Data tempValue(value);
1646    typeMatchLeft(tempValue);    typeMatchLeft(tempValue);
1647    typeMatchRight(tempValue);    typeMatchRight(tempValue);
# Line 1159  Data::setTaggedValue(int tagKey, Line 1698  Data::setTaggedValue(int tagKey,
1698  }  }
1699    
1700  void  void
1701    Data::setTaggedValueFromCPP(int tagKey,
1702                                const DataArrayView& value)
1703    {
1704      //
1705      // Ensure underlying data object is of type DataTagged
1706      tag();
1707    
1708      if (!isTagged()) {
1709        throw DataException("Error - DataTagged conversion failed!!");
1710      }
1711                                                                                                                  
1712      //
1713      // Call DataAbstract::setTaggedValue
1714      m_data->setTaggedValue(tagKey,value);
1715    }
1716    
1717    int
1718    Data::getTagNumber(int dpno)
1719    {
1720      return m_data->getTagNumber(dpno);
1721    }
1722    
1723    void
1724  Data::setRefValue(int ref,  Data::setRefValue(int ref,
1725                    const boost::python::numeric::array& value)                    const boost::python::numeric::array& value)
1726  {  {
# Line 1197  Data::getRefValue(int ref, Line 1759  Data::getRefValue(int ref,
1759    DataArrayView valueView = valueDataArray.getView();    DataArrayView valueView = valueDataArray.getView();
1760    
1761    if (rank==0) {    if (rank==0) {
1762      throw DataException("Data::getRefValue error: only rank 1 data handled for now.");        boost::python::numeric::array temp_numArray(valueView());
1763          value = temp_numArray;
1764    }    }
1765    if (rank==1) {    if (rank==1) {
1766      for (int i=0; i < shape[0]; i++) {      for (int i=0; i < shape[0]; i++) {
# Line 1205  Data::getRefValue(int ref, Line 1768  Data::getRefValue(int ref,
1768      }      }
1769    }    }
1770    if (rank==2) {    if (rank==2) {
1771      throw DataException("Data::getRefValue error: only rank 1 data handled for now.");      for (int i=0; i < shape[0]; i++) {
1772          for (int j=0; j < shape[1]; j++) {
1773            value[i][j] = valueView(i,j);
1774          }
1775        }
1776    }    }
1777    if (rank==3) {    if (rank==3) {
1778      throw DataException("Data::getRefValue error: only rank 1 data handled for now.");      for (int i=0; i < shape[0]; i++) {
1779          for (int j=0; j < shape[1]; j++) {
1780            for (int k=0; k < shape[2]; k++) {
1781              value[i][j][k] = valueView(i,j,k);
1782            }
1783          }
1784        }
1785    }    }
1786    if (rank==4) {    if (rank==4) {
1787      throw DataException("Data::getRefValue error: only rank 1 data handled for now.");      for (int i=0; i < shape[0]; i++) {
1788          for (int j=0; j < shape[1]; j++) {
1789            for (int k=0; k < shape[2]; k++) {
1790              for (int l=0; l < shape[3]; l++) {
1791                value[i][j][k][l] = valueView(i,j,k,l);
1792              }
1793            }
1794          }
1795        }
1796    }    }
1797    
1798  }  }
1799    
 /*  
 Note: this version removed for now. Not needed, and breaks escript.cpp  
 void  
 Data::setTaggedValue(int tagKey,  
                      const DataArrayView& value)  
 {  
   //  
   // Ensure underlying data object is of type DataTagged  
   tag();  
   
   if (!isTagged()) {  
     throw DataException("Error - DataTagged conversion failed!!");  
   }  
                                                                                                                 
   //  
   // Call DataAbstract::setTaggedValue  
   m_data->setTaggedValue(tagKey,value);  
 }  
 */  
   
1800  void  void
1801  Data::archiveData(const std::string fileName)  Data::archiveData(const std::string fileName)
1802  {  {
# Line 1261  Data::archiveData(const std::string file Line 1822  Data::archiveData(const std::string file
1822      dataType = 3;      dataType = 3;
1823      cout << "\tdataType: DataExpanded" << endl;      cout << "\tdataType: DataExpanded" << endl;
1824    }    }
1825    
1826    if (dataType == -1) {    if (dataType == -1) {
1827      throw DataException("archiveData Error: undefined dataType");      throw DataException("archiveData Error: undefined dataType");
1828    }    }
# Line 1300  Data::archiveData(const std::string file Line 1862  Data::archiveData(const std::string file
1862    cout << ">" << endl;    cout << ">" << endl;
1863    
1864    //    //
1865    // Write common data items to archive file    // Open archive file
1866    ofstream archiveFile;    ofstream archiveFile;
1867    archiveFile.open(fileName.data(), ios::out);    archiveFile.open(fileName.data(), ios::out);
1868    
# Line 1308  Data::archiveData(const std::string file Line 1870  Data::archiveData(const std::string file
1870      throw DataException("archiveData Error: problem opening archive file");      throw DataException("archiveData Error: problem opening archive file");
1871    }    }
1872    
1873      //
1874      // Write common data items to archive file
1875    archiveFile.write(reinterpret_cast<char *>(&dataType),sizeof(int));    archiveFile.write(reinterpret_cast<char *>(&dataType),sizeof(int));
1876    archiveFile.write(reinterpret_cast<char *>(&noSamples),sizeof(int));    archiveFile.write(reinterpret_cast<char *>(&noSamples),sizeof(int));
1877    archiveFile.write(reinterpret_cast<char *>(&noDPPSample),sizeof(int));    archiveFile.write(reinterpret_cast<char *>(&noDPPSample),sizeof(int));
# Line 1331  Data::archiveData(const std::string file Line 1895  Data::archiveData(const std::string file
1895      throw DataException("archiveData Error: problem writing to archive file");      throw DataException("archiveData Error: problem writing to archive file");
1896    }    }
1897    
   archiveFile.close();  
   
   if (!archiveFile.good()) {  
     throw DataException("archiveData Error: problem closing archive file");  
   }  
   
1898    //    //
1899    // Collect and archive underlying data values for each Data type    // Archive underlying data values for each Data type
1900      int noValues;
1901    switch (dataType) {    switch (dataType) {
1902      case 0:      case 0:
1903        // DataEmpty        // DataEmpty
1904          noValues = 0;
1905          archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
1906          cout << "\tnoValues: " << noValues << endl;
1907        break;        break;
1908      case 1:      case 1:
1909        // DataConstant        // DataConstant
1910          noValues = m_data->getLength();
1911          archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
1912          cout << "\tnoValues: " << noValues << endl;
1913          if (m_data->archiveData(archiveFile,noValues)) {
1914            throw DataException("archiveData Error: problem writing data to archive file");
1915          }
1916        break;        break;
1917      case 2:      case 2:
1918        // DataTagged        // DataTagged
1919          noValues = m_data->getLength();
1920          archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
1921          cout << "\tnoValues: " << noValues << endl;
1922          if (m_data->archiveData(archiveFile,noValues)) {
1923            throw DataException("archiveData Error: problem writing data to archive file");
1924          }
1925        break;        break;
1926      case 3:      case 3:
1927        // DataExpanded        // DataExpanded
1928          noValues = m_data->getLength();
1929          archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
1930          cout << "\tnoValues: " << noValues << endl;
1931          if (m_data->archiveData(archiveFile,noValues)) {
1932            throw DataException("archiveData Error: problem writing data to archive file");
1933          }
1934        break;        break;
1935    }    }
1936    
1937      if (!archiveFile.good()) {
1938        throw DataException("archiveData Error: problem writing data to archive file");
1939      }
1940    
1941      //
1942      // Close archive file
1943      archiveFile.close();
1944    
1945      if (!archiveFile.good()) {
1946        throw DataException("archiveData Error: problem closing archive file");
1947      }
1948    
1949  }  }
1950    
1951  void  void
# Line 1379  Data::extractData(const std::string file Line 1971  Data::extractData(const std::string file
1971    int flatShape[4];    int flatShape[4];
1972    
1973    //    //
1974    // Open the archive file and read common data items    // Open the archive file
1975    ifstream archiveFile;    ifstream archiveFile;
1976    archiveFile.open(fileName.data(), ios::in);    archiveFile.open(fileName.data(), ios::in);
1977    
# Line 1387  Data::extractData(const std::string file Line 1979  Data::extractData(const std::string file
1979      throw DataException("extractData Error: problem opening archive file");      throw DataException("extractData Error: problem opening archive file");
1980    }    }
1981    
1982      //
1983      // Read common data items from archive file
1984    archiveFile.read(reinterpret_cast<char *>(&dataType),sizeof(int));    archiveFile.read(reinterpret_cast<char *>(&dataType),sizeof(int));
1985    archiveFile.read(reinterpret_cast<char *>(&noSamples),sizeof(int));    archiveFile.read(reinterpret_cast<char *>(&noSamples),sizeof(int));
1986    archiveFile.read(reinterpret_cast<char *>(&noDPPSample),sizeof(int));    archiveFile.read(reinterpret_cast<char *>(&noDPPSample),sizeof(int));
# Line 1415  Data::extractData(const std::string file Line 2009  Data::extractData(const std::string file
2009      throw DataException("extractData Error: problem reading from archive file");      throw DataException("extractData Error: problem reading from archive file");
2010    }    }
2011    
2012    archiveFile.close();    //
2013      // Verify the values just read from the archive file
   if (!archiveFile.good()) {  
     throw DataException("extractData Error: problem closing archive file");  
   }  
   
2014    switch (dataType) {    switch (dataType) {
2015      case 0:      case 0:
2016        cout << "\tdataType: DataEmpty" << endl;        cout << "\tdataType: DataEmpty" << endl;
# Line 1475  Data::extractData(const std::string file Line 2065  Data::extractData(const std::string file
2065    
2066    //    //
2067    // Load this DataVector with the appropriate values    // Load this DataVector with the appropriate values
2068      int noValues;
2069      archiveFile.read(reinterpret_cast<char *>(&noValues),sizeof(int));
2070      cout << "\tnoValues: " << noValues << endl;
2071    switch (dataType) {    switch (dataType) {
2072      case 0:      case 0:
2073        // DataEmpty        // DataEmpty
2074          if (noValues != 0) {
2075            throw DataException("extractData Error: problem reading data from archive file");
2076          }
2077        break;        break;
2078      case 1:      case 1:
2079        // DataConstant        // DataConstant
2080          if (dataVec.extractData(archiveFile,noValues)) {
2081            throw DataException("extractData Error: problem reading data from archive file");
2082          }
2083        break;        break;
2084      case 2:      case 2:
2085        // DataTagged        // DataTagged
2086          if (dataVec.extractData(archiveFile,noValues)) {
2087            throw DataException("extractData Error: problem reading data from archive file");
2088          }
2089        break;        break;
2090      case 3:      case 3:
2091        // DataExpanded        // DataExpanded
2092          if (dataVec.extractData(archiveFile,noValues)) {
2093            throw DataException("extractData Error: problem reading data from archive file");
2094          }
2095        break;        break;
2096    }    }
2097    
2098      if (!archiveFile.good()) {
2099        throw DataException("extractData Error: problem reading from archive file");
2100      }
2101    
2102      //
2103      // Close archive file
2104      archiveFile.close();
2105    
2106      if (!archiveFile.good()) {
2107        throw DataException("extractData Error: problem closing archive file");
2108      }
2109    
2110    //    //
2111    // Construct an appropriate Data object    // Construct an appropriate Data object
2112    DataAbstract* tempData;    DataAbstract* tempData;
# Line 1513  Data::extractData(const std::string file Line 2130  Data::extractData(const std::string file
2130    }    }
2131    shared_ptr<DataAbstract> temp_data(tempData);    shared_ptr<DataAbstract> temp_data(tempData);
2132    m_data=temp_data;    m_data=temp_data;
   
2133  }  }
2134    
2135  ostream& escript::operator<<(ostream& o, const Data& data)  ostream& escript::operator<<(ostream& o, const Data& data)

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

  ViewVC Help
Powered by ViewVC 1.1.26