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

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

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

revision 757 by woo409, Mon Jun 26 13:12:56 2006 UTC revision 1693 by jfenwick, Fri Aug 8 04:22:58 2008 UTC
# Line 1  Line 1 
1  // $Id$  
2  /*  /* $Id$ */
3   ************************************************************  
4   *          Copyright 2006 by ACcESS MNRF                   *  /*******************************************************
5   *                                                          *   *
6   *              http://www.access.edu.au                    *   *           Copyright 2003-2007 by ACceSS MNRF
7   *       Primary Business: Queensland, Australia            *   *       Copyright 2007 by University of Queensland
8   *  Licensed under the Open Software License version 3.0    *   *
9   *     http://www.opensource.org/licenses/osl-3.0.php       *   *                http://esscc.uq.edu.au
10   *                                                          *   *        Primary Business: Queensland, Australia
11   ************************************************************   *  Licensed under the Open Software License version 3.0
12  */   *     http://www.opensource.org/licenses/osl-3.0.php
13     *
14     *******************************************************/
15    
16  /** \file Data.h */  /** \file Data.h */
17    
# Line 26  Line 28 
28    
29  extern "C" {  extern "C" {
30  #include "DataC.h"  #include "DataC.h"
31  #include "paso/Paso.h"  /* #include "paso/Paso.h" doesn't belong in this file...causes trouble for BruceFactory.cpp */
32  }  }
33    
34    #include "esysmpi.h"
35  #include <string>  #include <string>
36  #include <algorithm>  #include <algorithm>
37    #include <sstream>
38    
39    
40  #include <boost/shared_ptr.hpp>  #include <boost/shared_ptr.hpp>
41  #include <boost/python/object.hpp>  #include <boost/python/object.hpp>
# Line 47  class DataExpanded; Line 52  class DataExpanded;
52    
53  /**  /**
54     \brief     \brief
55     Data creates the appropriate Data object for the given construction     Data creates the appropriate Data object for the given construction
56     arguments.     arguments.
57    
58     Description:     Description:
59     Data is essentially a factory class which creates the appropriate Data     Data is essentially a factory class which creates the appropriate Data
# Line 66  class Data { Line 71  class Data {
71    typedef double (*UnaryDFunPtr)(double);    typedef double (*UnaryDFunPtr)(double);
72    typedef double (*BinaryDFunPtr)(double,double);    typedef double (*BinaryDFunPtr)(double,double);
73    
74    
75    /**    /**
76       Constructors.       Constructors.
77    */    */
# Line 209  class Data { Line 215  class Data {
215       Constructor which creates a DataConstant of "shape" with constant value.       Constructor which creates a DataConstant of "shape" with constant value.
216    */    */
217    ESCRIPT_DLL_API    ESCRIPT_DLL_API
218    Data(double value,    Data(double value,
219         const boost::python::tuple& shape=boost::python::make_tuple(),         const boost::python::tuple& shape=boost::python::make_tuple(),
220         const FunctionSpace& what=FunctionSpace(),         const FunctionSpace& what=FunctionSpace(),
221         bool expanded=false);         bool expanded=false);
222    /**    /**
# Line 234  class Data { Line 240  class Data {
240    
241    /**    /**
242       \brief       \brief
243       Return the values of all data-points as a single python numarray object.       switches on update protection
244    
245    */    */
246    ESCRIPT_DLL_API    ESCRIPT_DLL_API
247    const boost::python::numeric::array    void
248    convertToNumArray();    setProtection();
249    
250    /**    /**
251       \brief       \brief
252       Return the values of all data-points for the given sample as a single python numarray object.       Returns trueif the data object is protected against update
253    
254    */    */
255    ESCRIPT_DLL_API    ESCRIPT_DLL_API
256    const boost::python::numeric::array    bool
257    convertToNumArrayFromSampleNo(int sampleNo);    isProtected() const;
258    
259    /**    /**
260       \brief       \brief
261       Return the value of the specified data-point as a single python numarray object.       Return the values of a data point on this process
262    */    */
263    ESCRIPT_DLL_API    ESCRIPT_DLL_API
264    const boost::python::numeric::array    const boost::python::numeric::array
265    convertToNumArrayFromDPNo(int sampleNo,    getValueOfDataPoint(int dataPointNo);
266                              int dataPointNo);  
267      /**
268         \brief
269         sets the values of a data-point from a python object on this process
270      */
271      ESCRIPT_DLL_API
272      void
273      setValueOfDataPointToPyObject(int dataPointNo, const boost::python::object& py_object);
274    
275      /**
276         \brief
277         sets the values of a data-point from a numarray object on this process
278      */
279      ESCRIPT_DLL_API
280      void
281      setValueOfDataPointToArray(int dataPointNo, const boost::python::numeric::array&);
282    
283    /**    /**
284       \brief       \brief
285       Fills the expanded Data object from values of a python numarray object.       sets the values of a data-point on this process
286    */    */
287    ESCRIPT_DLL_API    ESCRIPT_DLL_API
288    void    void
289    fillFromNumArray(const boost::python::numeric::array);    setValueOfDataPoint(int dataPointNo, const double);
290    
291      /**
292         \brief
293         Return the value of the specified data-point across all processors
294      */
295      ESCRIPT_DLL_API
296      const boost::python::numeric::array
297      getValueOfGlobalDataPoint(int procNo, int dataPointNo);
298    
299    /**    /**
300       \brief       \brief
# Line 294  class Data { Line 325  class Data {
325    
326    /**    /**
327       \brief       \brief
328       Write the data as a string.       Write the data as a string. For large amounts of data, a summary is printed.
329    */    */
330    ESCRIPT_DLL_API    ESCRIPT_DLL_API
   inline  
331    std::string    std::string
332    toString() const    toString() const;
333    {  
     return m_data->toString();  
   }  
334    
335    /**    /**
336       \brief       \brief
# Line 423  class Data { Line 451  class Data {
451    
452    /**    /**
453       \brief       \brief
454         Return the number of data points
455      */
456      ESCRIPT_DLL_API
457      inline
458      int
459      getNumDataPoints() const
460      {
461        return getNumSamples() * getNumDataPointsPerSample();
462      }
463      /**
464         \brief
465       Return the number of samples.       Return the number of samples.
466    */    */
467    ESCRIPT_DLL_API    ESCRIPT_DLL_API
# Line 444  class Data { Line 483  class Data {
483    {    {
484      return m_data->getNumDPPSample();      return m_data->getNumDPPSample();
485    }    }
486      /**
487         \brief
488         dumps the object into a netCDF file
489      */
490      ESCRIPT_DLL_API
491      void
492      dump(const std::string fileName) const;
493    /**    /**
494       \brief       \brief
495       Return the sample data for the given sample no. This is not the       Return the sample data for the given sample no. This is not the
# Line 475  class Data { Line 520  class Data {
520    
521    /**    /**
522       \brief       \brief
      Assign the given value to the data-points referenced by the given  
      reference number.  
   
      The value supplied is a python numarray object.  The data from this numarray  
      is unpacked into a DataArray, and this is used to set the corresponding  
      data-points in the underlying Data object.  
   
      If the underlying Data object cannot be accessed via reference numbers, an  
      exception will be thrown.  
   
      \param ref - Input - reference number.  
      \param value - Input - value to assign to data-points associated with  
                             the given reference number.  
   */  
   ESCRIPT_DLL_API  
   void  
   setRefValue(int ref,  
               const boost::python::numeric::array& value);  
   
   /**  
      \brief  
      Return the values associated with the data-points referenced by the given  
      reference number.  
   
      The value supplied is a python numarray object. The data from the corresponding  
      data-points in this Data object are packed into the given numarray object.  
   
      If the underlying Data object cannot be accessed via reference numbers, an  
      exception will be thrown.  
   
      \param ref - Input - reference number.  
      \param value - Output - object to receive values from data-points  
                              associated with the given reference number.  
   */  
   ESCRIPT_DLL_API  
   void  
   getRefValue(int ref,  
               boost::python::numeric::array& value);  
   
   /**  
      \brief  
523       Return a view into the data for the data point specified.       Return a view into the data for the data point specified.
524       NOTE: Construction of the DataArrayView is a relatively expensive       NOTE: Construction of the DataArrayView is a relatively expensive
525       operation.       operation.
# Line 528  class Data { Line 532  class Data {
532    getDataPoint(int sampleNo,    getDataPoint(int sampleNo,
533                 int dataPointNo)                 int dataPointNo)
534    {    {
535      return m_data->getDataPoint(sampleNo,dataPointNo);                  return m_data->getDataPoint(sampleNo,dataPointNo);
536    }    }
537    
538    /**    /**
# Line 564  class Data { Line 568  class Data {
568    DataArrayView::ValueType::size_type    DataArrayView::ValueType::size_type
569    getLength() const;    getLength() const;
570    
571    
572    
573    /**    /**
574       \brief       \brief
575       Assign the given value to the tag. Implicitly converts this       Assign the given value to the tag assocciated with name. Implicitly converts this
576       object to type DataTagged. Throws an exception if this object       object to type DataTagged. Throws an exception if this object
577       cannot be converted to a DataTagged object.       cannot be converted to a DataTagged object or name cannot be mapped onto a tag key.
578         \param tagKey - Input - Integer key.
579         \param value - Input - Value to associate with given key.
580        ==>*
581      */
582      ESCRIPT_DLL_API
583      void
584      setTaggedValueByName(std::string name,
585                           const boost::python::object& value);
586    
587      /**
588         \brief
589         Assign the given value to the tag. Implicitly converts this
590         object to type DataTagged if it is constant.
591    
592       \param tagKey - Input - Integer key.       \param tagKey - Input - Integer key.
593       \param value - Input - Value to associate with given key.       \param value - Input - Value to associate with given key.
594      ==>*      ==>*
# Line 581  class Data { Line 601  class Data {
601    /**    /**
602       \brief       \brief
603       Assign the given value to the tag. Implicitly converts this       Assign the given value to the tag. Implicitly converts this
604       object to type DataTagged. Throws an exception if this object       object to type DataTagged if it is constant.
605       cannot be converted to a DataTagged object.  
606       \param tagKey - Input - Integer key.       \param tagKey - Input - Integer key.
607       \param value - Input - Value to associate with given key.       \param value - Input - Value to associate with given key.
608      ==>*      ==>*
# Line 607  class Data { Line 627  class Data {
627    
628    /**    /**
629       \brief       \brief
630         set all values to zero
631         *
632      */
633      ESCRIPT_DLL_API
634      void
635      setToZero();
636    
637      /**
638         \brief
639       Interpolates this onto the given functionspace and returns       Interpolates this onto the given functionspace and returns
640       the result as a Data object.       the result as a Data object.
641       *       *
# Line 614  class Data { Line 643  class Data {
643    ESCRIPT_DLL_API    ESCRIPT_DLL_API
644    Data    Data
645    interpolate(const FunctionSpace& functionspace) const;    interpolate(const FunctionSpace& functionspace) const;
   
646    /**    /**
647       \brief       \brief
648       Calculates the gradient of the data at the data points of functionspace.       Calculates the gradient of the data at the data points of functionspace.
# Line 640  class Data { Line 668  class Data {
668    
669    /**    /**
670       \brief       \brief
671         Returns 1./ Data object
672         *
673      */
674      ESCRIPT_DLL_API
675      Data
676      oneOver() const;
677      /**
678         \brief
679       Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.       Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.
680       *       *
681    */    */
# Line 703  class Data { Line 739  class Data {
739    
740    /**    /**
741       \brief       \brief
      Return the minimum absolute value of this Data object.  
      *  
   */  
   ESCRIPT_DLL_API  
   double  
   Linf() const;  
   
   /**  
      \brief  
742       Return the maximum value of this Data object.       Return the maximum value of this Data object.
743       *       *
744    */    */
# Line 762  class Data { Line 789  class Data {
789    */    */
790    ESCRIPT_DLL_API    ESCRIPT_DLL_API
791    const boost::python::tuple    const boost::python::tuple
792    mindp() const;    minGlobalDataPoint() const;
793    
794    ESCRIPT_DLL_API    ESCRIPT_DLL_API
795    void    void
796    calc_mindp(int& SampleNo,    calc_minGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;
              int& DataPointNo) const;  
   
797    /**    /**
798       \brief       \brief
799       Return the sign of each data point of this Data object.       Return the sign of each data point of this Data object.
# Line 781  class Data { Line 806  class Data {
806    
807    /**    /**
808       \brief       \brief
809         Return the symmetric part of a matrix which is half the matrix plus its transpose.
810         *
811      */
812      ESCRIPT_DLL_API
813      Data
814      symmetric() const;
815    
816      /**
817         \brief
818         Return the nonsymmetric part of a matrix which is half the matrix minus its transpose.
819         *
820      */
821      ESCRIPT_DLL_API
822      Data
823      nonsymmetric() const;
824    
825      /**
826         \brief
827         Return the trace of a matrix
828         *
829      */
830      ESCRIPT_DLL_API
831      Data
832      trace(int axis_offset) const;
833    
834      /**
835         \brief
836         Transpose each data point of this Data object around the given axis.
837         *
838      */
839      ESCRIPT_DLL_API
840      Data
841      transpose(int axis_offset) const;
842    
843      /**
844         \brief
845       Return the eigenvalues of the symmetric part at each data point of this Data object in increasing values.       Return the eigenvalues of the symmetric part at each data point of this Data object in increasing values.
846       Currently this function is restricted to rank 2, square shape, and dimension 3.       Currently this function is restricted to rank 2, square shape, and dimension 3.
847       *       *
# Line 792  class Data { Line 853  class Data {
853    /**    /**
854       \brief       \brief
855       Return the eigenvalues and corresponding eigenvcetors of the symmetric part at each data point of this Data object.       Return the eigenvalues and corresponding eigenvcetors of the symmetric part at each data point of this Data object.
856       the eigenvalues are ordered in increasing size where eigenvalues with relative difference less than       the eigenvalues are ordered in increasing size where eigenvalues with relative difference less than
857       tol are treated as equal. The eigenvectors are orthogonal, normalized and the sclaed such that the       tol are treated as equal. The eigenvectors are orthogonal, normalized and the sclaed such that the
858       first non-zero entry is positive.       first non-zero entry is positive.
859       Currently this function is restricted to rank 2, square shape, and dimension 3       Currently this function is restricted to rank 2, square shape, and dimension 3
860       *       *
# Line 804  class Data { Line 865  class Data {
865    
866    /**    /**
867       \brief       \brief
868       Transpose each data point of this Data object around the given axis.       swaps the components axis0 and axis1
      --* not implemented yet *--  
869       *       *
870    */    */
871    ESCRIPT_DLL_API    ESCRIPT_DLL_API
872    Data    Data
873    transpose(int axis) const;    swapaxes(const int axis0, const int axis1) const;
874    
875    /**    /**
876       \brief       \brief
877       Calculate the trace of each data point of this Data object.       Return the error function erf of each data point of this Data object.
878       *       *
879    */    */
880    ESCRIPT_DLL_API    ESCRIPT_DLL_API
881    Data    Data
882    trace() const;    erf() const;
883    
884    /**    /**
885       \brief       \brief
# Line 998  class Data { Line 1058  class Data {
1058    /**    /**
1059       \brief       \brief
1060       Return the given power of each data point of this boost python object.       Return the given power of each data point of this boost python object.
1061        
1062       \param right Input - the power to raise the object to.       \param right Input - the power to raise the object to.
1063       *       *
1064     */     */
# Line 1009  class Data { Line 1069  class Data {
1069    /**    /**
1070       \brief       \brief
1071       Return the given power of each data point of this boost python object.       Return the given power of each data point of this boost python object.
1072        
1073       \param left Input - the bases       \param left Input - the bases
1074       *       *
1075     */     */
# Line 1045  class Data { Line 1105  class Data {
1105    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1106    Data& operator+=(const boost::python::object& right);    Data& operator+=(const boost::python::object& right);
1107    
1108      ESCRIPT_DLL_API
1109      Data& operator=(const Data& other);
1110    
1111    /**    /**
1112       \brief       \brief
1113       Overloaded operator -=       Overloaded operator -=
# Line 1137  class Data { Line 1200  class Data {
1200    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1201    inline    inline
1202    void    void
1203    unaryOp(UnaryFunction operation);    unaryOp2(UnaryFunction operation);
1204    
1205    /**    /**
1206       \brief       \brief
# Line 1188  class Data { Line 1251  class Data {
1251    
1252    /**    /**
1253       \brief       \brief
1254       print the data values to stdout. Used for debugging       print the data values to stdout. Used for debugging
1255    */    */
1256    void print();    ESCRIPT_DLL_API
1257      void
1258            print(void);
1259    
1260      /**
1261         \brief
1262         return the MPI rank number of the local data
1263                     MPI_COMM_WORLD is assumed and the result of MPI_Comm_size()
1264                     is returned
1265      */
1266      ESCRIPT_DLL_API
1267            int
1268            get_MPIRank(void) const;
1269    
1270      /**
1271         \brief
1272         return the MPI rank number of the local data
1273                     MPI_COMM_WORLD is assumed and the result of MPI_Comm_rank()
1274                     is returned
1275      */
1276      ESCRIPT_DLL_API
1277            int
1278            get_MPISize(void) const;
1279    
1280      /**
1281         \brief
1282         return the MPI rank number of the local data
1283                     MPI_COMM_WORLD is assumed and returned.
1284      */
1285      ESCRIPT_DLL_API
1286            MPI_Comm
1287            get_MPIComm(void) const;
1288    
1289      /**
1290         \brief
1291         return the object produced by the factory, which is a DataConstant or DataExpanded
1292      */
1293      ESCRIPT_DLL_API
1294            DataAbstract*
1295            borrowData(void) const;
1296    
1297   protected:   protected:
1298    
# Line 1249  class Data { Line 1351  class Data {
1351    
1352    /**    /**
1353       \brief       \brief
      Perform the given binary operation on all of the data's elements.  
      RHS is a boost::python object.  
   */  
   template <class BinaryFunction>  
   inline  
   void  
   binaryOp(const boost::python::object& right,  
            BinaryFunction operation);  
   
   /**  
      \brief  
1354       Convert the data type of the RHS to match this.       Convert the data type of the RHS to match this.
1355       \param right - Input - data type to match.       \param right - Input - data type to match.
1356    */    */
# Line 1284  class Data { Line 1375  class Data {
1375               const FunctionSpace& what,               const FunctionSpace& what,
1376               bool expanded);               bool expanded);
1377    
1378    /**    //
1379       \brief    // flag to protect the data object against any update
1380       Reshape the data point if the data point is currently rank 0.    bool m_protected;
      Will throw an exception if the data points are not rank 0.  
      The original data point value is used for all values of the new  
      data point.  
   */  
   void  
   reshapeDataPoint(const DataArrayView::ShapeType& shape);  
1381    
1382    //    //
1383    // pointer to the actual data object    // pointer to the actual data object
1384    boost::shared_ptr<DataAbstract> m_data;    boost::shared_ptr<DataAbstract> m_data;
1385    
   //  
   // pointer to the internal profiling data  
   struct profDataEntry *profData;  
   
1386  };  };
1387    
1388  template <class IValueType>  template <class IValueType>
# Line 1329  Data::initialise(const IValueType& value Line 1410  Data::initialise(const IValueType& value
1410  /**  /**
1411     Binary Data object operators.     Binary Data object operators.
1412  */  */
1413    inline double rpow(double x,double y)
1414    {
1415        return pow(y,x);
1416    }
1417    
1418  /**  /**
1419    \brief    \brief
# Line 1422  ESCRIPT_DLL_API Data operator*(const boo Line 1507  ESCRIPT_DLL_API Data operator*(const boo
1507  */  */
1508  ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);  ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);
1509    
1510    
1511    
1512  /**  /**
1513    \brief    \brief
1514    Output operator    Output operator
# Line 1430  ESCRIPT_DLL_API std::ostream& operator<< Line 1517  ESCRIPT_DLL_API std::ostream& operator<<
1517    
1518  /**  /**
1519    \brief    \brief
1520      Compute a tensor product of two Data objects
1521      \param arg0 - Input - Data object
1522      \param arg1 - Input - Data object
1523      \param axis_offset - Input - axis offset
1524      \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1
1525    */
1526    ESCRIPT_DLL_API
1527    Data
1528    C_GeneralTensorProduct(Data& arg0,
1529                         Data& arg1,
1530                         int axis_offset=0,
1531                         int transpose=0);
1532    
1533    
1534    
1535    /**
1536      \brief
1537    Return true if operands are equivalent, else return false.    Return true if operands are equivalent, else return false.
1538    NB: this operator does very little at this point, and isn't to    NB: this operator does very little at this point, and isn't to
1539    be relied on. Requires further implementation.    be relied on. Requires further implementation.
1540  */  */
1541  //ESCRIPT_DLL_API bool operator==(const Data& left, const Data& right);  // ESCRIPT_DLL_API bool operator==(const Data& left, const Data& right);
1542    
1543  /**  /**
1544    \brief    \brief
# Line 1450  Data::binaryOp(const Data& right, Line 1554  Data::binaryOp(const Data& right,
1554     //     //
1555     // if this has a rank of zero promote it to the rank of the RHS     // if this has a rank of zero promote it to the rank of the RHS
1556     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {
1557       reshapeDataPoint(right.getPointDataView().getShape());       throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");
1558     }     }
1559     //     //
1560     // initially make the temporary a shallow copy     // initially make the temporary a shallow copy
# Line 1458  Data::binaryOp(const Data& right, Line 1562  Data::binaryOp(const Data& right,
1562     if (getFunctionSpace()!=right.getFunctionSpace()) {     if (getFunctionSpace()!=right.getFunctionSpace()) {
1563       if (right.probeInterpolation(getFunctionSpace())) {       if (right.probeInterpolation(getFunctionSpace())) {
1564         //         //
1565         // an interpolation is required so create a new Data         // an interpolation is required so create a new Data
1566         tempRight=Data(right,this->getFunctionSpace());         tempRight=Data(right,this->getFunctionSpace());
1567       } else if (probeInterpolation(right.getFunctionSpace())) {       } else if (probeInterpolation(right.getFunctionSpace())) {
1568         //         //
# Line 1506  Data::binaryOp(const Data& right, Line 1610  Data::binaryOp(const Data& right,
1610    
1611  /**  /**
1612    \brief    \brief
   Perform the given binary operation with this and right as operands.  
   Right is a boost::python object.  
 */  
 template <class BinaryFunction>  
 inline  
 void  
 Data::binaryOp(const boost::python::object& right,  
                BinaryFunction operation)  
 {  
    DataArray temp(right);  
    //  
    // if this has a rank of zero promote it to the rank of the RHS.  
    if (getPointDataView().getRank()==0 && temp.getView().getRank()!=0) {  
       reshapeDataPoint(temp.getView().getShape());  
    }  
    //  
    // Always allow scalar values for the RHS but check other shapes  
    if (temp.getView().getRank()!=0) {  
      if (!getPointDataView().checkShape(temp.getView().getShape())) {  
        throw DataException(getPointDataView().createShapeErrorMessage(  
                   "Error - RHS shape doesn't match LHS shape.",temp.getView().getShape()));  
      }  
    }  
    if (isExpanded()) {  
      DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());  
      EsysAssert((leftC!=0),"Programming error - casting to DataExpanded.");  
      escript::binaryOp(*leftC,temp.getView(),operation);  
    } else if (isTagged()) {  
      DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());  
      EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");  
      escript::binaryOp(*leftC,temp.getView(),operation);  
    } else if (isConstant()) {  
      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());  
      EsysAssert((leftC!=0),"Programming error - casting to DataConstant.");  
      escript::binaryOp(*leftC,temp.getView(),operation);  
    }  
 }  
   
 /**  
   \brief  
   Perform the given unary operation on other and return the result.  
   Given operation is performed on each element of each data point, thus  
   argument object is a rank n Data object, and returned object is a rank n  
   Data object.  
   Calls Data::unaryOp.  
 */  
 template <class UnaryFunction>  
 inline  
 Data  
 unaryOp(const Data& other,  
         UnaryFunction operation)  
 {  
   Data result;  
   result.copy(other);  
   result.unaryOp(operation);  
   return result;  
 }  
   
 /**  
   \brief  
   Perform the given unary operation on this.  
   Given operation is performed on each element of each data point.  
   Calls escript::unaryOp.  
 */  
 template <class UnaryFunction>  
 inline  
 void  
 Data::unaryOp(UnaryFunction operation)  
 {  
   if (isExpanded()) {  
     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());  
     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");  
     escript::unaryOp(*leftC,operation);  
   } else if (isTagged()) {  
     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());  
     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");  
     escript::unaryOp(*leftC,operation);  
   } else if (isConstant()) {  
     DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());  
     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");  
     escript::unaryOp(*leftC,operation);  
   }  
 }  
   
 /**  
   \brief  
1613    Perform the given Data object reduction algorithm on this and return the result.    Perform the given Data object reduction algorithm on this and return the result.
1614    Given operation combines each element of each data point, thus argument    Given operation combines each element of each data point, thus argument
1615    object (*this) is a rank n Data object, and returned object is a scalar.    object (*this) is a rank n Data object, and returned object is a scalar.
# Line 1622  Data::algorithm(BinaryFunction operation Line 1640  Data::algorithm(BinaryFunction operation
1640    \brief    \brief
1641    Perform the given data point reduction algorithm on data and return the result.    Perform the given data point reduction algorithm on data and return the result.
1642    Given operation combines each element within each data point into a scalar,    Given operation combines each element within each data point into a scalar,
1643    thus argument object is a rank n Data object, and returned object is a    thus argument object is a rank n Data object, and returned object is a
1644    rank 0 Data object.    rank 0 Data object.
1645    Calls escript::dp_algorithm.    Calls escript::dp_algorithm.
1646  */  */
# Line 1671  Data::dp_algorithm(BinaryFunction operat Line 1689  Data::dp_algorithm(BinaryFunction operat
1689    return falseRetVal;    return falseRetVal;
1690  }  }
1691    
1692    /**
1693      \brief
1694      Compute a tensor operation with two Data objects
1695      \param arg0 - Input - Data object
1696      \param arg1 - Input - Data object
1697      \param operation - Input - Binary op functor
1698    */
1699    template <typename BinaryFunction>
1700    inline
1701    Data
1702    C_TensorBinaryOperation(Data const &arg_0,
1703                            Data const &arg_1,
1704                            BinaryFunction operation)
1705    {
1706      // Interpolate if necessary and find an appropriate function space
1707      Data arg_0_Z, arg_1_Z;
1708      if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {
1709        if (arg_0.probeInterpolation(arg_1.getFunctionSpace())) {
1710          arg_0_Z = arg_0.interpolate(arg_1.getFunctionSpace());
1711          arg_1_Z = Data(arg_1);
1712        }
1713        else if (arg_1.probeInterpolation(arg_0.getFunctionSpace())) {
1714          arg_1_Z=arg_1.interpolate(arg_0.getFunctionSpace());
1715          arg_0_Z =Data(arg_0);
1716        }
1717        else {
1718          throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible function spaces.");
1719        }
1720      } else {
1721          arg_0_Z = Data(arg_0);
1722          arg_1_Z = Data(arg_1);
1723      }
1724      // Get rank and shape of inputs
1725      int rank0 = arg_0_Z.getDataPointRank();
1726      int rank1 = arg_1_Z.getDataPointRank();
1727      DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();
1728      DataArrayView::ShapeType shape1 = arg_1_Z.getDataPointShape();
1729      int size0 = arg_0_Z.getDataPointSize();
1730      int size1 = arg_1_Z.getDataPointSize();
1731    
1732      // Declare output Data object
1733      Data res;
1734    
1735      if (shape0 == shape1) {
1736    
1737        if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
1738          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
1739          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
1740          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
1741          double *ptr_2 = &((res.getPointDataView().getData())[0]);
1742          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1743        }
1744        else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
1745    
1746          // Prepare the DataConstant input
1747          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
1748    
1749          // Borrow DataTagged input from Data object
1750          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1751    
1752          // Prepare a DataTagged output 2
1753          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataTagged output
1754          res.tag();
1755          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1756    
1757          // Prepare offset into DataConstant
1758          int offset_0 = tmp_0->getPointOffset(0,0);
1759          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1760          // Get the views
1761          DataArrayView view_1 = tmp_1->getDefaultValue();
1762          DataArrayView view_2 = tmp_2->getDefaultValue();
1763          // Get the pointers to the actual data
1764          double *ptr_1 = &((view_1.getData())[0]);
1765          double *ptr_2 = &((view_2.getData())[0]);
1766          // Compute a result for the default
1767          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1768          // Compute a result for each tag
1769          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
1770          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1771          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
1772            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1773            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
1774            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1775            double *ptr_1 = &view_1.getData(0);
1776            double *ptr_2 = &view_2.getData(0);
1777            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1778          }
1779    
1780        }
1781        else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
1782    
1783          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1784          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
1785          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1786          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1787    
1788          int sampleNo_1,dataPointNo_1;
1789          int numSamples_1 = arg_1_Z.getNumSamples();
1790          int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
1791          int offset_0 = tmp_0->getPointOffset(0,0);
1792          #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
1793          for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
1794            for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
1795              int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
1796              int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
1797              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1798              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1799              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1800              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1801            }
1802          }
1803    
1804        }
1805        else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
1806    
1807          // Borrow DataTagged input from Data object
1808          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1809    
1810          // Prepare the DataConstant input
1811          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
1812    
1813          // Prepare a DataTagged output 2
1814          res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataTagged output
1815          res.tag();
1816          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1817    
1818          // Prepare offset into DataConstant
1819          int offset_1 = tmp_1->getPointOffset(0,0);
1820          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1821          // Get the views
1822          DataArrayView view_0 = tmp_0->getDefaultValue();
1823          DataArrayView view_2 = tmp_2->getDefaultValue();
1824          // Get the pointers to the actual data
1825          double *ptr_0 = &((view_0.getData())[0]);
1826          double *ptr_2 = &((view_2.getData())[0]);
1827          // Compute a result for the default
1828          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1829          // Compute a result for each tag
1830          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
1831          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1832          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
1833            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1834            DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
1835            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1836            double *ptr_0 = &view_0.getData(0);
1837            double *ptr_2 = &view_2.getData(0);
1838            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1839          }
1840    
1841        }
1842        else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
1843    
1844          // Borrow DataTagged input from Data object
1845          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1846    
1847          // Borrow DataTagged input from Data object
1848          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1849    
1850          // Prepare a DataTagged output 2
1851          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
1852          res.tag();        // DataTagged output
1853          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1854    
1855          // Get the views
1856          DataArrayView view_0 = tmp_0->getDefaultValue();
1857          DataArrayView view_1 = tmp_1->getDefaultValue();
1858          DataArrayView view_2 = tmp_2->getDefaultValue();
1859          // Get the pointers to the actual data
1860          double *ptr_0 = &((view_0.getData())[0]);
1861          double *ptr_1 = &((view_1.getData())[0]);
1862          double *ptr_2 = &((view_2.getData())[0]);
1863          // Compute a result for the default
1864          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1865          // Merge the tags
1866          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1867          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
1868          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
1869          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
1870            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
1871          }
1872          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
1873            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1874          }
1875          // Compute a result for each tag
1876          const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
1877          for (i=lookup_2.begin();i!=lookup_2.end();i++) {
1878            DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
1879            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
1880            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1881            double *ptr_0 = &view_0.getData(0);
1882            double *ptr_1 = &view_1.getData(0);
1883            double *ptr_2 = &view_2.getData(0);
1884            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1885          }
1886    
1887        }
1888        else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
1889    
1890          // After finding a common function space above the two inputs have the same numSamples and num DPPS
1891          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1892          DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1893          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1894          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1895    
1896          int sampleNo_0,dataPointNo_0;
1897          int numSamples_0 = arg_0_Z.getNumSamples();
1898          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1899          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1900          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1901            int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
1902            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1903            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1904              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
1905              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1906              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1907              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1908              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1909            }
1910          }
1911    
1912        }
1913        else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
1914    
1915          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1916          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1917          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
1918          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1919    
1920          int sampleNo_0,dataPointNo_0;
1921          int numSamples_0 = arg_0_Z.getNumSamples();
1922          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1923          int offset_1 = tmp_1->getPointOffset(0,0);
1924          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1925          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1926            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1927              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1928              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1929              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1930              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1931              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1932              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1933            }
1934          }
1935    
1936        }
1937        else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
1938    
1939          // After finding a common function space above the two inputs have the same numSamples and num DPPS
1940          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1941          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1942          DataTagged*   tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1943          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1944    
1945          int sampleNo_0,dataPointNo_0;
1946          int numSamples_0 = arg_0_Z.getNumSamples();
1947          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1948          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1949          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1950            int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
1951            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1952            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1953              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1954              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1955              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1956              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1957              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1958            }
1959          }
1960    
1961        }
1962        else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
1963    
1964          // After finding a common function space above the two inputs have the same numSamples and num DPPS
1965          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1966          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1967          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1968          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1969    
1970          int sampleNo_0,dataPointNo_0;
1971          int numSamples_0 = arg_0_Z.getNumSamples();
1972          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1973          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1974          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1975            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1976              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1977              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
1978              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1979              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1980              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1981              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1982              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1983            }
1984          }
1985    
1986        }
1987        else {
1988          throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
1989        }
1990    
1991      } else if (0 == rank0) {
1992    
1993        if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
1994          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataConstant output
1995          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
1996          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
1997          double *ptr_2 = &((res.getPointDataView().getData())[0]);
1998          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
1999        }
2000        else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
2001    
2002          // Prepare the DataConstant input
2003          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2004    
2005          // Borrow DataTagged input from Data object
2006          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2007    
2008          // Prepare a DataTagged output 2
2009          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataTagged output
2010          res.tag();
2011          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2012    
2013          // Prepare offset into DataConstant
2014          int offset_0 = tmp_0->getPointOffset(0,0);
2015          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2016          // Get the views
2017          DataArrayView view_1 = tmp_1->getDefaultValue();
2018          DataArrayView view_2 = tmp_2->getDefaultValue();
2019          // Get the pointers to the actual data
2020          double *ptr_1 = &((view_1.getData())[0]);
2021          double *ptr_2 = &((view_2.getData())[0]);
2022          // Compute a result for the default
2023          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2024          // Compute a result for each tag
2025          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2026          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2027          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2028            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2029            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2030            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2031            double *ptr_1 = &view_1.getData(0);
2032            double *ptr_2 = &view_2.getData(0);
2033            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2034          }
2035    
2036        }
2037        else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2038    
2039          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2040          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2041          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2042          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2043    
2044          int sampleNo_1,dataPointNo_1;
2045          int numSamples_1 = arg_1_Z.getNumSamples();
2046          int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2047          int offset_0 = tmp_0->getPointOffset(0,0);
2048          #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2049          for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2050            for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2051              int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2052              int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2053              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2054              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2055              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2056              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2057    
2058            }
2059          }
2060    
2061        }
2062        else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
2063    
2064          // Borrow DataTagged input from Data object
2065          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2066    
2067          // Prepare the DataConstant input
2068          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2069    
2070          // Prepare a DataTagged output 2
2071          res = Data(0.0, shape1, arg_0_Z.getFunctionSpace());      // DataTagged output
2072          res.tag();
2073          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2074    
2075          // Prepare offset into DataConstant
2076          int offset_1 = tmp_1->getPointOffset(0,0);
2077          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2078          // Get the views
2079          DataArrayView view_0 = tmp_0->getDefaultValue();
2080          DataArrayView view_2 = tmp_2->getDefaultValue();
2081          // Get the pointers to the actual data
2082          double *ptr_0 = &((view_0.getData())[0]);
2083          double *ptr_2 = &((view_2.getData())[0]);
2084          // Compute a result for the default
2085          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2086          // Compute a result for each tag
2087          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2088          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2089          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2090            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2091            DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2092            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2093            double *ptr_0 = &view_0.getData(0);
2094            double *ptr_2 = &view_2.getData(0);
2095            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2096          }
2097    
2098        }
2099        else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
2100    
2101          // Borrow DataTagged input from Data object
2102          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2103    
2104          // Borrow DataTagged input from Data object
2105          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2106    
2107          // Prepare a DataTagged output 2
2108          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());
2109          res.tag();        // DataTagged output
2110          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2111    
2112          // Get the views
2113          DataArrayView view_0 = tmp_0->getDefaultValue();
2114          DataArrayView view_1 = tmp_1->getDefaultValue();
2115          DataArrayView view_2 = tmp_2->getDefaultValue();
2116          // Get the pointers to the actual data
2117          double *ptr_0 = &((view_0.getData())[0]);
2118          double *ptr_1 = &((view_1.getData())[0]);
2119          double *ptr_2 = &((view_2.getData())[0]);
2120          // Compute a result for the default
2121          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2122          // Merge the tags
2123          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2124          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2125          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2126          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2127            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
2128          }
2129          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2130            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2131          }
2132          // Compute a result for each tag
2133          const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2134          for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2135            DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2136            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2137            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2138            double *ptr_0 = &view_0.getData(0);
2139            double *ptr_1 = &view_1.getData(0);
2140            double *ptr_2 = &view_2.getData(0);
2141            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2142          }
2143    
2144        }
2145        else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
2146    
2147          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2148          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2149          DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2150          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2151          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2152    
2153          int sampleNo_0,dataPointNo_0;
2154          int numSamples_0 = arg_0_Z.getNumSamples();
2155          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2156          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2157          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2158            int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2159            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2160            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2161              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2162              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2163              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2164              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2165              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2166            }
2167          }
2168    
2169        }
2170        else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
2171    
2172          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2173          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2174          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2175          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2176    
2177          int sampleNo_0,dataPointNo_0;
2178          int numSamples_0 = arg_0_Z.getNumSamples();
2179          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2180          int offset_1 = tmp_1->getPointOffset(0,0);
2181          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2182          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2183            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2184              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2185              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2186              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2187              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2188              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2189              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2190            }
2191          }
2192    
2193    
2194        }
2195        else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
2196    
2197          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2198          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2199          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2200          DataTagged*   tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2201          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2202    
2203          int sampleNo_0,dataPointNo_0;
2204          int numSamples_0 = arg_0_Z.getNumSamples();
2205          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2206          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2207          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2208            int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2209            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2210            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2211              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2212              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2213              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2214              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2215              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2216            }
2217          }
2218    
2219        }
2220        else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
2221    
2222          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2223          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2224          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2225          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2226          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2227    
2228          int sampleNo_0,dataPointNo_0;
2229          int numSamples_0 = arg_0_Z.getNumSamples();
2230          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2231          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2232          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2233            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2234              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2235              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2236              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2237              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2238              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2239              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2240              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2241            }
2242          }
2243    
2244        }
2245        else {
2246          throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2247        }
2248    
2249      } else if (0 == rank1) {
2250    
2251        if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2252          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2253          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2254          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
2255          double *ptr_2 = &((res.getPointDataView().getData())[0]);
2256          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2257        }
2258        else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
2259    
2260          // Prepare the DataConstant input
2261          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2262    
2263          // Borrow DataTagged input from Data object
2264          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2265    
2266          // Prepare a DataTagged output 2
2267          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataTagged output
2268          res.tag();
2269          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2270    
2271          // Prepare offset into DataConstant
2272          int offset_0 = tmp_0->getPointOffset(0,0);
2273          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2274          // Get the views
2275          DataArrayView view_1 = tmp_1->getDefaultValue();
2276          DataArrayView view_2 = tmp_2->getDefaultValue();
2277          // Get the pointers to the actual data
2278          double *ptr_1 = &((view_1.getData())[0]);
2279          double *ptr_2 = &((view_2.getData())[0]);
2280          // Compute a result for the default
2281          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2282          // Compute a result for each tag
2283          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2284          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2285          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2286            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2287            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2288            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2289            double *ptr_1 = &view_1.getData(0);
2290            double *ptr_2 = &view_2.getData(0);
2291            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2292          }
2293    
2294        }
2295        else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2296    
2297          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2298          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2299          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2300          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2301    
2302          int sampleNo_1,dataPointNo_1;
2303          int numSamples_1 = arg_1_Z.getNumSamples();
2304          int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2305          int offset_0 = tmp_0->getPointOffset(0,0);
2306          #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2307          for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2308            for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2309              int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2310              int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2311              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2312              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2313              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2314              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2315            }
2316          }
2317    
2318        }
2319        else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
2320    
2321          // Borrow DataTagged input from Data object
2322          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2323    
2324          // Prepare the DataConstant input
2325          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2326    
2327          // Prepare a DataTagged output 2
2328          res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataTagged output
2329          res.tag();
2330          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2331    
2332          // Prepare offset into DataConstant
2333          int offset_1 = tmp_1->getPointOffset(0,0);
2334          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2335          // Get the views
2336          DataArrayView view_0 = tmp_0->getDefaultValue();
2337          DataArrayView view_2 = tmp_2->getDefaultValue();
2338          // Get the pointers to the actual data
2339          double *ptr_0 = &((view_0.getData())[0]);
2340          double *ptr_2 = &((view_2.getData())[0]);
2341          // Compute a result for the default
2342          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2343          // Compute a result for each tag
2344          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2345          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2346          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2347            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2348            DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2349            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2350            double *ptr_0 = &view_0.getData(0);
2351            double *ptr_2 = &view_2.getData(0);
2352            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2353          }
2354    
2355        }
2356        else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
2357    
2358          // Borrow DataTagged input from Data object
2359          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2360    
2361          // Borrow DataTagged input from Data object
2362          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2363    
2364          // Prepare a DataTagged output 2
2365          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2366          res.tag();        // DataTagged output
2367          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2368    
2369          // Get the views
2370          DataArrayView view_0 = tmp_0->getDefaultValue();
2371          DataArrayView view_1 = tmp_1->getDefaultValue();
2372          DataArrayView view_2 = tmp_2->getDefaultValue();
2373          // Get the pointers to the actual data
2374          double *ptr_0 = &((view_0.getData())[0]);
2375          double *ptr_1 = &((view_1.getData())[0]);
2376          double *ptr_2 = &((view_2.getData())[0]);
2377          // Compute a result for the default
2378          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2379          // Merge the tags
2380          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2381          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2382          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2383          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2384            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
2385          }
2386          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2387            tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2388          }
2389          // Compute a result for each tag
2390          const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2391          for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2392            DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2393            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2394            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2395            double *ptr_0 = &view_0.getData(0);
2396            double *ptr_1 = &view_1.getData(0);
2397            double *ptr_2 = &view_2.getData(0);
2398            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2399          }
2400    
2401        }
2402        else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
2403    
2404          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2405          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2406          DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2407          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2408          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2409    
2410          int sampleNo_0,dataPointNo_0;
2411          int numSamples_0 = arg_0_Z.getNumSamples();
2412          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2413          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2414          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2415            int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2416            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2417            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2418              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2419              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2420              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2421              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2422              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2423            }
2424          }
2425    
2426        }
2427        else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
2428    
2429          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2430          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2431          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2432          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2433    
2434          int sampleNo_0,dataPointNo_0;
2435          int numSamples_0 = arg_0_Z.getNumSamples();
2436          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2437          int offset_1 = tmp_1->getPointOffset(0,0);
2438          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2439          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2440            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2441              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2442              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2443              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2444              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2445              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2446              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2447            }
2448          }
2449    
2450    
2451        }
2452        else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
2453    
2454          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2455          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2456          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2457          DataTagged*   tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2458          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2459    
2460          int sampleNo_0,dataPointNo_0;
2461          int numSamples_0 = arg_0_Z.getNumSamples();
2462          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2463          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2464          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2465            int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2466            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2467            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2468              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2469              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2470              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2471              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2472              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2473            }
2474          }
2475    
2476        }
2477        else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
2478    
2479          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2480          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2481          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2482          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2483          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2484    
2485          int sampleNo_0,dataPointNo_0;
2486          int numSamples_0 = arg_0_Z.getNumSamples();
2487          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2488          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2489          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2490            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2491              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2492              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2493              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2494              double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2495              double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2496              double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2497              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2498            }
2499          }
2500    
2501        }
2502        else {
2503          throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2504        }
2505    
2506      } else {
2507        throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible shapes");
2508      }
2509    
2510      return res;
2511    }
2512    
2513    template <typename UnaryFunction>
2514    Data
2515    C_TensorUnaryOperation(Data const &arg_0,
2516                           UnaryFunction operation)
2517    {
2518      // Interpolate if necessary and find an appropriate function space
2519      Data arg_0_Z = Data(arg_0);
2520    
2521      // Get rank and shape of inputs
2522      int rank0 = arg_0_Z.getDataPointRank();
2523      DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();
2524      int size0 = arg_0_Z.getDataPointSize();
2525    
2526      // Declare output Data object
2527      Data res;
2528    
2529      if (arg_0_Z.isConstant()) {
2530        res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataConstant output
2531        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2532        double *ptr_2 = &((res.getPointDataView().getData())[0]);
2533        tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2534      }
2535      else if (arg_0_Z.isTagged()) {
2536    
2537        // Borrow DataTagged input from Data object
2538        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2539    
2540        // Prepare a DataTagged output 2
2541        res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());   // DataTagged output
2542        res.tag();
2543        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2544    
2545        // Get the views
2546        DataArrayView view_0 = tmp_0->getDefaultValue();
2547        DataArrayView view_2 = tmp_2->getDefaultValue();
2548        // Get the pointers to the actual data
2549        double *ptr_0 = &((view_0.getData())[0]);
2550        double *ptr_2 = &((view_2.getData())[0]);
2551        // Compute a result for the default
2552        tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2553        // Compute a result for each tag
2554        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2555        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2556        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2557          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2558          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2559          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2560          double *ptr_0 = &view_0.getData(0);
2561          double *ptr_2 = &view_2.getData(0);
2562          tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2563        }
2564    
2565      }
2566      else if (arg_0_Z.isExpanded()) {
2567    
2568        res = Data(0.0, shape0, arg_0_Z.getFunctionSpace(),true); // DataExpanded output
2569        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2570        DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2571    
2572        int sampleNo_0,dataPointNo_0;
2573        int numSamples_0 = arg_0_Z.getNumSamples();
2574        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2575        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2576        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2577          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2578            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2579            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2580            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2581            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2582            tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2583          }
2584        }
2585    
2586      }
2587      else {
2588        throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");
2589      }
2590    
2591      return res;
2592    }
2593    
2594  }  }
2595  #endif  #endif

Legend:
Removed from v.757  
changed lines
  Added in v.1693

  ViewVC Help
Powered by ViewVC 1.1.26