escript  Revision_Unversioneddirectory
Data.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2016 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16 
17 
20 #ifndef DATA_H
21 #define DATA_H
22 #include "system_dep.h"
23 
24 #include "DataTypes.h"
25 #include "DataAbstract.h"
26 #include "DataAlgorithm.h"
27 #include "FunctionSpace.h"
28 #include "BinaryOp.h"
29 #include "UnaryOp.h"
30 #include "DataException.h"
31 
32 #ifdef _OPENMP
33 #include <omp.h>
34 #endif
35 
36 #include "esysUtils/Esys_MPI.h"
37 #include <string>
38 #include <algorithm>
39 #include <sstream>
40 
41 #include <boost/shared_ptr.hpp>
42 #include <boost/python/object.hpp>
43 #include <boost/python/tuple.hpp>
44 #include <boost/math/special_functions/bessel.hpp>
45 
46 namespace escript {
47 
48 //
49 // Forward declaration for various implementations of Data.
50 class DataConstant;
51 class DataTagged;
52 class DataExpanded;
53 class DataLazy;
54 
68 class Data {
69 
70  public:
71 
72  // These typedefs allow function names to be cast to pointers
73  // to functions of the appropriate type when calling unaryOp etc.
74  typedef double (*UnaryDFunPtr)(double);
75  typedef double (*BinaryDFunPtr)(double,double);
76 
77 
88  Data();
89 
96  Data(const Data& inData);
97 
105  Data(const Data& inData,
106  const FunctionSpace& what);
107 
113  Data(const DataTypes::ValueType& value,
114  const DataTypes::ShapeType& shape,
115  const FunctionSpace& what=FunctionSpace(),
116  bool expanded=false);
117 
130  Data(double value,
131  const DataTypes::ShapeType& dataPointShape=DataTypes::ShapeType(),
132  const FunctionSpace& what=FunctionSpace(),
133  bool expanded=false);
134 
143  Data(const Data& inData,
144  const DataTypes::RegionType& region);
145 
157  Data(const boost::python::object& value,
158  const FunctionSpace& what=FunctionSpace(),
159  bool expanded=false);
160 
172  Data(const WrappedArray& w, const FunctionSpace& what,
173  bool expanded=false);
174 
175 
186  Data(const boost::python::object& value,
187  const Data& other);
188 
194  Data(double value,
195  const boost::python::tuple& shape=boost::python::make_tuple(),
196  const FunctionSpace& what=FunctionSpace(),
197  bool expanded=false);
198 
204  explicit Data(DataAbstract* underlyingdata);
205 
210  explicit Data(DataAbstract_ptr underlyingdata);
211 
217  ~Data();
218 
223  void
224  copy(const Data& other);
225 
230  Data
231  copySelf();
232 
233 
238  Data
239  delay();
240 
245  void
246  delaySelf();
247 
248 
259  void
260  setProtection();
261 
268  bool
269  isProtected() const;
270 
271 
277  const boost::python::object
278  getValueOfDataPointAsTuple(int dataPointNo);
279 
285  void
286  setValueOfDataPointToPyObject(int dataPointNo, const boost::python::object& py_object);
287 
293  void
294  setValueOfDataPointToArray(int dataPointNo, const boost::python::object&);
295 
301  void
302  setValueOfDataPoint(int dataPointNo, const double);
303 
308  const boost::python::object
309  getValueOfGlobalDataPointAsTuple(int procNo, int dataPointNo);
310 
311 
316  void
317  setTupleForGlobalDataPoint(int id, int proc, boost::python::object);
318 
325  int
326  getTagNumber(int dpno);
327 
328 
334  std::string
335  toString() const;
336 
342  void
343  expand();
344 
352  void
353  tag();
354 
360  void
361  resolve();
362 
369  bool
370  hasNaN();
371 
376  void
377  replaceNaN(double value);
378 
387  void
388  requireWrite();
389 
396  bool
397  isExpanded() const;
398 
405  bool
406  actsExpanded() const;
407 
408 
414  bool
415  isTagged() const;
416 
422  bool
423  isConstant() const;
424 
429  bool
430  isLazy() const;
431 
436  bool
437  isReady() const;
438 
445  bool
446  isEmpty() const;
447 
453  inline
454  const FunctionSpace&
456  {
457  return m_data->getFunctionSpace();
458  }
459 
465  inline
466 // const AbstractDomain&
468  getDomain() const
469  {
470  return getFunctionSpace().getDomain();
471  }
472 
473 
480  inline
481 // const AbstractDomain&
482  Domain_ptr
484  {
486  }
487 
493  inline
494  unsigned int
496  {
497  return m_data->getRank();
498  }
499 
505  inline
506  int
508  {
510  }
516  inline
517  int
519  {
520  return m_data->getNumSamples();
521  }
522 
528  inline
529  int
531  {
532  return m_data->getNumDPPSample();
533  }
534 
541  inline
542  bool numSamplesEqual(int numDataPointsPerSample, int numSamples) const
543  {
544  return (isEmpty() ||
545  (numDataPointsPerSample==getNumDataPointsPerSample() && numSamples==getNumSamples()));
546  }
547 
553  inline
554  bool isDataPointShapeEqual(int rank, const int* dimensions) const
555  {
556  if (isEmpty())
557  return true;
558  const DataTypes::ShapeType givenShape(&dimensions[0],&dimensions[rank]);
559  return (getDataPointShape()==givenShape);
560  }
561 
567  int
568  getNoValues() const
569  {
570  return m_data->getNoValues();
571  }
572 
573 
579  void
580  dump(const std::string fileName) const;
581 
589  const boost::python::object
590  toListOfTuples(bool scalarastuple=true);
591 
592 
601  inline
604 
605 
614  inline
617 
618 
626  inline
628  getDataRO() const;
629 
637  inline
640  {
641  return m_data->getSampleDataByTag(tag);
642  }
643 
652  getDataPointRO(int sampleNo, int dataPointNo);
653 
662  getDataPointRW(int sampleNo, int dataPointNo);
663 
664 
665 
671  inline
673  getDataOffset(int sampleNo,
674  int dataPointNo)
675  {
676  return m_data->getPointOffset(sampleNo,dataPointNo);
677  }
678 
684  inline
685  const DataTypes::ShapeType&
687  {
688  return m_data->getShape();
689  }
690 
696  const boost::python::tuple
697  getShapeTuple() const;
698 
705  int
706  getDataPointSize() const;
707 
714  getLength() const;
715 
721  bool
722  hasNoSamples() const
723  {
724  return getLength()==0;
725  }
726 
736  void
737  setTaggedValueByName(std::string name,
738  const boost::python::object& value);
739 
750  void
751  setTaggedValue(int tagKey,
752  const boost::python::object& value);
753 
765  void
766  setTaggedValueFromCPP(int tagKey,
767  const DataTypes::ShapeType& pointshape,
768  const DataTypes::ValueType& value,
769  int dataOffset=0);
770 
771 
772 
778  void
779  copyWithMask(const Data& other,
780  const Data& mask);
781 
792  void
793  setToZero();
794 
802  Data
803  interpolate(const FunctionSpace& functionspace) const;
804 
806  Data
807  interpolateFromTable3D(const WrappedArray& table, double Amin, double Astep,
808  double undef, Data& B, double Bmin, double Bstep, Data& C,
809  double Cmin, double Cstep, bool check_boundaries);
810 
812  Data
813  interpolateFromTable2D(const WrappedArray& table, double Amin, double Astep,
814  double undef, Data& B, double Bmin, double Bstep,bool check_boundaries);
815 
817  Data
818  interpolateFromTable1D(const WrappedArray& table, double Amin, double Astep,
819  double undef,bool check_boundaries);
820 
821 
823  Data
824  interpolateFromTable3DP(boost::python::object table, double Amin, double Astep,
825  Data& B, double Bmin, double Bstep, Data& C, double Cmin, double Cstep, double undef,bool check_boundaries);
826 
827 
829  Data
830  interpolateFromTable2DP(boost::python::object table, double Amin, double Astep,
831  Data& B, double Bmin, double Bstep, double undef,bool check_boundaries);
832 
834  Data
835  interpolateFromTable1DP(boost::python::object table, double Amin, double Astep,
836  double undef,bool check_boundaries);
837 
839  Data
840  nonuniforminterp(boost::python::object in, boost::python::object out, bool check_boundaries);
841 
843  Data
844  nonuniformslope(boost::python::object in, boost::python::object out, bool check_boundaries);
845 
853  Data
854  gradOn(const FunctionSpace& functionspace) const;
855 
857  Data
858  grad() const;
859 
865  boost::python::object
866  integrateToTuple_const() const;
867 
868 
874  boost::python::object
876 
877 
878 
885  Data
886  oneOver() const;
893  Data
894  wherePositive() const;
895 
902  Data
903  whereNegative() const;
904 
911  Data
912  whereNonNegative() const;
913 
920  Data
921  whereNonPositive() const;
922 
929  Data
930  whereZero(double tol=0.0) const;
931 
938  Data
939  whereNonZero(double tol=0.0) const;
940 
953  double
954  Lsup();
955 
957  double
958  Lsup_const() const;
959 
960 
973  double
974  sup();
975 
977  double
978  sup_const() const;
979 
980 
993  double
994  inf();
995 
997  double
998  inf_const() const;
999 
1000 
1001 
1008  Data
1009  abs() const;
1010 
1017  Data
1018  maxval() const;
1019 
1026  Data
1027  minval() const;
1028 
1037  const boost::python::tuple
1038  minGlobalDataPoint() const;
1039 
1048  const boost::python::tuple
1049  maxGlobalDataPoint() const;
1050 
1051 
1052 
1060  Data
1061  sign() const;
1062 
1069  Data
1070  symmetric() const;
1071 
1078  Data
1079  nonsymmetric() const;
1080 
1087  Data
1088  trace(int axis_offset) const;
1089 
1096  Data
1097  transpose(int axis_offset) const;
1098 
1106  Data
1107  eigenvalues() const;
1108 
1119  const boost::python::tuple
1120  eigenvalues_and_eigenvectors(const double tol=1.e-12) const;
1121 
1128  Data
1129  swapaxes(const int axis0, const int axis1) const;
1130 
1137  Data
1138  erf() const;
1139 
1146  Data
1147  sin() const;
1148 
1155  Data
1156  cos() const;
1157 
1164  Data
1165  bessel(int order, double (*besselfunc) (int,double) );
1166 
1167 
1174  Data
1175  besselFirstKind(int order);
1176 
1183  Data
1184  besselSecondKind(int order);
1185 
1186 
1193  Data
1194  tan() const;
1195 
1202  Data
1203  asin() const;
1204 
1211  Data
1212  acos() const;
1213 
1220  Data
1221  atan() const;
1222 
1229  Data
1230  sinh() const;
1231 
1238  Data
1239  cosh() const;
1240 
1247  Data
1248  tanh() const;
1249 
1256  Data
1257  asinh() const;
1258 
1265  Data
1266  acosh() const;
1267 
1274  Data
1275  atanh() const;
1276 
1283  Data
1284  log10() const;
1285 
1292  Data
1293  log() const;
1294 
1301  Data
1302  exp() const;
1303 
1310  Data
1311  sqrt() const;
1312 
1319  Data
1320  neg() const;
1321 
1329  Data
1330  pos() const;
1331 
1340  Data
1341  powD(const Data& right) const;
1342 
1351  Data
1352  powO(const boost::python::object& right) const;
1353 
1363  Data
1364  rpowO(const boost::python::object& left) const;
1365 
1373  Data& operator+=(const Data& right);
1375  Data& operator+=(const boost::python::object& right);
1376 
1378  Data& operator=(const Data& other);
1379 
1387  Data& operator-=(const Data& right);
1389  Data& operator-=(const boost::python::object& right);
1390 
1398  Data& operator*=(const Data& right);
1400  Data& operator*=(const boost::python::object& right);
1401 
1409  Data& operator/=(const Data& right);
1411  Data& operator/=(const boost::python::object& right);
1412 
1418  Data truedivD(const Data& right);
1419 
1425  Data truedivO(const boost::python::object& right);
1426 
1432  Data rtruedivO(const boost::python::object& left);
1433 
1439  boost::python::object __add__(const boost::python::object& right);
1440 
1441 
1447  boost::python::object __sub__(const boost::python::object& right);
1448 
1454  boost::python::object __rsub__(const boost::python::object& right);
1455 
1461  boost::python::object __mul__(const boost::python::object& right);
1462 
1468  boost::python::object __div__(const boost::python::object& right);
1469 
1475  boost::python::object __rdiv__(const boost::python::object& right);
1476 
1481  Data
1482  matrixInverse() const;
1483 
1489  bool
1490  probeInterpolation(const FunctionSpace& functionspace) const;
1491 
1508  Data
1509  getItem(const boost::python::object& key) const;
1510 
1523  void
1524  setItemD(const boost::python::object& key,
1525  const Data& value);
1526 
1528  void
1529  setItemO(const boost::python::object& key,
1530  const boost::python::object& value);
1531 
1532  // These following public methods should be treated as private.
1533 
1539  template <class UnaryFunction>
1541  inline
1542  void
1543  unaryOp2(UnaryFunction operation);
1544 
1553  Data
1554  getSlice(const DataTypes::RegionType& region) const;
1555 
1565  void
1566  setSlice(const Data& value,
1567  const DataTypes::RegionType& region);
1568 
1574  void
1575  print(void);
1576 
1584  int
1585  get_MPIRank(void) const;
1586 
1594  int
1595  get_MPISize(void) const;
1596 
1603  MPI_Comm
1604  get_MPIComm(void) const;
1605 
1612  DataAbstract*
1613  borrowData(void) const;
1614 
1617  borrowDataPtr(void) const;
1618 
1621  borrowReadyPtr(void) const;
1622 
1623 
1624 
1635 
1636 
1640 
1653 
1654 
1660  size_t
1661  getNumberOfTaggedValues() const;
1662 
1663  protected:
1664 
1665  private:
1666 
1667 template <class BinaryOp>
1668  double
1669 #ifdef ESYS_MPI
1670  lazyAlgWorker(double init, MPI_Op mpiop_type);
1671 #else
1672  lazyAlgWorker(double init);
1673 #endif
1674 
1675  double
1676  LsupWorker() const;
1677 
1678  double
1679  supWorker() const;
1680 
1681  double
1682  infWorker() const;
1683 
1684  boost::python::object
1685  integrateWorker() const;
1686 
1687  void
1688  calc_minGlobalDataPoint(int& ProcNo, int& DataPointNo) const;
1689 
1690  void
1691  calc_maxGlobalDataPoint(int& ProcNo, int& DataPointNo) const;
1692 
1693  // For internal use in Data.cpp only!
1694  // other uses should call the main entry points and allow laziness
1695  Data
1696  minval_nonlazy() const;
1697 
1698  // For internal use in Data.cpp only!
1699  Data
1700  maxval_nonlazy() const;
1701 
1702 
1709  inline
1710  void
1711  operandCheck(const Data& right) const
1712  {
1713  return m_data->operandCheck(*(right.m_data.get()));
1714  }
1715 
1721  template <class BinaryFunction>
1722  inline
1723  double
1724  algorithm(BinaryFunction operation,
1725  double initial_value) const;
1726 
1734  template <class BinaryFunction>
1735  inline
1736  Data
1737  dp_algorithm(BinaryFunction operation,
1738  double initial_value) const;
1739 
1748  template <class BinaryFunction>
1749  inline
1750  void
1751  binaryOp(const Data& right,
1752  BinaryFunction operation);
1753 
1759  void
1760  typeMatchLeft(Data& right) const;
1761 
1767  void
1768  typeMatchRight(const Data& right);
1769 
1775  void
1776  initialise(const DataTypes::ValueType& value,
1777  const DataTypes::ShapeType& shape,
1778  const FunctionSpace& what,
1779  bool expanded);
1780 
1781  void
1782  initialise(const WrappedArray& value,
1783  const FunctionSpace& what,
1784  bool expanded);
1785 
1786  void
1787  initialise(const double value,
1788  const DataTypes::ShapeType& shape,
1789  const FunctionSpace& what,
1790  bool expanded);
1791 
1792  //
1793  // flag to protect the data object against any update
1795  mutable bool m_shared;
1796  bool m_lazy;
1797 
1798  //
1799  // pointer to the actual data object
1800 // boost::shared_ptr<DataAbstract> m_data;
1802 
1803 // If possible please use getReadyPtr instead.
1804 // But see warning below.
1805  const DataReady*
1806  getReady() const;
1807 
1808  DataReady*
1809  getReady();
1810 
1811 
1812 // Be wary of using this for local operations since it (temporarily) increases reference count.
1813 // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1814 // getReady() instead
1816  getReadyPtr();
1817 
1819  getReadyPtr() const;
1820 
1821 
1827  void updateShareStatus(bool nowshared) const
1828  {
1829  m_shared=nowshared; // m_shared is mutable
1830  }
1831 
1832  // In the isShared() method below:
1833  // A problem would occur if m_data (the address pointed to) were being modified
1834  // while the call m_data->is_shared is being executed.
1835  //
1836  // Q: So why do I think this code can be thread safe/correct?
1837  // A: We need to make some assumptions.
1838  // 1. We assume it is acceptable to return true under some conditions when we aren't shared.
1839  // 2. We assume that no constructions or assignments which will share previously unshared
1840  // will occur while this call is executing. This is consistent with the way Data:: and C are written.
1841  //
1842  // This means that the only transition we need to consider, is when a previously shared object is
1843  // not shared anymore. ie. the other objects have been destroyed or a deep copy has been made.
1844  // In those cases the m_shared flag changes to false after m_data has completed changing.
1845  // For any threads executing before the flag switches they will assume the object is still shared.
1846  bool isShared() const
1847  {
1848  return m_shared;
1849 /* if (m_shared) return true;
1850  if (m_data->isShared())
1851  {
1852  updateShareStatus(true);
1853  return true;
1854  }
1855  return false;*/
1856  }
1857 
1859  {
1860  if (isLazy())
1861  {
1862  #ifdef _OPENMP
1863  if (omp_in_parallel())
1864  { // Yes this is throwing an exception out of an omp thread which is forbidden.
1865  throw DataException("Please do not call forceResolve() in a parallel region.");
1866  }
1867  #endif
1868  resolve();
1869  }
1870  }
1871 
1877  {
1878 #ifdef _OPENMP
1879  if (omp_in_parallel())
1880  {
1881 // *((int*)0)=17;
1882  throw DataException("Programming error. Please do not run exclusiveWrite() in multi-threaded sections.");
1883  }
1884 #endif
1885  forceResolve();
1886  if (isShared())
1887  {
1888  DataAbstract* t=m_data->deepCopy();
1890  }
1891 #ifdef EXWRITECHK
1892  m_data->exclusivewritecalled=true;
1893 #endif
1894  }
1895 
1900  {
1901  if (isLazy() || isShared())
1902  {
1903  throw DataException("Programming error. ExclusiveWrite required - please call requireWrite()");
1904  }
1905  }
1906 
1913  void set_m_data(DataAbstract_ptr p);
1914 
1915  friend class DataAbstract; // To allow calls to updateShareStatus
1916  friend class TestDomain; // so its getX will work quickly
1917 #ifdef IKNOWWHATIMDOING
1918  friend ESCRIPT_DLL_API Data applyBinaryCFunction(boost::python::object cfunc, boost::python::tuple shape, escript::Data& d, escript::Data& e);
1919 #endif
1920  friend ESCRIPT_DLL_API Data condEval(escript::Data& mask, escript::Data& trueval, escript::Data& falseval);
1921  friend ESCRIPT_DLL_API Data randomData(const boost::python::tuple& shape, const FunctionSpace& what, long seed, const boost::python::tuple& filter);
1922 
1923 };
1924 
1925 
1926 #ifdef IKNOWWHATIMDOING
1928 Data
1929 applyBinaryCFunction(boost::python::object func, boost::python::tuple shape, escript::Data& d, escript::Data& e);
1930 #endif
1931 
1932 
1934 Data
1935 condEval(escript::Data& mask, escript::Data& trueval, escript::Data& falseval);
1936 
1937 
1938 
1943 Data randomData(const boost::python::tuple& shape,
1944  const FunctionSpace& what,
1945  long seed, const boost::python::tuple& filter);
1946 
1947 
1948 } // end namespace escript
1949 
1950 
1951 // No, this is not supposed to be at the top of the file
1952 // DataAbstact needs to be declared first, then DataReady needs to be fully declared
1953 // so that I can dynamic cast between them below.
1954 #include "DataReady.h"
1955 #include "DataLazy.h"
1956 
1957 namespace escript
1958 {
1959 
1960 inline
1961 const DataReady*
1963 {
1964  const DataReady* dr=dynamic_cast<const DataReady*>(m_data.get());
1965  EsysAssert((dr!=0), "Error - casting to DataReady.");
1966  return dr;
1967 }
1968 
1969 inline
1970 DataReady*
1972 {
1973  DataReady* dr=dynamic_cast<DataReady*>(m_data.get());
1974  EsysAssert((dr!=0), "Error - casting to DataReady.");
1975  return dr;
1976 }
1977 
1978 // Be wary of using this for local operations since it (temporarily) increases reference count.
1979 // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1980 // getReady() instead
1981 inline
1984 {
1985  DataReady_ptr dr=boost::dynamic_pointer_cast<DataReady>(m_data);
1986  EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1987  return dr;
1988 }
1989 
1990 
1991 inline
1994 {
1995  const_DataReady_ptr dr=boost::dynamic_pointer_cast<const DataReady>(m_data);
1996  EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1997  return dr;
1998 }
1999 
2000 inline
2003 {
2004  if (isLazy())
2005  {
2006  throw DataException("Error, attempt to acquire RW access to lazy data. Please call requireWrite() first.");
2007  }
2008 #ifdef EXWRITECHK
2009  if (!getReady()->exclusivewritecalled)
2010  {
2011  throw DataException("Error, call to Data::getSampleDataRW without a preceeding call to requireWrite/exclusiveWrite.");
2012  }
2013 #endif
2014  return getReady()->getSampleDataRW(sampleNo);
2015 }
2016 
2017 inline
2020 {
2021  DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
2022  if (l!=0)
2023  {
2024  size_t offset=0;
2025  const DataTypes::ValueType* res=l->resolveSample(sampleNo,offset);
2026  return &((*res)[offset]);
2027  }
2028  return getReady()->getSampleDataRO(sampleNo);
2029 }
2030 
2031 inline
2034 {
2035  if (isLazy())
2036  {
2037  throw DataException("Programmer error - getDataRO must not be called on Lazy Data.");
2038  }
2039  if (getNumSamples()==0)
2040  {
2041  return 0;
2042  }
2043  else
2044  {
2045  return &(getReady()->getVectorRO()[0]);
2046  }
2047 }
2048 
2049 
2053 inline double rpow(double x,double y)
2054 {
2055  return pow(y,x);
2056 }
2057 
2063 ESCRIPT_DLL_API Data operator+(const Data& left, const Data& right);
2064 
2070 ESCRIPT_DLL_API Data operator-(const Data& left, const Data& right);
2071 
2077 ESCRIPT_DLL_API Data operator*(const Data& left, const Data& right);
2078 
2084 ESCRIPT_DLL_API Data operator/(const Data& left, const Data& right);
2085 
2092 ESCRIPT_DLL_API Data operator+(const Data& left, const boost::python::object& right);
2093 
2100 ESCRIPT_DLL_API Data operator-(const Data& left, const boost::python::object& right);
2101 
2108 ESCRIPT_DLL_API Data operator*(const Data& left, const boost::python::object& right);
2109 
2116 ESCRIPT_DLL_API Data operator/(const Data& left, const boost::python::object& right);
2117 
2124 ESCRIPT_DLL_API Data operator+(const boost::python::object& left, const Data& right);
2125 
2132 ESCRIPT_DLL_API Data operator-(const boost::python::object& left, const Data& right);
2133 
2140 ESCRIPT_DLL_API Data operator*(const boost::python::object& left, const Data& right);
2141 
2148 ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);
2149 
2150 
2151 
2156 ESCRIPT_DLL_API std::ostream& operator<<(std::ostream& o, const Data& data);
2157 
2167 Data
2169  Data& arg_1,
2170  int axis_offset=0,
2171  int transpose=0);
2172 
2178 inline
2179 Data
2180 Data::truedivD(const Data& right)
2181 {
2182  return *this / right;
2183 }
2184 
2190 inline
2191 Data
2192 Data::truedivO(const boost::python::object& right)
2193 {
2194  Data tmp(right, getFunctionSpace(), false);
2195  return truedivD(tmp);
2196 }
2197 
2203 inline
2204 Data
2205 Data::rtruedivO(const boost::python::object& left)
2206 {
2207  Data tmp(left, getFunctionSpace(), false);
2208  return tmp.truedivD(*this);
2209 }
2210 
2216 template <class BinaryFunction>
2217 inline
2218 void
2219 Data::binaryOp(const Data& right,
2220  BinaryFunction operation)
2221 {
2222  //
2223  // if this has a rank of zero promote it to the rank of the RHS
2224  if (getDataPointRank()==0 && right.getDataPointRank()!=0) {
2225  throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");
2226  }
2227 
2228  if (isLazy() || right.isLazy())
2229  {
2230  throw DataException("Programmer error - attempt to call binaryOp with Lazy Data.");
2231  }
2232  //
2233  // initially make the temporary a shallow copy
2234  Data tempRight(right);
2236  FunctionSpace fsr=right.getFunctionSpace();
2237  if (fsl!=fsr) {
2238  signed char intres=fsl.getDomain()->preferredInterpolationOnDomain(fsr.getTypeCode(), fsl.getTypeCode());
2239  if (intres==0)
2240  {
2241  std::string msg="Error - attempt to combine incompatible FunctionSpaces.";
2242  msg+=fsl.toString();
2243  msg+=" ";
2244  msg+=fsr.toString();
2245  throw DataException(msg.c_str());
2246  }
2247  else if (intres==1)
2248  {
2249  // an interpolation is required so create a new Data
2250  tempRight=Data(right,fsl);
2251  }
2252  else // reverse interpolation preferred
2253  {
2254  // interpolate onto the RHS function space
2255  Data tempLeft(*this,fsr);
2256  set_m_data(tempLeft.m_data);
2257  }
2258  }
2259  operandCheck(tempRight);
2260  //
2261  // ensure this has the right type for the RHS
2262  typeMatchRight(tempRight);
2263  //
2264  // Need to cast to the concrete types so that the correct binaryOp
2265  // is called.
2266  if (isExpanded()) {
2267  //
2268  // Expanded data will be done in parallel, the right hand side can be
2269  // of any data type
2270  DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2271  EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2272  escript::binaryOp(*leftC,*(tempRight.getReady()),operation);
2273  } else if (isTagged()) {
2274  //
2275  // Tagged data is operated on serially, the right hand side can be
2276  // either DataConstant or DataTagged
2277  DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
2278  EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
2279  if (right.isTagged()) {
2280  DataTagged* rightC=dynamic_cast<DataTagged*>(tempRight.m_data.get());
2281  EsysAssert((rightC!=0), "Programming error - casting to DataTagged.");
2282  escript::binaryOp(*leftC,*rightC,operation);
2283  } else {
2284  DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
2285  EsysAssert((rightC!=0), "Programming error - casting to DataConstant.");
2286  escript::binaryOp(*leftC,*rightC,operation);
2287  }
2288  } else if (isConstant()) {
2289  DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2290  DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
2291  EsysAssert((leftC!=0 && rightC!=0), "Programming error - casting to DataConstant.");
2292  escript::binaryOp(*leftC,*rightC,operation);
2293  }
2294 }
2295 
2303 template <class BinaryFunction>
2304 inline
2305 double
2306 Data::algorithm(BinaryFunction operation, double initial_value) const
2307 {
2308  if (isExpanded()) {
2309  DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2310  EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2311  return escript::algorithm(*leftC,operation,initial_value);
2312  } else if (isTagged()) {
2313  DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
2314  EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
2315  return escript::algorithm(*leftC,operation,initial_value);
2316  } else if (isConstant()) {
2317  DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2318  EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
2319  return escript::algorithm(*leftC,operation,initial_value);
2320  } else if (isEmpty()) {
2321  throw DataException("Error - Operations (algorithm) not permitted on instances of DataEmpty.");
2322  } else if (isLazy()) {
2323  throw DataException("Error - Operations not permitted on instances of DataLazy.");
2324  } else {
2325  throw DataException("Error - Data encapsulates an unknown type.");
2326  }
2327 }
2328 
2337 template <class BinaryFunction>
2338 inline
2339 Data
2340 Data::dp_algorithm(BinaryFunction operation, double initial_value) const
2341 {
2342  if (isEmpty()) {
2343  throw DataException("Error - Operations (dp_algorithm) not permitted on instances of DataEmpty.");
2344  }
2345  else if (isExpanded()) {
2347  DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
2348  DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
2349  EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
2350  EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
2351  escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
2352  return result;
2353  }
2354  else if (isTagged()) {
2355  DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
2356  EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
2357  DataTypes::ValueType defval(1);
2358  defval[0]=0;
2359  DataTagged* resultT=new DataTagged(getFunctionSpace(), DataTypes::scalarShape, defval, dataT);
2360  escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
2361  return Data(resultT); // note: the Data object now owns the resultT pointer
2362  }
2363  else if (isConstant()) {
2365  DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
2366  DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
2367  EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
2368  EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
2369  escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
2370  return result;
2371  } else if (isLazy()) {
2372  throw DataException("Error - Operations not permitted on instances of DataLazy.");
2373  } else {
2374  throw DataException("Error - Data encapsulates an unknown type.");
2375  }
2376 }
2377 
2385 template <typename BinaryFunction>
2386 inline
2387 Data
2389  Data const &arg_1,
2390  BinaryFunction operation)
2391 {
2392  if (arg_0.isEmpty() || arg_1.isEmpty())
2393  {
2394  throw DataException("Error - Operations (C_TensorBinaryOperation) not permitted on instances of DataEmpty.");
2395  }
2396  if (arg_0.isLazy() || arg_1.isLazy())
2397  {
2398  throw DataException("Error - Operations not permitted on lazy data.");
2399  }
2400  // Interpolate if necessary and find an appropriate function space
2401  Data arg_0_Z, arg_1_Z;
2402  FunctionSpace fsl=arg_0.getFunctionSpace();
2403  FunctionSpace fsr=arg_1.getFunctionSpace();
2404  if (fsl!=fsr) {
2405  signed char intres=fsl.getDomain()->preferredInterpolationOnDomain(fsr.getTypeCode(), fsl.getTypeCode());
2406  if (intres==0)
2407  {
2408  std::string msg="Error - C_TensorBinaryOperation: arguments have incompatible function spaces.";
2409  msg+=fsl.toString();
2410  msg+=" ";
2411  msg+=fsr.toString();
2412  throw DataException(msg.c_str());
2413  }
2414  else if (intres==1)
2415  {
2416  arg_1_Z=arg_1.interpolate(arg_0.getFunctionSpace());
2417  arg_0_Z =Data(arg_0);
2418  }
2419  else // reverse interpolation preferred
2420  {
2421  arg_0_Z = arg_0.interpolate(arg_1.getFunctionSpace());
2422  arg_1_Z = Data(arg_1);
2423  }
2424  } else {
2425  arg_0_Z = Data(arg_0);
2426  arg_1_Z = Data(arg_1);
2427  }
2428  // Get rank and shape of inputs
2429  int rank0 = arg_0_Z.getDataPointRank();
2430  int rank1 = arg_1_Z.getDataPointRank();
2431  DataTypes::ShapeType shape0 = arg_0_Z.getDataPointShape();
2432  DataTypes::ShapeType shape1 = arg_1_Z.getDataPointShape();
2433  int size0 = arg_0_Z.getDataPointSize();
2434  int size1 = arg_1_Z.getDataPointSize();
2435  // Declare output Data object
2436  Data res;
2437 
2438  if (shape0 == shape1) {
2439  if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2440  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
2441  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2442  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2443  double *ptr_2 = &(res.getDataAtOffsetRW(0));
2444 
2445  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2446  }
2447  else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2448 
2449  // Prepare the DataConstant input
2450  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2451 
2452  // Borrow DataTagged input from Data object
2453  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2454 
2455  // Prepare a DataTagged output 2
2456  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
2457  res.tag();
2458  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2459 
2460  // Prepare offset into DataConstant
2461  int offset_0 = tmp_0->getPointOffset(0,0);
2462  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2463 
2464  // Get the pointers to the actual data
2465  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2466  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2467 
2468  // Compute a result for the default
2469  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2470  // Compute a result for each tag
2471  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2472  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2473  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2474  tmp_2->addTag(i->first);
2475  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2476  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2477 
2478  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2479  }
2480 
2481  }
2482  else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2483  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2484  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2485  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2486  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2487 
2488  int sampleNo_1,dataPointNo_1;
2489  int numSamples_1 = arg_1_Z.getNumSamples();
2490  int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2491  int offset_0 = tmp_0->getPointOffset(0,0);
2492  res.requireWrite();
2493  #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2494  for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2495  for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2496  int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2497  int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2498  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2499  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2500  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2501  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2502  }
2503  }
2504 
2505  }
2506  else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2507  // Borrow DataTagged input from Data object
2508  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2509 
2510  // Prepare the DataConstant input
2511  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2512 
2513  // Prepare a DataTagged output 2
2514  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
2515  res.tag();
2516  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2517 
2518  // Prepare offset into DataConstant
2519  int offset_1 = tmp_1->getPointOffset(0,0);
2520 
2521  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2522  // Get the pointers to the actual data
2523  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2524  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2525  // Compute a result for the default
2526  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2527  // Compute a result for each tag
2528  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2529  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2530  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2531  tmp_2->addTag(i->first);
2532  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2533  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2534  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2535  }
2536 
2537  }
2538  else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2539  // Borrow DataTagged input from Data object
2540  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2541 
2542  // Borrow DataTagged input from Data object
2543  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2544 
2545  // Prepare a DataTagged output 2
2546  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2547  res.tag(); // DataTagged output
2548  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2549 
2550  // Get the pointers to the actual data
2551  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2552  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2553  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2554 
2555  // Compute a result for the default
2556  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2557  // Merge the tags
2558  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2559  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2560  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2561  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2562  tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2563  }
2564  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2565  tmp_2->addTag(i->first);
2566  }
2567  // Compute a result for each tag
2568  const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2569  for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2570 
2571  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2572  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2573  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2574 
2575  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2576  }
2577 
2578  }
2579  else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2580  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2581  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2582  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2583  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2584  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2585 
2586  int sampleNo_0,dataPointNo_0;
2587  int numSamples_0 = arg_0_Z.getNumSamples();
2588  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2589  res.requireWrite();
2590  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2591  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2592  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2593  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2594  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2595  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2596  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2597  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2598  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2599  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2600  }
2601  }
2602 
2603  }
2604  else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2605  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2606  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2607  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2608  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2609 
2610  int sampleNo_0,dataPointNo_0;
2611  int numSamples_0 = arg_0_Z.getNumSamples();
2612  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2613  int offset_1 = tmp_1->getPointOffset(0,0);
2614  res.requireWrite();
2615  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2616  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2617  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2618  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2619  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2620 
2621  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2622  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2623  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2624 
2625 
2626  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2627  }
2628  }
2629 
2630  }
2631  else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2632  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2633  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2634  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2635  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2636  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2637 
2638  int sampleNo_0,dataPointNo_0;
2639  int numSamples_0 = arg_0_Z.getNumSamples();
2640  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2641  res.requireWrite();
2642  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2643  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2644  int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2645  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2646  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2647  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2648  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2649  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2650  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2651  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2652  }
2653  }
2654 
2655  }
2656  else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2657  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2658  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2659  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2660  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2661  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2662 
2663  int sampleNo_0,dataPointNo_0;
2664  int numSamples_0 = arg_0_Z.getNumSamples();
2665  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2666  res.requireWrite();
2667  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2668  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2669  dataPointNo_0=0;
2670 // for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2671  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2672  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2673  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2674  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2675  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2676  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2677  tensor_binary_operation(size0*numDataPointsPerSample_0, ptr_0, ptr_1, ptr_2, operation);
2678 // }
2679  }
2680 
2681  }
2682  else {
2683  throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2684  }
2685 
2686  } else if (0 == rank0) {
2687  if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2688  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataConstant output
2689  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2690  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2691  double *ptr_2 = &(res.getDataAtOffsetRW(0));
2692  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2693  }
2694  else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2695 
2696  // Prepare the DataConstant input
2697  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2698 
2699  // Borrow DataTagged input from Data object
2700  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2701 
2702  // Prepare a DataTagged output 2
2703  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataTagged output
2704  res.tag();
2705  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2706 
2707  // Prepare offset into DataConstant
2708  int offset_0 = tmp_0->getPointOffset(0,0);
2709  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2710 
2711  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2712  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2713 
2714  // Compute a result for the default
2715  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2716  // Compute a result for each tag
2717  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2718  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2719  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2720  tmp_2->addTag(i->first);
2721  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2722  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2723  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2724  }
2725 
2726  }
2727  else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2728 
2729  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2730  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2731  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2732  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2733 
2734  int sampleNo_1;
2735  int numSamples_1 = arg_1_Z.getNumSamples();
2736  int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2737  int offset_0 = tmp_0->getPointOffset(0,0);
2738  const double *ptr_src = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2739  double ptr_0 = ptr_src[0];
2740  int size = size1*numDataPointsPerSample_1;
2741  res.requireWrite();
2742  #pragma omp parallel for private(sampleNo_1) schedule(static)
2743  for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2744 // for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2745  int offset_1 = tmp_1->getPointOffset(sampleNo_1,0);
2746  int offset_2 = tmp_2->getPointOffset(sampleNo_1,0);
2747 // const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2748  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2749  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2750  tensor_binary_operation(size, ptr_0, ptr_1, ptr_2, operation);
2751 
2752 // }
2753  }
2754 
2755  }
2756  else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2757 
2758  // Borrow DataTagged input from Data object
2759  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2760 
2761  // Prepare the DataConstant input
2762  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2763 
2764  // Prepare a DataTagged output 2
2765  res = Data(0.0, shape1, arg_0_Z.getFunctionSpace()); // DataTagged output
2766  res.tag();
2767  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2768 
2769  // Prepare offset into DataConstant
2770  int offset_1 = tmp_1->getPointOffset(0,0);
2771  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2772 
2773  // Get the pointers to the actual data
2774  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2775  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2776 
2777 
2778  // Compute a result for the default
2779  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2780  // Compute a result for each tag
2781  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2782  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2783  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2784  tmp_2->addTag(i->first);
2785  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2786  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2787 
2788  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2789  }
2790 
2791  }
2792  else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2793 
2794  // Borrow DataTagged input from Data object
2795  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2796 
2797  // Borrow DataTagged input from Data object
2798  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2799 
2800  // Prepare a DataTagged output 2
2801  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());
2802  res.tag(); // DataTagged output
2803  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2804 
2805  // Get the pointers to the actual data
2806  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2807  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2808  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2809 
2810  // Compute a result for the default
2811  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2812  // Merge the tags
2813  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2814  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2815  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2816  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2817  tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2818  }
2819  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2820  tmp_2->addTag(i->first);
2821  }
2822  // Compute a result for each tag
2823  const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2824  for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2825  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2826  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2827  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2828 
2829  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2830  }
2831 
2832  }
2833  else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2834 
2835  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2836  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2837  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2838  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2839  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2840 
2841  int sampleNo_0,dataPointNo_0;
2842  int numSamples_0 = arg_0_Z.getNumSamples();
2843  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2844  res.requireWrite();
2845  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2846  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2847  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2848  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2849  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2850  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2851  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2852  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2853  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2854  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2855  }
2856  }
2857 
2858  }
2859  else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2860  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2861  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2862  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2863  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2864 
2865  int sampleNo_0,dataPointNo_0;
2866  int numSamples_0 = arg_0_Z.getNumSamples();
2867  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2868  int offset_1 = tmp_1->getPointOffset(0,0);
2869  res.requireWrite();
2870  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2871  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2872  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2873  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2874  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2875  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2876  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2877  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2878  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2879  }
2880  }
2881 
2882 
2883  }
2884  else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2885 
2886  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2887  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2888  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2889  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2890  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2891 
2892  int sampleNo_0,dataPointNo_0;
2893  int numSamples_0 = arg_0_Z.getNumSamples();
2894  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2895  res.requireWrite();
2896  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2897  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2898  int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2899  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2900  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2901  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2902  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2903  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2904  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2905  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2906  }
2907  }
2908 
2909  }
2910  else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2911 
2912  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2913  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2914  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2915  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2916  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2917 
2918  int sampleNo_0,dataPointNo_0;
2919  int numSamples_0 = arg_0_Z.getNumSamples();
2920  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2921  res.requireWrite();
2922  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2923  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2924  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2925  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2926  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2927  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2928  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2929  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2930  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2931  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2932  }
2933  }
2934 
2935  }
2936  else {
2937  throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2938  }
2939 
2940  } else if (0 == rank1) {
2941  if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2942  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
2943  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2944  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2945  double *ptr_2 = &(res.getDataAtOffsetRW(0));
2946  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2947  }
2948  else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2949 
2950  // Prepare the DataConstant input
2951  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2952 
2953  // Borrow DataTagged input from Data object
2954  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2955 
2956  // Prepare a DataTagged output 2
2957  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
2958  res.tag();
2959  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2960 
2961  // Prepare offset into DataConstant
2962  int offset_0 = tmp_0->getPointOffset(0,0);
2963  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2964 
2965  //Get the pointers to the actual data
2966  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2967  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2968 
2969  // Compute a result for the default
2970  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2971  // Compute a result for each tag
2972  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2973  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2974  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2975  tmp_2->addTag(i->first);
2976  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2977  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2978  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2979  }
2980  }
2981  else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2982 
2983  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2984  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2985  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2986  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2987 
2988  int sampleNo_1,dataPointNo_1;
2989  int numSamples_1 = arg_1_Z.getNumSamples();
2990  int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2991  int offset_0 = tmp_0->getPointOffset(0,0);
2992  res.requireWrite();
2993  #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2994  for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2995  for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2996  int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2997  int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2998  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2999  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3000  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3001  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3002  }
3003  }
3004 
3005  }
3006  else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
3007 
3008  // Borrow DataTagged input from Data object
3009  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
3010 
3011  // Prepare the DataConstant input
3012  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
3013 
3014  // Prepare a DataTagged output 2
3015  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
3016  res.tag();
3017  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
3018 
3019  // Prepare offset into DataConstant
3020  int offset_1 = tmp_1->getPointOffset(0,0);
3021  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3022  // Get the pointers to the actual data
3023  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
3024  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
3025  // Compute a result for the default
3026  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3027  // Compute a result for each tag
3028  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
3029  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
3030  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
3031  tmp_2->addTag(i->first);
3032  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
3033  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
3034  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3035  }
3036 
3037  }
3038  else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
3039 
3040  // Borrow DataTagged input from Data object
3041  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
3042 
3043  // Borrow DataTagged input from Data object
3044  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
3045 
3046  // Prepare a DataTagged output 2
3047  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
3048  res.tag(); // DataTagged output
3049  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
3050 
3051  // Get the pointers to the actual data
3052  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
3053  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
3054  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
3055 
3056  // Compute a result for the default
3057  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3058  // Merge the tags
3059  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
3060  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
3061  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
3062  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
3063  tmp_2->addTag(i->first); // use tmp_2 to get correct shape
3064  }
3065  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
3066  tmp_2->addTag(i->first);
3067  }
3068  // Compute a result for each tag
3069  const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
3070  for (i=lookup_2.begin();i!=lookup_2.end();i++) {
3071  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
3072  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
3073  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
3074  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3075  }
3076 
3077  }
3078  else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
3079 
3080  // After finding a common function space above the two inputs have the same numSamples and num DPPS
3081  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3082  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
3083  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
3084  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3085 
3086  int sampleNo_0,dataPointNo_0;
3087  int numSamples_0 = arg_0_Z.getNumSamples();
3088  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3089  res.requireWrite();
3090  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3091  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3092  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
3093  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3094  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3095  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
3096  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3097  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3098  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3099  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3100  }
3101  }
3102 
3103  }
3104  else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
3105  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3106  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3107  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
3108  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3109 
3110  int sampleNo_0;
3111  int numSamples_0 = arg_0_Z.getNumSamples();
3112  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3113  int offset_1 = tmp_1->getPointOffset(0,0);
3114  const double *ptr_src = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3115  double ptr_1 = ptr_src[0];
3116  int size = size0 * numDataPointsPerSample_0;
3117  res.requireWrite();
3118  #pragma omp parallel for private(sampleNo_0) schedule(static)
3119  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3120 // for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3121  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0);
3122  int offset_2 = tmp_2->getPointOffset(sampleNo_0,0);
3123  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3124 // const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3125  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3126  tensor_binary_operation(size, ptr_0, ptr_1, ptr_2, operation);
3127 // }
3128  }
3129 
3130 
3131  }
3132  else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
3133 
3134  // After finding a common function space above the two inputs have the same numSamples and num DPPS
3135  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3136  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3137  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
3138  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3139 
3140  int sampleNo_0,dataPointNo_0;
3141  int numSamples_0 = arg_0_Z.getNumSamples();
3142  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3143  res.requireWrite();
3144  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3145  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3146  int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
3147  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3148  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3149  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
3150  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3151  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3152  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3153  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3154  }
3155  }
3156 
3157  }
3158  else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
3159 
3160  // After finding a common function space above the two inputs have the same numSamples and num DPPS
3161  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3162  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3163  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
3164  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3165 
3166  int sampleNo_0,dataPointNo_0;
3167  int numSamples_0 = arg_0_Z.getNumSamples();
3168  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3169  res.requireWrite();
3170  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3171  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3172  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3173  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
3174  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
3175  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3176  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3177  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3178  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3179  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3180  }
3181  }
3182 
3183  }
3184  else {
3185  throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
3186  }
3187 
3188  } else {
3189  throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible shapes");
3190  }
3191 
3192  return res;
3193 }
3194 
3195 template <typename UnaryFunction>
3196 Data
3198  UnaryFunction operation)
3199 {
3200  if (arg_0.isEmpty()) // do this before we attempt to interpolate
3201  {
3202  throw DataException("Error - Operations (C_TensorUnaryOperation) not permitted on instances of DataEmpty.");
3203  }
3204  if (arg_0.isLazy())
3205  {
3206  throw DataException("Error - Operations not permitted on lazy data.");
3207  }
3208  // Interpolate if necessary and find an appropriate function space
3209  Data arg_0_Z = Data(arg_0);
3210 
3211  // Get rank and shape of inputs
3212  const DataTypes::ShapeType& shape0 = arg_0_Z.getDataPointShape();
3213  int size0 = arg_0_Z.getDataPointSize();
3214 
3215  // Declare output Data object
3216  Data res;
3217 
3218  if (arg_0_Z.isConstant()) {
3219  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataConstant output
3220  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
3221  double *ptr_2 = &(res.getDataAtOffsetRW(0));
3222  tensor_unary_operation(size0, ptr_0, ptr_2, operation);
3223  }
3224  else if (arg_0_Z.isTagged()) {
3225 
3226  // Borrow DataTagged input from Data object
3227  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
3228 
3229  // Prepare a DataTagged output 2
3230  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
3231  res.tag();
3232  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
3233 
3234  // Get the pointers to the actual data
3235  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
3236  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
3237  // Compute a result for the default
3238  tensor_unary_operation(size0, ptr_0, ptr_2, operation);
3239  // Compute a result for each tag
3240  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
3241  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
3242  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
3243  tmp_2->addTag(i->first);
3244  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
3245  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
3246  tensor_unary_operation(size0, ptr_0, ptr_2, operation);
3247  }
3248 
3249  }
3250  else if (arg_0_Z.isExpanded()) {
3251 
3252  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace(),true); // DataExpanded output
3253  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3254  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3255 
3256  int sampleNo_0,dataPointNo_0;
3257  int numSamples_0 = arg_0_Z.getNumSamples();
3258  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3259  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3260  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3261  dataPointNo_0=0;
3262 // for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3263  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
3264  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3265  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3266  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3267  tensor_unary_operation(size0*numDataPointsPerSample_0, ptr_0, ptr_2, operation);
3268 // }
3269  }
3270  }
3271  else {
3272  throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");
3273  }
3274 
3275  return res;
3276 }
3277 
3278 }
3279 #endif
double sup_const() const
Definition: Data.cpp:1703
Definition: FunctionSpace.h:34
DataVector implements an arbitrarily long vector of data values. DataVector is the underlying data co...
Definition: DataVector.h:44
std::ostream & operator<<(std::ostream &o, const Data &data)
Output operator.
Definition: Data.cpp:2859
DataAbstract_ptr m_data
Definition: Data.h:1801
const boost::python::object toListOfTuples(bool scalarastuple=true)
returns the values of the object as a list of tuples (one for each datapoint).
Definition: Data.cpp:1072
int getNumDataPointsPerSample() const
Return the number of data points per sample.
Definition: Data.h:530
DataTypes::ValueType::reference getDataAtOffsetRW(DataTypes::ValueType::size_type i)
Definition: Data.cpp:3277
Data asin() const
Return the asin of each data point of this Data object.
Definition: Data.cpp:1524
const DataTypes::ShapeType & getDataPointShape() const
Return a reference to the data point shape.
Definition: Data.h:686
double LsupWorker() const
Definition: Data.cpp:1851
Data wherePositive() const
Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.
Definition: Data.cpp:968
void addTag(int tagKey)
addTag - does not modify the default value for this object. ** Not unit tested ** ...
Definition: DataTagged.cpp:432
void setItemO(const boost::python::object &key, const boost::python::object &value)
Definition: Data.cpp:2713
MPI_Comm get_MPIComm(void) const
return the MPI rank number of the local data MPI_COMM_WORLD is assumed and returned.
Definition: Data.cpp:4101
Data whereNonNegative() const
Return a Data with a 1 for +ive or 0 values and a 0 for -ive values.
Definition: Data.cpp:982
int get_MPIRank(void) const
return the MPI rank number of the local data MPI_COMM_WORLD is assumed and the result of MPI_Comm_siz...
Definition: Data.cpp:4089
unsigned int getDataPointRank() const
Return the rank of the point data.
Definition: Data.h:495
int getNumSamples() const
Return the number of samples.
Definition: Data.h:518
virtual DataTypes::ValueType::size_type getPointOffset(int sampleNo, int dataPointNo) const
Return the offset for the given given data point. This returns the offset in bytes for the given poin...
Definition: DataExpanded.cpp:305
boost::python::object __add__(const boost::python::object &right)
wrapper for python add operation
Definition: Data.cpp:4536
void dp_algorithm(const DataExpanded &data, DataExpanded &result, BinaryFunction operation, double initial_value)
Perform the given data-point reduction operation on all data-points in data, storing results in corre...
Definition: DataAlgorithm.h:274
void checkExclusiveWrite()
checks if caller can have exclusive write to the object
Definition: Data.h:1899
boost::python::object __rsub__(const boost::python::object &right)
wrapper for python reverse subtract operation
Definition: Data.cpp:4584
DataReady_ptr borrowReadyPtr(void) const
Definition: Data.cpp:3240
void setToZero()
set all values to zero
Definition: Data.cpp:580
Data getItem(const boost::python::object &key) const
Returns a slice from this Data object.
Definition: Data.cpp:2691
double(* UnaryDFunPtr)(double)
Definition: Data.h:74
Data dp_algorithm(BinaryFunction operation, double initial_value) const
Reduce each data-point in this Data object using the given operation. Return a Data object with the s...
Definition: Data.h:2340
bool isExpanded() const
Return true if this Data is expanded.
Definition: Data.cpp:832
double Lsup()
Return the maximum absolute value of this Data object.
Definition: Data.cpp:1682
Data tanh() const
Return the tanh of each data point of this Data object.
Definition: Data.cpp:1560
Data delay()
produce a delayed evaluation version of this Data.
Definition: Data.cpp:551
Data & operator*=(const Data &right)
Overloaded operator *=.
Definition: Data.cpp:2474
void tag()
If possible convert this Data to DataTagged. This will only allow Constant data to be converted to ta...
Definition: Data.cpp:919
void set_m_data(DataAbstract_ptr p)
Modify the data abstract hosted by this Data object For internal use only. Passing a pointer to null ...
Definition: Data.cpp:433
boost::python::object integrateToTuple()
Calculate the integral over the function space domain as a python tuple.
Definition: Data.cpp:1353
virtual DataAbstract * deepCopy()=0
Return a deep copy of the current object.
int getDataPointSize() const
Return the size of the data point. It is the product of the data point shape dimensions.
Definition: Data.cpp:1054
boost::shared_ptr< AbstractDomain > Domain_ptr
Definition: AbstractDomain.h:36
double lazyAlgWorker(double init)
Definition: Data.cpp:1769
Definition: DataReady.h:35
double inf_const() const
Definition: Data.cpp:1734
bool m_shared
Definition: Data.h:1795
std::vector< std::pair< int, int > > RegionType
Definition: DataTypes.h:39
Definition: AbstractContinuousDomain.cpp:24
bool probeInterpolation(const FunctionSpace &functionspace) const
Returns true if this can be interpolated to functionspace.
Definition: Data.cpp:1020
boost::shared_ptr< const DataReady > const_DataReady_ptr
Definition: DataAbstract.h:59
DataTypes::ValueType::const_reference getDataAtOffsetRO(DataTypes::ValueType::size_type i)
Return a pointer to the beginning of the datapoint at the specified offset. TODO Eventually these sho...
Definition: Data.cpp:3285
double infWorker() const
Definition: Data.cpp:1925
void replaceNaN(double value)
replaces all NaN values with value
Definition: Data.cpp:1837
void setItemD(const boost::python::object &key, const Data &value)
Copies slice from value into this Data object.
Definition: Data.cpp:2721
void tensor_unary_operation(const int size, const double *arg1, double *argRes, UnaryFunction operation)
Definition: LocalOps.h:512
virtual ValueType::size_type getPointOffset(int sampleNo, int dataPointNo) const
getPointOffset
Definition: DataTagged.cpp:526
Data atan() const
Return the atan of each data point of this Data object.
Definition: Data.cpp:1539
Data operator-(const Data &left, const Data &right)
Operator- Takes two Data objects.
Definition: Data.cpp:2583
Data rpowO(const boost::python::object &left) const
Return the given power of each data point of this boost python object.
Definition: Data.cpp:2550
Data besselSecondKind(int order)
Return the Bessel function of the second kind for each data point of this Data object.
Definition: Data.cpp:1407
void print(void)
print the data values to stdout. Used for debugging
Definition: Data.cpp:4041
void typeMatchRight(const Data &right)
Convert the data type of this to match the RHS.
Definition: Data.cpp:2768
void setSlice(const Data &value, const DataTypes::RegionType &region)
Copy the specified slice from the given value into this Data object.
Definition: Data.cpp:2737
Data neg() const
Return the negation of each data point of this Data object.
Definition: Data.cpp:1640
DataConstant stores a single data point which represents the entire function space.
Definition: DataConstant.h:37
static const ShapeType scalarShape
Use this instead of creating empty shape objects for scalars.
Definition: DataTypes.h:42
std::string toString() const
Write the data as a string. For large amounts of data, a summary is printed.
Definition: Data.cpp:3248
void setValueOfDataPointToPyObject(int dataPointNo, const boost::python::object &py_object)
sets the values of a data-point from a python object on this process
Definition: Data.cpp:1180
Data log10() const
Return the log to base 10 of each data point of this Data object.
Definition: Data.cpp:1612
Data interpolateFromTable2DP(boost::python::object table, double Amin, double Astep, Data &B, double Bmin, double Bstep, double undef, bool check_boundaries)
Definition: Data.cpp:3338
Data getSlice(const DataTypes::RegionType &region) const
Return a Data object containing the specified slice of this Data object.
Definition: Data.cpp:2705
bool isLazy() const
Return true if this Data is lazy.
Definition: Data.cpp:867
Data swapaxes(const int axis0, const int axis1) const
swaps the components axis0 and axis1
Definition: Data.cpp:1999
void setTaggedValueByName(std::string name, const boost::python::object &value)
Assign the given value to the tag assocciated with name. Implicitly converts this object to type Data...
Definition: Data.cpp:2792
Give a short description of what DataExpanded does.
Definition: DataExpanded.h:44
bool numSamplesEqual(int numDataPointsPerSample, int numSamples) const
Returns true if the number of data points per sample and the number of samples match the respective a...
Definition: Data.h:542
double Lsup_const() const
Definition: Data.cpp:1672
const DataMapType & getTagLookup() const
getTagLookup
Definition: DataTagged.h:629
DataReady_ptr getReadyPtr()
Definition: Data.h:1983
double * getSampleDataRW(ValueType::size_type sampleNo)
Return the sample data for the given sample number.
Definition: DataReady.h:120
DataTypes::ValueType::const_reference getDataByTagRO(int tag, DataTypes::ValueType::size_type i) const
Definition: DataTagged.cpp:563
virtual DataTypes::ValueType::size_type getPointOffset(int sampleNo, int dataPointNo) const
Return the offset for the given sample. This is a somewhat artificial notion but returns the offset i...
Definition: DataConstant.cpp:137
Data trace(int axis_offset) const
Return the trace of a matrix.
Definition: Data.cpp:2114
friend Data randomData(const boost::python::tuple &shape, const FunctionSpace &what, long seed, const boost::python::tuple &filter)
Create a new Expanded Data object filled with pseudo-random data.
Data maxval_nonlazy() const
Definition: Data.cpp:1973
Data symmetric() const
Return the symmetric part of a matrix which is half the matrix plus its transpose.
Definition: Data.cpp:2056
(Testing use only) Provides a domain to wrap a collection of values.
Definition: TestDomain.h:41
Data nonuniformslope(boost::python::object in, boost::python::object out, bool check_boundaries)
Definition: Data.cpp:3962
Domain_ptr getDomainPython() const
Return the function space domain. Internal use only! This gets around some python difficulties by cas...
Definition: FunctionSpace.cpp:113
std::vector< int > ShapeType
The shape of a single datapoint.
Definition: DataTypes.h:38
boost::python::object __sub__(const boost::python::object &right)
wrapper for python subtract operation
Definition: Data.cpp:4560
DataAbstract_ptr borrowDataPtr(void) const
Definition: Data.cpp:3233
Simulates a full dataset accessible via sampleNo and dataPointNo.
Definition: DataTagged.h:44
boost::shared_ptr< DataAbstract > DataAbstract_ptr
Definition: DataAbstract.h:51
void operandCheck(const Data &right) const
Check *this and the right operand are compatible. Throws an exception if they aren&#39;t.
Definition: Data.h:1711
Data interpolateFromTable1D(const WrappedArray &table, double Amin, double Astep, double undef, bool check_boundaries)
Definition: Data.cpp:3356
Data sqrt() const
Return the square root of each data point of this Data object.
Definition: Data.cpp:1665
DataAbstract::ValueType::value_type * getSampleDataByTag(int tag)
Return the sample data for the given tag. If an attempt is made to access data that isn&#39;t tagged an e...
Definition: Data.h:639
size_t getNumberOfTaggedValues() const
For tagged Data returns the number of tags with values. For non-tagged data will return 0 (even Data ...
Definition: Data.cpp:4467
bool isConstant() const
Return true if this Data is constant.
Definition: Data.cpp:860
DataTypes::ValueType & getExpandedVectorReference()
Ensures that the Data is expanded and returns its underlying vector Does not check for exclusive writ...
Definition: Data.cpp:4458
void binaryOp(const Data &right, BinaryFunction operation)
Perform the given binary operation on all of the data&#39;s elements. The underlying type of the right ha...
Definition: Data.h:2219
Data cos() const
Return the cos of each data point of this Data object.
Definition: Data.cpp:1510
Data minval() const
Return the minimum value of each data point of this Data object.
Definition: Data.cpp:1991
Data matrixInverse() const
return inverse of matricies.
Definition: Data.cpp:2524
int MPI_Comm
Definition: Esys_MPI.h:38
Data cosh() const
Return the cosh of each data point of this Data object.
Definition: Data.cpp:1553
Data C_GeneralTensorProduct(Data &arg_0, Data &arg_1, int axis_offset=0, int transpose=0)
Compute a tensor product of two Data objects.
Definition: Data.cpp:2866
bool m_protected
Definition: Data.h:1794
bool isTagged() const
Return true if this Data is tagged.
Definition: Data.cpp:846
Data gradOn(const FunctionSpace &functionspace) const
Calculates the gradient of the data at the data points of functionspace. If functionspace is not pres...
Definition: Data.cpp:1026
Data transpose(int axis_offset) const
Transpose each data point of this Data object around the given axis.
Definition: Data.cpp:2169
Data & operator/=(const Data &right)
Overloaded operator /=.
Definition: Data.cpp:2497
Data pos() const
Return the identity of each data point of this Data object. Simply returns this object unmodified...
Definition: Data.cpp:1647
void resolve()
If this data is lazy, then convert it to ready data. What type of ready data depends on the expressio...
Definition: Data.cpp:945
Data powD(const Data &right) const
Return the given power of each data point of this Data object.
Definition: Data.cpp:2564
bool isReady() const
Return true if this data is ready.
Definition: Data.cpp:874
DataTypes::ValueType::reference getDataByTagRW(int tag, DataTypes::ValueType::size_type i)
getDataByTag
Definition: DataTagged.cpp:574
Data sin() const
Return the sin of each data point of this Data object.
Definition: Data.cpp:1503
void calc_minGlobalDataPoint(int &ProcNo, int &DataPointNo) const
Definition: Data.cpp:2254
double(* BinaryDFunPtr)(double, double)
Definition: Data.h:75
Data minval_nonlazy() const
Definition: Data.cpp:1963
const boost::python::tuple getShapeTuple() const
Return the data point shape as a tuple of integers.
Definition: Data.cpp:510
Data C_TensorUnaryOperation(Data const &arg_0, UnaryFunction operation)
Definition: Data.h:3197
bool isShared() const
Definition: Data.h:1846
Data acos() const
Return the acos of each data point of this Data object.
Definition: Data.cpp:1531
bool isProtected() const
Returns true, if the data object is protected against update.
Definition: Data.cpp:887
Data rtruedivO(const boost::python::object &left)
Newer style division operator for python.
Definition: Data.h:2205
Data interpolateFromTable1DP(boost::python::object table, double Amin, double Astep, double undef, bool check_boundaries)
Definition: Data.cpp:3347
std::map< int, int > DataMapType
Definition: DataTagged.h:57
std::string toString() const
Returns a text description of the function space.
Definition: FunctionSpace.cpp:122
Data log() const
Return the natural log of each data point of this Data object.
Definition: Data.cpp:1619
Data nonsymmetric() const
Return the nonsymmetric part of a matrix which is half the matrix minus its transpose.
Definition: Data.cpp:2079
Data interpolateFromTable3DP(boost::python::object table, double Amin, double Astep, Data &B, double Bmin, double Bstep, Data &C, double Cmin, double Cstep, double undef, bool check_boundaries)
Definition: Data.cpp:3328
void delaySelf()
convert the current data into lazy data.
Definition: Data.cpp:562
bool hasNoSamples() const
Return true if this object contains no samples. This is not the same as isEmpty() ...
Definition: Data.h:722
double rpow(double x, double y)
Definition: Data.h:2053
DataTypes::ValueType::size_type getDataOffset(int sampleNo, int dataPointNo)
Return the offset for the given sample and point within the sample.
Definition: Data.h:673
void setTaggedValueFromCPP(int tagKey, const DataTypes::ShapeType &pointshape, const DataTypes::ValueType &value, int dataOffset=0)
Assign the given value to the tag. Implicitly converts this object to type DataTagged if it is consta...
Definition: Data.cpp:2830
void expand()
Whatever the current Data type make this into a DataExpanded.
Definition: Data.cpp:895
const FunctionSpace & getFunctionSpace() const
Return the function space.
Definition: Data.h:455
Data operator+(const Data &left, const Data &right)
Operator+ Takes two Data objects.
Definition: Data.cpp:2574
boost::python::object integrateToTuple_const() const
Calculate the integral over the function space domain as a python tuple.
Definition: Data.cpp:1343
DataTypes::ValueType::size_type getLength() const
Return the number of doubles stored for this Data.
Definition: Data.cpp:1061
Data represents a collection of datapoints.
Definition: Data.h:68
void updateShareStatus(bool nowshared) const
Update the Data&#39;s shared flag This indicates that the DataAbstract used by this object is now shared ...
Definition: Data.h:1827
double algorithm(BinaryFunction operation, double initial_value) const
Perform the specified reduction algorithm on every element of every data point in this Data object ac...
Definition: Data.h:2306
const DataAbstract::ValueType::value_type * getSampleDataRO(DataAbstract::ValueType::size_type sampleNo) const
Return the sample data for the given sample no. Please do not use this unless you NEED to access samp...
Definition: Data.h:2019
Data sinh() const
Return the sinh of each data point of this Data object.
Definition: Data.cpp:1546
double supWorker() const
Definition: Data.cpp:1889
Data maxval() const
Return the maximum value of each data point of this Data object.
Definition: Data.cpp:1983
Data interpolateFromTable3D(const WrappedArray &table, double Amin, double Astep, double undef, Data &B, double Bmin, double Bstep, Data &C, double Cmin, double Cstep, bool check_boundaries)
Definition: Data.cpp:3643
void setValueOfDataPointToArray(int dataPointNo, const boost::python::object &)
sets the values of a data-point from a array-like object on this process
Definition: Data.cpp:1232
#define EsysAssert(AssertTest, AssertMessage)
EsysAssert is a MACRO that will throw an exception if the boolean condition specified is false...
Definition: EsysAssert.h:96
int get_MPISize(void) const
return the MPI rank number of the local data MPI_COMM_WORLD is assumed and the result of MPI_Comm_ran...
Definition: Data.cpp:4077
const ElementType & const_reference
Definition: DataVector.h:62
const double * getSampleDataRO(ValueType::size_type sampleNo) const
Definition: DataReady.h:126
DataAbstract::ValueType::value_type * getSampleDataRW(DataAbstract::ValueType::size_type sampleNo)
Return the sample data for the given sample no. Please do not use this unless you NEED to access samp...
Definition: Data.h:2002
boost::python::object integrateWorker() const
Definition: Data.cpp:1364
Data & operator=(const Data &other)
Definition: Data.cpp:2443
void forceResolve()
Definition: Data.h:1858
Data C_TensorBinaryOperation(Data const &arg_0, Data const &arg_1, BinaryFunction operation)
Compute a tensor operation with two Data objects.
Definition: Data.h:2388
int getTypeCode() const
Returns the function space type code.
Definition: FunctionSpace.cpp:99
int getNumDataPoints() const
Return the number of data points.
Definition: Data.h:507
void requireWrite()
Ensures data is ready for write access. This means that the data will be resolved if lazy and will be...
Definition: Data.cpp:954
virtual const DataTypes::ValueType & getVectorRO() const =0
Data whereNonPositive() const
Return a Data with a 1 for -ive or 0 values and a 0 for +ive values.
Definition: Data.cpp:989
Data atanh() const
Return the atanh of each data point of this Data object.
Definition: Data.cpp:1601
Data copySelf()
Return a pointer to a deep copy of this object.
Definition: Data.cpp:535
DataTypes::ValueType::const_reference getDataPointRO(int sampleNo, int dataPointNo)
Return a reference into the DataVector which points to the specified data point.
Definition: Data.cpp:3304
double algorithm(const DataExpanded &data, BinaryFunction operation, double initial_value)
Perform the given operation upon all values in all data-points in the given Data object and return th...
Definition: DataAlgorithm.h:186
DataTypes::ValueType::reference getDefaultValueRW(DataTypes::ValueType::size_type i)
getDefaultValue
Definition: DataTagged.h:615
bool m_lazy
Definition: Data.h:1796
bool isEmpty() const
Return true if this Data holds an instance of DataEmpty. This is _not_ the same as asking if the obj...
Definition: Data.cpp:853
void setTupleForGlobalDataPoint(int id, int proc, boost::python::object)
Set the value of a global data point.
Definition: Data.cpp:1188
DataTypes::ValueType::const_reference getDefaultValueRO(DataTypes::ValueType::size_type i) const
Definition: DataTagged.h:622
boost::python::object __div__(const boost::python::object &right)
wrapper for python divide operation
Definition: Data.cpp:4633
DataAbstract * borrowData(void) const
return the object produced by the factory, which is a DataConstant or DataExpanded TODO Ownership of ...
Definition: Data.cpp:3226
Data grad() const
Definition: Data.cpp:1044
const DataAbstract::ValueType::value_type * getDataRO() const
Return a pointer to the beginning of the underlying data.
Definition: Data.h:2033
void setValueOfDataPoint(int dataPointNo, const double)
sets the values of a data-point on this process
Definition: Data.cpp:1269
friend Data condEval(escript::Data &mask, escript::Data &trueval, escript::Data &falseval)
void calc_maxGlobalDataPoint(int &ProcNo, int &DataPointNo) const
Definition: Data.cpp:2340
void copy(const Data &other)
Make this object a deep copy of "other".
Definition: Data.cpp:542
Wraps an expression tree of other DataObjects. The data will be evaluated when required.
Definition: DataLazy.h:102
int MPI_Op
Definition: Esys_MPI.h:40
Data bessel(int order, double(*besselfunc)(int, double))
Bessel worker function.
Definition: Data.cpp:1413
int getTagNumber(int dpno)
Return the tag number associated with the given data-point.
Definition: Data.cpp:2849
double inf()
Return the minimum value of this Data object.
Definition: Data.cpp:1744
DataException exception class.
Definition: DataException.h:35
#define ESCRIPT_DLL_API
Definition: escriptcore/src/system_dep.h:54
Data nonuniforminterp(boost::python::object in, boost::python::object out, bool check_boundaries)
Definition: Data.cpp:3886
Data erf() const
Return the error function erf of each data point of this Data object.
Definition: Data.cpp:1568
Data acosh() const
Return the acosh of each data point of this Data object.
Definition: Data.cpp:1590
ElementType & reference
Definition: DataVector.h:61
Describes binary operations performed on instances of DataAbstract.
ElementType value_type
Definition: DataVector.h:59
Data operator*(const AbstractSystemMatrix &left, const Data &right)
Definition: AbstractSystemMatrix.cpp:43
Data sign() const
Return the sign of each data point of this Data object. -1 for negative values, zero for zero values...
Definition: Data.cpp:1626
void exclusiveWrite()
if another object is sharing out member data make a copy to work with instead. This code should only ...
Definition: Data.h:1876
Data oneOver() const
Returns 1./ Data object.
Definition: Data.cpp:961
Data abs() const
Return the absolute value of each data point of this Data object.
Definition: Data.cpp:1633
Definition: DataAbstract.h:61
const boost::python::tuple minGlobalDataPoint() const
Return the (sample number, data-point number) of the data point with the minimum component value in t...
Definition: Data.cpp:2241
const ValueType * resolveSample(int sampleNo, size_t &roffset) const
Compute the value of the expression for the given sample.
Definition: DataLazy.cpp:1644
bool actsExpanded() const
Return true if this Data is expanded or resolves to expanded. That is, if it has a separate value for...
Definition: Data.cpp:839
void setProtection()
switches on update protection
Definition: Data.cpp:881
Data besselFirstKind(int order)
Return the Bessel function of the first kind for each data point of this Data object.
Definition: Data.cpp:1401
boost::shared_ptr< DataReady > DataReady_ptr
Definition: DataAbstract.h:56
boost::python::object __rdiv__(const boost::python::object &right)
wrapper for python reverse divide operation
Definition: Data.cpp:4657
bool hasNaN()
returns return true if data contains NaN.
Definition: Data.cpp:1826
Data powO(const boost::python::object &right) const
Return the given power of each data point of this boost python object.
Definition: Data.cpp:2557
Data interpolateFromTable2D(const WrappedArray &table, double Amin, double Astep, double undef, Data &B, double Bmin, double Bstep, bool check_boundaries)
Definition: Data.cpp:3484
const boost::python::tuple maxGlobalDataPoint() const
Return the (sample number, data-point number) of the data point with the minimum component value in t...
Definition: Data.cpp:2331
Data eigenvalues() const
Return the eigenvalues of the symmetric part at each data point of this Data object in increasing val...
Definition: Data.cpp:2193
Data tan() const
Return the tan of each data point of this Data object.
Definition: Data.cpp:1517
Data operator/(const Data &left, const Data &right)
Operator/ Takes two Data objects.
Definition: Data.cpp:2601
Data()
Default constructor. Creates a DataEmpty object.
Definition: Data.cpp:241
void initialise(const DataTypes::ValueType &value, const DataTypes::ShapeType &shape, const FunctionSpace &what, bool expanded)
Construct a Data object of the appropriate type.
Definition: Data.cpp:468
void copyWithMask(const Data &other, const Data &mask)
Copy other Data object into this Data object where mask is positive.
Definition: Data.cpp:602
boost::python::object __mul__(const boost::python::object &right)
wrapper for python multiply operation
Definition: Data.cpp:4609
void typeMatchLeft(Data &right) const
Convert the data type of the RHS to match this.
Definition: Data.cpp:2752
long size_type
Definition: DataVector.h:60
int getNoValues() const
Return the number of values in the shape for this object.
Definition: Data.h:568
Data whereZero(double tol=0.0) const
Return a Data with a 1 for 0 values and a 0 for +ive or -ive values.
Definition: Data.cpp:996
Data asinh() const
Return the asinh of each data point of this Data object.
Definition: Data.cpp:1579
DataTypes::ValueType::reference getDataPointRW(int sampleNo, int dataPointNo)
Return a reference into the DataVector which points to the specified data point.
Definition: Data.cpp:3320
Data interpolate(const FunctionSpace &functionspace) const
Interpolates this onto the given functionspace and returns the result as a Data object.
Definition: Data.cpp:1014
Data truedivO(const boost::python::object &right)
Newer style division operator for python.
Definition: Data.h:2192
Definition: WrappedArray.h:29
double sup()
Return the maximum value of this Data object.
Definition: Data.cpp:1713
Data exp() const
Return the exponential function of each data point of this Data object.
Definition: Data.cpp:1658
~Data()
Destructor.
Definition: Data.cpp:425
void unaryOp2(UnaryFunction operation)
Perform the given unary operation on every element of every data point in this Data object...
void tensor_binary_operation(const int size, const double *arg1, const double *arg2, double *argRes, BinaryFunction operation)
Definition: LocalOps.h:524
const boost::python::object getValueOfDataPointAsTuple(int dataPointNo)
Return the value of a data point as a python tuple.
Definition: Data.cpp:1149
const_Domain_ptr getDomain() const
Returns the function space domain.
Definition: FunctionSpace.cpp:107
Data truedivD(const Data &right)
Newer style division operator for python.
Definition: Data.h:2180
boost::shared_ptr< const AbstractDomain > const_Domain_ptr
Definition: AbstractDomain.h:39
Domain_ptr getDomainPython() const
Return the domain. TODO: For internal use only. This should be removed.
Definition: Data.h:483
Data & operator-=(const Data &right)
Overloaded operator -=.
Definition: Data.cpp:2451
Data whereNonZero(double tol=0.0) const
Return a Data with a 0 for 0 values and a 1 for +ive or -ive values.
Definition: Data.cpp:1005
void setTaggedValue(int tagKey, const boost::python::object &value)
Assign the given value to the tag. Implicitly converts this object to type DataTagged if it is consta...
Definition: Data.cpp:2809
const boost::python::object getValueOfGlobalDataPointAsTuple(int procNo, int dataPointNo)
Return a data point across all processors as a python tuple.
Definition: Data.cpp:1291
void dump(const std::string fileName) const
dumps the object into a netCDF file
Definition: Data.cpp:4055
const DataReady * getReady() const
Definition: Data.h:1962
const boost::python::tuple eigenvalues_and_eigenvectors(const double tol=1.e-12) const
Return the eigenvalues and corresponding eigenvcetors of the symmetric part at each data point of thi...
Definition: Data.cpp:2216
Data whereNegative() const
Return a Data with a 1 for -ive values and a 0 for +ive or 0 values.
Definition: Data.cpp:975
bool isDataPointShapeEqual(int rank, const int *dimensions) const
Returns true if the shape matches the vector (dimensions[0],..., dimensions[rank-1]). DataEmpty always returns true.
Definition: Data.h:554
const_Domain_ptr getDomain() const
Return the domain.
Definition: Data.h:468
Data & operator+=(const Data &right)
Overloaded operator +=.
Definition: Data.cpp:2419
void binaryOp(DataTagged &left, const DataConstant &right, BinaryFunction operation)
Perform the given binary operation.
Definition: BinaryOp.h:45