escript  Revision_4925
Data.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2014 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
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 
33 
34 #include "DataC.h"
35 
36 #ifdef _OPENMP
37 #include <omp.h>
38 #endif
39 
40 #include "esysUtils/Esys_MPI.h"
41 #include <string>
42 #include <algorithm>
43 #include <sstream>
44 
45 #include <boost/shared_ptr.hpp>
46 #include <boost/python/object.hpp>
47 #include <boost/python/tuple.hpp>
48 
49 namespace escript {
50 
51 //
52 // Forward declaration for various implementations of Data.
53 class DataConstant;
54 class DataTagged;
55 class DataExpanded;
56 class DataLazy;
57 
71 class Data {
72 
73  public:
74 
75  // These typedefs allow function names to be cast to pointers
76  // to functions of the appropriate type when calling unaryOp etc.
77  typedef double (*UnaryDFunPtr)(double);
78  typedef double (*BinaryDFunPtr)(double,double);
79 
80 
91  Data();
92 
99  Data(const Data& inData);
100 
108  Data(const Data& inData,
109  const FunctionSpace& what);
110 
116  Data(const DataTypes::ValueType& value,
117  const DataTypes::ShapeType& shape,
118  const FunctionSpace& what=FunctionSpace(),
119  bool expanded=false);
120 
133  Data(double value,
134  const DataTypes::ShapeType& dataPointShape=DataTypes::ShapeType(),
135  const FunctionSpace& what=FunctionSpace(),
136  bool expanded=false);
137 
146  Data(const Data& inData,
147  const DataTypes::RegionType& region);
148 
160  Data(const boost::python::object& value,
161  const FunctionSpace& what=FunctionSpace(),
162  bool expanded=false);
163 
175  Data(const WrappedArray& w, const FunctionSpace& what,
176  bool expanded=false);
177 
178 
189  Data(const boost::python::object& value,
190  const Data& other);
191 
197  Data(double value,
198  const boost::python::tuple& shape=boost::python::make_tuple(),
199  const FunctionSpace& what=FunctionSpace(),
200  bool expanded=false);
201 
207  explicit Data(DataAbstract* underlyingdata);
208 
213  explicit Data(DataAbstract_ptr underlyingdata);
214 
220  ~Data();
221 
226  void
227  copy(const Data& other);
228 
233  Data
234  copySelf();
235 
236 
241  Data
242  delay();
243 
248  void
249  delaySelf();
250 
251 
262  void
263  setProtection();
264 
271  bool
272  isProtected() const;
273 
274 
280  const boost::python::object
281  getValueOfDataPointAsTuple(int dataPointNo);
282 
288  void
289  setValueOfDataPointToPyObject(int dataPointNo, const boost::python::object& py_object);
290 
296  void
297  setValueOfDataPointToArray(int dataPointNo, const boost::python::object&);
298 
304  void
305  setValueOfDataPoint(int dataPointNo, const double);
306 
311  const boost::python::object
312  getValueOfGlobalDataPointAsTuple(int procNo, int dataPointNo);
313 
314 
319  void
320  setTupleForGlobalDataPoint(int id, int proc, boost::python::object);
321 
328  int
329  getTagNumber(int dpno);
330 
337  getDataC();
338 
339 
340 
347  getDataC() const;
348 
349 
355  std::string
356  toString() const;
357 
363  void
364  expand();
365 
373  void
374  tag();
375 
381  void
382  resolve();
383 
384 
393  void
394  requireWrite();
395 
402  bool
403  isExpanded() const;
404 
411  bool
412  actsExpanded() const;
413 
414 
420  bool
421  isTagged() const;
422 
428  bool
429  isConstant() const;
430 
435  bool
436  isLazy() const;
437 
442  bool
443  isReady() const;
444 
451  bool
452  isEmpty() const;
453 
459  inline
460  const FunctionSpace&
462  {
463  return m_data->getFunctionSpace();
464  }
465 
471  inline
472 // const AbstractDomain&
474  getDomain() const
475  {
476  return getFunctionSpace().getDomain();
477  }
478 
479 
486  inline
487 // const AbstractDomain&
488  Domain_ptr
490  {
492  }
493 
499  inline
500  unsigned int
502  {
503  return m_data->getRank();
504  }
505 
511  inline
512  int
514  {
516  }
522  inline
523  int
525  {
526  return m_data->getNumSamples();
527  }
528 
534  inline
535  int
537  {
538  return m_data->getNumDPPSample();
539  }
540 
547  inline
548  bool numSamplesEqual(int numDataPointsPerSample, int numSamples) const
549  {
550  return (isEmpty() ||
551  (numDataPointsPerSample==getNumDataPointsPerSample() && numSamples==getNumSamples()));
552  }
553 
559  inline
560  bool isDataPointShapeEqual(int rank, const int* dimensions) const
561  {
562  if (isEmpty())
563  return true;
564  const DataTypes::ShapeType givenShape(&dimensions[0],&dimensions[rank]);
565  return (getDataPointShape()==givenShape);
566  }
567 
573  int
574  getNoValues() const
575  {
576  return m_data->getNoValues();
577  }
578 
579 
585  void
586  dump(const std::string fileName) const;
587 
595  const boost::python::object
596  toListOfTuples(bool scalarastuple=true);
597 
598 
607  inline
610 
611 
620  inline
623 
624 
632  inline
635  {
636  return m_data->getSampleDataByTag(tag);
637  }
638 
647  getDataPointRO(int sampleNo, int dataPointNo);
648 
657  getDataPointRW(int sampleNo, int dataPointNo);
658 
659 
660 
666  inline
668  getDataOffset(int sampleNo,
669  int dataPointNo)
670  {
671  return m_data->getPointOffset(sampleNo,dataPointNo);
672  }
673 
679  inline
680  const DataTypes::ShapeType&
682  {
683  return m_data->getShape();
684  }
685 
691  const boost::python::tuple
692  getShapeTuple() const;
693 
700  int
701  getDataPointSize() const;
702 
709  getLength() const;
710 
716  bool
717  hasNoSamples() const
718  {
719  return getLength()==0;
720  }
721 
731  void
732  setTaggedValueByName(std::string name,
733  const boost::python::object& value);
734 
745  void
746  setTaggedValue(int tagKey,
747  const boost::python::object& value);
748 
760  void
761  setTaggedValueFromCPP(int tagKey,
762  const DataTypes::ShapeType& pointshape,
763  const DataTypes::ValueType& value,
764  int dataOffset=0);
765 
766 
767 
773  void
774  copyWithMask(const Data& other,
775  const Data& mask);
776 
787  void
788  setToZero();
789 
797  Data
798  interpolate(const FunctionSpace& functionspace) const;
799 
801  Data
802  interpolateFromTable3D(const WrappedArray& table, double Amin, double Astep,
803  double undef, Data& B, double Bmin, double Bstep, Data& C,
804  double Cmin, double Cstep, bool check_boundaries);
805 
807  Data
808  interpolateFromTable2D(const WrappedArray& table, double Amin, double Astep,
809  double undef, Data& B, double Bmin, double Bstep,bool check_boundaries);
810 
812  Data
813  interpolateFromTable1D(const WrappedArray& table, double Amin, double Astep,
814  double undef,bool check_boundaries);
815 
816 
818  Data
819  interpolateFromTable3DP(boost::python::object table, double Amin, double Astep,
820  Data& B, double Bmin, double Bstep, Data& C, double Cmin, double Cstep, double undef,bool check_boundaries);
821 
822 
824  Data
825  interpolateFromTable2DP(boost::python::object table, double Amin, double Astep,
826  Data& B, double Bmin, double Bstep, double undef,bool check_boundaries);
827 
829  Data
830  interpolateFromTable1DP(boost::python::object table, double Amin, double Astep,
831  double undef,bool check_boundaries);
832 
834  Data
835  nonuniforminterp(boost::python::object in, boost::python::object out, bool check_boundaries);
836 
838  Data
839  nonuniformslope(boost::python::object in, boost::python::object out, bool check_boundaries);
840 
848  Data
849  gradOn(const FunctionSpace& functionspace) const;
850 
852  Data
853  grad() const;
854 
860  boost::python::object
861  integrateToTuple_const() const;
862 
863 
869  boost::python::object
871 
872 
873 
880  Data
881  oneOver() const;
888  Data
889  wherePositive() const;
890 
897  Data
898  whereNegative() const;
899 
906  Data
907  whereNonNegative() const;
908 
915  Data
916  whereNonPositive() const;
917 
924  Data
925  whereZero(double tol=0.0) const;
926 
933  Data
934  whereNonZero(double tol=0.0) const;
935 
948  double
949  Lsup();
950 
952  double
953  Lsup_const() const;
954 
955 
968  double
969  sup();
970 
972  double
973  sup_const() const;
974 
975 
988  double
989  inf();
990 
992  double
993  inf_const() const;
994 
995 
996 
1003  Data
1004  abs() const;
1005 
1012  Data
1013  maxval() const;
1014 
1021  Data
1022  minval() const;
1023 
1032  const boost::python::tuple
1033  minGlobalDataPoint() const;
1034 
1043  const boost::python::tuple
1044  maxGlobalDataPoint() const;
1045 
1046 
1047 
1055  Data
1056  sign() const;
1057 
1064  Data
1065  symmetric() const;
1066 
1073  Data
1074  nonsymmetric() const;
1075 
1082  Data
1083  trace(int axis_offset) const;
1084 
1091  Data
1092  transpose(int axis_offset) const;
1093 
1101  Data
1102  eigenvalues() const;
1103 
1114  const boost::python::tuple
1115  eigenvalues_and_eigenvectors(const double tol=1.e-12) const;
1116 
1123  Data
1124  swapaxes(const int axis0, const int axis1) const;
1125 
1132  Data
1133  erf() const;
1134 
1141  Data
1142  sin() const;
1143 
1150  Data
1151  cos() const;
1152 
1159  Data
1160  tan() const;
1161 
1168  Data
1169  asin() const;
1170 
1177  Data
1178  acos() const;
1179 
1186  Data
1187  atan() const;
1188 
1195  Data
1196  sinh() const;
1197 
1204  Data
1205  cosh() const;
1206 
1213  Data
1214  tanh() const;
1215 
1222  Data
1223  asinh() const;
1224 
1231  Data
1232  acosh() const;
1233 
1240  Data
1241  atanh() const;
1242 
1249  Data
1250  log10() const;
1251 
1258  Data
1259  log() const;
1260 
1267  Data
1268  exp() const;
1269 
1276  Data
1277  sqrt() const;
1278 
1285  Data
1286  neg() const;
1287 
1295  Data
1296  pos() const;
1297 
1306  Data
1307  powD(const Data& right) const;
1308 
1317  Data
1318  powO(const boost::python::object& right) const;
1319 
1329  Data
1330  rpowO(const boost::python::object& left) const;
1331 
1339  Data& operator+=(const Data& right);
1341  Data& operator+=(const boost::python::object& right);
1342 
1344  Data& operator=(const Data& other);
1345 
1353  Data& operator-=(const Data& right);
1355  Data& operator-=(const boost::python::object& right);
1356 
1364  Data& operator*=(const Data& right);
1366  Data& operator*=(const boost::python::object& right);
1367 
1375  Data& operator/=(const Data& right);
1377  Data& operator/=(const boost::python::object& right);
1378 
1384  Data truedivD(const Data& right);
1385 
1391  Data truedivO(const boost::python::object& right);
1392 
1398  Data rtruedivO(const boost::python::object& left);
1399 
1405  boost::python::object __add__(const boost::python::object& right);
1406 
1407 
1413  boost::python::object __sub__(const boost::python::object& right);
1414 
1420  boost::python::object __rsub__(const boost::python::object& right);
1421 
1427  boost::python::object __mul__(const boost::python::object& right);
1428 
1434  boost::python::object __div__(const boost::python::object& right);
1435 
1441  boost::python::object __rdiv__(const boost::python::object& right);
1442 
1447  Data
1448  matrixInverse() const;
1449 
1455  bool
1456  probeInterpolation(const FunctionSpace& functionspace) const;
1457 
1474  Data
1475  getItem(const boost::python::object& key) const;
1476 
1489  void
1490  setItemD(const boost::python::object& key,
1491  const Data& value);
1492 
1494  void
1495  setItemO(const boost::python::object& key,
1496  const boost::python::object& value);
1497 
1498  // These following public methods should be treated as private.
1499 
1505  template <class UnaryFunction>
1507  inline
1508  void
1509  unaryOp2(UnaryFunction operation);
1510 
1519  Data
1520  getSlice(const DataTypes::RegionType& region) const;
1521 
1531  void
1532  setSlice(const Data& value,
1533  const DataTypes::RegionType& region);
1534 
1540  void
1541  print(void);
1542 
1550  int
1551  get_MPIRank(void) const;
1552 
1560  int
1561  get_MPISize(void) const;
1562 
1569  MPI_Comm
1570  get_MPIComm(void) const;
1571 
1578  DataAbstract*
1579  borrowData(void) const;
1580 
1583  borrowDataPtr(void) const;
1584 
1587  borrowReadyPtr(void) const;
1588 
1589 
1590 
1601 
1602 
1606 
1619 
1620  protected:
1621 
1622  private:
1623 
1624 template <class BinaryOp>
1625  double
1626 #ifdef ESYS_MPI
1627  lazyAlgWorker(double init, MPI_Op mpiop_type);
1628 #else
1629  lazyAlgWorker(double init);
1630 #endif
1631 
1632  double
1633  LsupWorker() const;
1634 
1635  double
1636  supWorker() const;
1637 
1638  double
1639  infWorker() const;
1640 
1641  boost::python::object
1642  integrateWorker() const;
1643 
1644  void
1645  calc_minGlobalDataPoint(int& ProcNo, int& DataPointNo) const;
1646 
1647  void
1648  calc_maxGlobalDataPoint(int& ProcNo, int& DataPointNo) const;
1649 
1650  // For internal use in Data.cpp only!
1651  // other uses should call the main entry points and allow laziness
1652  Data
1653  minval_nonlazy() const;
1654 
1655  // For internal use in Data.cpp only!
1656  Data
1657  maxval_nonlazy() const;
1658 
1659 
1666  inline
1667  void
1668  operandCheck(const Data& right) const
1669  {
1670  return m_data->operandCheck(*(right.m_data.get()));
1671  }
1672 
1678  template <class BinaryFunction>
1679  inline
1680  double
1681  algorithm(BinaryFunction operation,
1682  double initial_value) const;
1683 
1691  template <class BinaryFunction>
1692  inline
1693  Data
1694  dp_algorithm(BinaryFunction operation,
1695  double initial_value) const;
1696 
1705  template <class BinaryFunction>
1706  inline
1707  void
1708  binaryOp(const Data& right,
1709  BinaryFunction operation);
1710 
1716  void
1717  typeMatchLeft(Data& right) const;
1718 
1724  void
1725  typeMatchRight(const Data& right);
1726 
1732  void
1733  initialise(const DataTypes::ValueType& value,
1734  const DataTypes::ShapeType& shape,
1735  const FunctionSpace& what,
1736  bool expanded);
1737 
1738  void
1739  initialise(const WrappedArray& value,
1740  const FunctionSpace& what,
1741  bool expanded);
1742 
1743  void
1744  initialise(const double value,
1745  const DataTypes::ShapeType& shape,
1746  const FunctionSpace& what,
1747  bool expanded);
1748 
1749  //
1750  // flag to protect the data object against any update
1752  mutable bool m_shared;
1753  bool m_lazy;
1754 
1755  //
1756  // pointer to the actual data object
1757 // boost::shared_ptr<DataAbstract> m_data;
1759 
1760 // If possible please use getReadyPtr instead.
1761 // But see warning below.
1762  const DataReady*
1763  getReady() const;
1764 
1765  DataReady*
1766  getReady();
1767 
1768 
1769 // Be wary of using this for local operations since it (temporarily) increases reference count.
1770 // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1771 // getReady() instead
1773  getReadyPtr();
1774 
1776  getReadyPtr() const;
1777 
1778 
1784  void updateShareStatus(bool nowshared) const
1785  {
1786  m_shared=nowshared; // m_shared is mutable
1787  }
1788 
1789  // In the isShared() method below:
1790  // A problem would occur if m_data (the address pointed to) were being modified
1791  // while the call m_data->is_shared is being executed.
1792  //
1793  // Q: So why do I think this code can be thread safe/correct?
1794  // A: We need to make some assumptions.
1795  // 1. We assume it is acceptable to return true under some conditions when we aren't shared.
1796  // 2. We assume that no constructions or assignments which will share previously unshared
1797  // will occur while this call is executing. This is consistent with the way Data:: and C are written.
1798  //
1799  // This means that the only transition we need to consider, is when a previously shared object is
1800  // not shared anymore. ie. the other objects have been destroyed or a deep copy has been made.
1801  // In those cases the m_shared flag changes to false after m_data has completed changing.
1802  // For any threads executing before the flag switches they will assume the object is still shared.
1803  bool isShared() const
1804  {
1805  return m_shared;
1806 /* if (m_shared) return true;
1807  if (m_data->isShared())
1808  {
1809  updateShareStatus(true);
1810  return true;
1811  }
1812  return false;*/
1813  }
1814 
1816  {
1817  if (isLazy())
1818  {
1819  #ifdef _OPENMP
1820  if (omp_in_parallel())
1821  { // Yes this is throwing an exception out of an omp thread which is forbidden.
1822  throw DataException("Please do not call forceResolve() in a parallel region.");
1823  }
1824  #endif
1825  resolve();
1826  }
1827  }
1828 
1834  {
1835 #ifdef _OPENMP
1836  if (omp_in_parallel())
1837  {
1838 // *((int*)0)=17;
1839  throw DataException("Programming error. Please do not run exclusiveWrite() in multi-threaded sections.");
1840  }
1841 #endif
1842  forceResolve();
1843  if (isShared())
1844  {
1845  DataAbstract* t=m_data->deepCopy();
1847  }
1848  }
1849 
1854  {
1855  if (isLazy() || isShared())
1856  {
1857  throw DataException("Programming error. ExclusiveWrite required - please call requireWrite()");
1858  }
1859  }
1860 
1867  void set_m_data(DataAbstract_ptr p);
1868 
1869  friend class DataAbstract; // To allow calls to updateShareStatus
1870  friend class TestDomain; // so its getX will work quickly
1871 #ifdef IKNOWWHATIMDOING
1872  friend ESCRIPT_DLL_API Data applyBinaryCFunction(boost::python::object cfunc, boost::python::tuple shape, escript::Data& d, escript::Data& e);
1873 #endif
1874  friend ESCRIPT_DLL_API Data condEval(escript::Data& mask, escript::Data& trueval, escript::Data& falseval);
1875  friend ESCRIPT_DLL_API Data randomData(const boost::python::tuple& shape, const FunctionSpace& what, long seed, const boost::python::tuple& filter);
1876 
1877 };
1878 
1879 
1880 #ifdef IKNOWWHATIMDOING
1882 Data
1883 applyBinaryCFunction(boost::python::object func, boost::python::tuple shape, escript::Data& d, escript::Data& e);
1884 #endif
1885 
1886 
1888 Data
1889 condEval(escript::Data& mask, escript::Data& trueval, escript::Data& falseval);
1890 
1891 
1892 
1897 Data randomData(const boost::python::tuple& shape,
1898  const FunctionSpace& what,
1899  long seed, const boost::python::tuple& filter);
1900 
1901 
1902 } // end namespace escript
1903 
1904 
1905 // No, this is not supposed to be at the top of the file
1906 // DataAbstact needs to be declared first, then DataReady needs to be fully declared
1907 // so that I can dynamic cast between them below.
1908 #include "DataReady.h"
1909 #include "DataLazy.h"
1910 
1911 namespace escript
1912 {
1913 
1914 inline
1915 const DataReady*
1917 {
1918  const DataReady* dr=dynamic_cast<const DataReady*>(m_data.get());
1919  EsysAssert((dr!=0), "Error - casting to DataReady.");
1920  return dr;
1921 }
1922 
1923 inline
1924 DataReady*
1926 {
1927  DataReady* dr=dynamic_cast<DataReady*>(m_data.get());
1928  EsysAssert((dr!=0), "Error - casting to DataReady.");
1929  return dr;
1930 }
1931 
1932 // Be wary of using this for local operations since it (temporarily) increases reference count.
1933 // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1934 // getReady() instead
1935 inline
1938 {
1939  DataReady_ptr dr=boost::dynamic_pointer_cast<DataReady>(m_data);
1940  EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1941  return dr;
1942 }
1943 
1944 
1945 inline
1948 {
1949  const_DataReady_ptr dr=boost::dynamic_pointer_cast<const DataReady>(m_data);
1950  EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1951  return dr;
1952 }
1953 
1954 inline
1957 {
1958  if (isLazy())
1959  {
1960  throw DataException("Error, attempt to acquire RW access to lazy data. Please call requireWrite() first.");
1961  }
1962  return getReady()->getSampleDataRW(sampleNo);
1963 }
1964 
1965 inline
1968 {
1969  DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
1970  if (l!=0)
1971  {
1972  size_t offset=0;
1973  const DataTypes::ValueType* res=l->resolveSample(sampleNo,offset);
1974  return &((*res)[offset]);
1975  }
1976  return getReady()->getSampleDataRO(sampleNo);
1977 }
1978 
1982 inline double rpow(double x,double y)
1983 {
1984  return pow(y,x);
1985 }
1986 
1992 ESCRIPT_DLL_API Data operator+(const Data& left, const Data& right);
1993 
1999 ESCRIPT_DLL_API Data operator-(const Data& left, const Data& right);
2000 
2006 ESCRIPT_DLL_API Data operator*(const Data& left, const Data& right);
2007 
2013 ESCRIPT_DLL_API Data operator/(const Data& left, const Data& right);
2014 
2021 ESCRIPT_DLL_API Data operator+(const Data& left, const boost::python::object& right);
2022 
2029 ESCRIPT_DLL_API Data operator-(const Data& left, const boost::python::object& right);
2030 
2037 ESCRIPT_DLL_API Data operator*(const Data& left, const boost::python::object& right);
2038 
2045 ESCRIPT_DLL_API Data operator/(const Data& left, const boost::python::object& right);
2046 
2053 ESCRIPT_DLL_API Data operator+(const boost::python::object& left, const Data& right);
2054 
2061 ESCRIPT_DLL_API Data operator-(const boost::python::object& left, const Data& right);
2062 
2069 ESCRIPT_DLL_API Data operator*(const boost::python::object& left, const Data& right);
2070 
2077 ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);
2078 
2079 
2080 
2085 ESCRIPT_DLL_API std::ostream& operator<<(std::ostream& o, const Data& data);
2086 
2096 Data
2097 C_GeneralTensorProduct(Data& arg_0,
2098  Data& arg_1,
2099  int axis_offset=0,
2100  int transpose=0);
2101 
2107 inline
2108 Data
2109 Data::truedivD(const Data& right)
2110 {
2111  return *this / right;
2112 }
2113 
2119 inline
2120 Data
2121 Data::truedivO(const boost::python::object& right)
2122 {
2123  Data tmp(right, getFunctionSpace(), false);
2124  return truedivD(tmp);
2125 }
2126 
2132 inline
2133 Data
2134 Data::rtruedivO(const boost::python::object& left)
2135 {
2136  Data tmp(left, getFunctionSpace(), false);
2137  return tmp.truedivD(*this);
2138 }
2139 
2145 template <class BinaryFunction>
2146 inline
2147 void
2148 Data::binaryOp(const Data& right,
2149  BinaryFunction operation)
2150 {
2151  //
2152  // if this has a rank of zero promote it to the rank of the RHS
2153  if (getDataPointRank()==0 && right.getDataPointRank()!=0) {
2154  throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");
2155  }
2156 
2157  if (isLazy() || right.isLazy())
2158  {
2159  throw DataException("Programmer error - attempt to call binaryOp with Lazy Data.");
2160  }
2161  //
2162  // initially make the temporary a shallow copy
2163  Data tempRight(right);
2165  FunctionSpace fsr=right.getFunctionSpace();
2166  if (fsl!=fsr) {
2167  signed char intres=fsl.getDomain()->preferredInterpolationOnDomain(fsr.getTypeCode(), fsl.getTypeCode());
2168  if (intres==0)
2169  {
2170  std::string msg="Error - attempt to combine incompatible FunctionSpaces.";
2171  msg+=fsl.toString();
2172  msg+=" ";
2173  msg+=fsr.toString();
2174  throw DataException(msg.c_str());
2175  }
2176  else if (intres==1)
2177  {
2178  // an interpolation is required so create a new Data
2179  tempRight=Data(right,fsl);
2180  }
2181  else // reverse interpolation preferred
2182  {
2183  // interpolate onto the RHS function space
2184  Data tempLeft(*this,fsr);
2185  set_m_data(tempLeft.m_data);
2186  }
2187  }
2188  operandCheck(tempRight);
2189  //
2190  // ensure this has the right type for the RHS
2191  typeMatchRight(tempRight);
2192  //
2193  // Need to cast to the concrete types so that the correct binaryOp
2194  // is called.
2195  if (isExpanded()) {
2196  //
2197  // Expanded data will be done in parallel, the right hand side can be
2198  // of any data type
2199  DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2200  EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2201  escript::binaryOp(*leftC,*(tempRight.getReady()),operation);
2202  } else if (isTagged()) {
2203  //
2204  // Tagged data is operated on serially, the right hand side can be
2205  // either DataConstant or DataTagged
2206  DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
2207  EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
2208  if (right.isTagged()) {
2209  DataTagged* rightC=dynamic_cast<DataTagged*>(tempRight.m_data.get());
2210  EsysAssert((rightC!=0), "Programming error - casting to DataTagged.");
2211  escript::binaryOp(*leftC,*rightC,operation);
2212  } else {
2213  DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
2214  EsysAssert((rightC!=0), "Programming error - casting to DataConstant.");
2215  escript::binaryOp(*leftC,*rightC,operation);
2216  }
2217  } else if (isConstant()) {
2218  DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2219  DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
2220  EsysAssert((leftC!=0 && rightC!=0), "Programming error - casting to DataConstant.");
2221  escript::binaryOp(*leftC,*rightC,operation);
2222  }
2223 }
2224 
2232 template <class BinaryFunction>
2233 inline
2234 double
2235 Data::algorithm(BinaryFunction operation, double initial_value) const
2236 {
2237  if (isExpanded()) {
2238  DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2239  EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2240  return escript::algorithm(*leftC,operation,initial_value);
2241  } else if (isTagged()) {
2242  DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
2243  EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
2244  return escript::algorithm(*leftC,operation,initial_value);
2245  } else if (isConstant()) {
2246  DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2247  EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
2248  return escript::algorithm(*leftC,operation,initial_value);
2249  } else if (isEmpty()) {
2250  throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2251  } else if (isLazy()) {
2252  throw DataException("Error - Operations not permitted on instances of DataLazy.");
2253  } else {
2254  throw DataException("Error - Data encapsulates an unknown type.");
2255  }
2256 }
2257 
2266 template <class BinaryFunction>
2267 inline
2268 Data
2269 Data::dp_algorithm(BinaryFunction operation, double initial_value) const
2270 {
2271  if (isEmpty()) {
2272  throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2273  }
2274  else if (isExpanded()) {
2276  DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
2277  DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
2278  EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
2279  EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
2280  escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
2281  return result;
2282  }
2283  else if (isTagged()) {
2284  DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
2285  EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
2286  DataTypes::ValueType defval(1);
2287  defval[0]=0;
2288  DataTagged* resultT=new DataTagged(getFunctionSpace(), DataTypes::scalarShape, defval, dataT);
2289  escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
2290  return Data(resultT); // note: the Data object now owns the resultT pointer
2291  }
2292  else if (isConstant()) {
2294  DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
2295  DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
2296  EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
2297  EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
2298  escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
2299  return result;
2300  } else if (isLazy()) {
2301  throw DataException("Error - Operations not permitted on instances of DataLazy.");
2302  } else {
2303  throw DataException("Error - Data encapsulates an unknown type.");
2304  }
2305 }
2306 
2314 template <typename BinaryFunction>
2315 inline
2316 Data
2318  Data const &arg_1,
2319  BinaryFunction operation)
2320 {
2321  if (arg_0.isEmpty() || arg_1.isEmpty())
2322  {
2323  throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2324  }
2325  if (arg_0.isLazy() || arg_1.isLazy())
2326  {
2327  throw DataException("Error - Operations not permitted on lazy data.");
2328  }
2329  // Interpolate if necessary and find an appropriate function space
2330  Data arg_0_Z, arg_1_Z;
2331  FunctionSpace fsl=arg_0.getFunctionSpace();
2332  FunctionSpace fsr=arg_1.getFunctionSpace();
2333  if (fsl!=fsr) {
2334  signed char intres=fsl.getDomain()->preferredInterpolationOnDomain(fsr.getTypeCode(), fsl.getTypeCode());
2335  if (intres==0)
2336  {
2337  std::string msg="Error - C_TensorBinaryOperation: arguments have incompatible function spaces.";
2338  msg+=fsl.toString();
2339  msg+=" ";
2340  msg+=fsr.toString();
2341  throw DataException(msg.c_str());
2342  }
2343  else if (intres==1)
2344  {
2345  arg_1_Z=arg_1.interpolate(arg_0.getFunctionSpace());
2346  arg_0_Z =Data(arg_0);
2347  }
2348  else // reverse interpolation preferred
2349  {
2350  arg_0_Z = arg_0.interpolate(arg_1.getFunctionSpace());
2351  arg_1_Z = Data(arg_1);
2352  }
2353  } else {
2354  arg_0_Z = Data(arg_0);
2355  arg_1_Z = Data(arg_1);
2356  }
2357  // Get rank and shape of inputs
2358  int rank0 = arg_0_Z.getDataPointRank();
2359  int rank1 = arg_1_Z.getDataPointRank();
2360  DataTypes::ShapeType shape0 = arg_0_Z.getDataPointShape();
2361  DataTypes::ShapeType shape1 = arg_1_Z.getDataPointShape();
2362  int size0 = arg_0_Z.getDataPointSize();
2363  int size1 = arg_1_Z.getDataPointSize();
2364  // Declare output Data object
2365  Data res;
2366 
2367  if (shape0 == shape1) {
2368  if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2369  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
2370  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2371  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2372  double *ptr_2 = &(res.getDataAtOffsetRW(0));
2373 
2374  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2375  }
2376  else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2377 
2378  // Prepare the DataConstant input
2379  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2380 
2381  // Borrow DataTagged input from Data object
2382  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2383 
2384  // Prepare a DataTagged output 2
2385  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
2386  res.tag();
2387  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2388 
2389  // Prepare offset into DataConstant
2390  int offset_0 = tmp_0->getPointOffset(0,0);
2391  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2392 
2393  // Get the pointers to the actual data
2394  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2395  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2396 
2397  // Compute a result for the default
2398  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2399  // Compute a result for each tag
2400  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2401  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2402  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2403  tmp_2->addTag(i->first);
2404  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2405  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2406 
2407  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2408  }
2409 
2410  }
2411  else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2412  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2413  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2414  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2415  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2416 
2417  int sampleNo_1,dataPointNo_1;
2418  int numSamples_1 = arg_1_Z.getNumSamples();
2419  int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2420  int offset_0 = tmp_0->getPointOffset(0,0);
2421  res.requireWrite();
2422  #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2423  for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2424  for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2425  int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2426  int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2427  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2428  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2429  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2430  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2431  }
2432  }
2433 
2434  }
2435  else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2436  // Borrow DataTagged input from Data object
2437  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2438 
2439  // Prepare the DataConstant input
2440  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2441 
2442  // Prepare a DataTagged output 2
2443  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
2444  res.tag();
2445  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2446 
2447  // Prepare offset into DataConstant
2448  int offset_1 = tmp_1->getPointOffset(0,0);
2449 
2450  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2451  // Get the pointers to the actual data
2452  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2453  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2454  // Compute a result for the default
2455  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2456  // Compute a result for each tag
2457  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2458  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2459  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2460  tmp_2->addTag(i->first);
2461  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2462  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2463  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2464  }
2465 
2466  }
2467  else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2468  // Borrow DataTagged input from Data object
2469  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2470 
2471  // Borrow DataTagged input from Data object
2472  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2473 
2474  // Prepare a DataTagged output 2
2475  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2476  res.tag(); // DataTagged output
2477  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2478 
2479  // Get the pointers to the actual data
2480  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2481  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2482  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2483 
2484  // Compute a result for the default
2485  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2486  // Merge the tags
2487  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2488  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2489  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2490  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2491  tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2492  }
2493  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2494  tmp_2->addTag(i->first);
2495  }
2496  // Compute a result for each tag
2497  const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2498  for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2499 
2500  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2501  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2502  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2503 
2504  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2505  }
2506 
2507  }
2508  else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2509  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2510  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2511  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2512  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2513  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2514 
2515  int sampleNo_0,dataPointNo_0;
2516  int numSamples_0 = arg_0_Z.getNumSamples();
2517  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2518  res.requireWrite();
2519  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2520  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2521  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2522  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2523  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2524  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2525  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2526  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2527  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2528  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2529  }
2530  }
2531 
2532  }
2533  else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2534  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2535  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2536  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2537  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2538 
2539  int sampleNo_0,dataPointNo_0;
2540  int numSamples_0 = arg_0_Z.getNumSamples();
2541  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2542  int offset_1 = tmp_1->getPointOffset(0,0);
2543  res.requireWrite();
2544  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2545  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2546  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2547  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2548  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2549 
2550  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2551  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2552  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2553 
2554 
2555  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2556  }
2557  }
2558 
2559  }
2560  else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2561  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2562  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2563  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2564  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2565  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2566 
2567  int sampleNo_0,dataPointNo_0;
2568  int numSamples_0 = arg_0_Z.getNumSamples();
2569  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2570  res.requireWrite();
2571  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2572  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2573  int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2574  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2575  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2576  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2577  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2578  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2579  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2580  tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2581  }
2582  }
2583 
2584  }
2585  else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2586  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2587  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2588  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2589  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2590  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2591 
2592  int sampleNo_0,dataPointNo_0;
2593  int numSamples_0 = arg_0_Z.getNumSamples();
2594  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2595  res.requireWrite();
2596  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2597  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2598  dataPointNo_0=0;
2599 // for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2600  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2601  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2602  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2603  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2604  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2605  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2606  tensor_binary_operation(size0*numDataPointsPerSample_0, ptr_0, ptr_1, ptr_2, operation);
2607 // }
2608  }
2609 
2610  }
2611  else {
2612  throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2613  }
2614 
2615  } else if (0 == rank0) {
2616  if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2617  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataConstant output
2618  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2619  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2620  double *ptr_2 = &(res.getDataAtOffsetRW(0));
2621  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2622  }
2623  else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2624 
2625  // Prepare the DataConstant input
2626  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2627 
2628  // Borrow DataTagged input from Data object
2629  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2630 
2631  // Prepare a DataTagged output 2
2632  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataTagged output
2633  res.tag();
2634  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2635 
2636  // Prepare offset into DataConstant
2637  int offset_0 = tmp_0->getPointOffset(0,0);
2638  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2639 
2640  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2641  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2642 
2643  // Compute a result for the default
2644  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2645  // Compute a result for each tag
2646  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2647  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2648  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2649  tmp_2->addTag(i->first);
2650  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2651  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2652  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2653  }
2654 
2655  }
2656  else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2657 
2658  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2659  DataConstant* tmp_0=dynamic_cast<DataConstant*>(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_1;
2664  int numSamples_1 = arg_1_Z.getNumSamples();
2665  int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2666  int offset_0 = tmp_0->getPointOffset(0,0);
2667  const double *ptr_src = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2668  double ptr_0 = ptr_src[0];
2669  int size = size1*numDataPointsPerSample_1;
2670  res.requireWrite();
2671  #pragma omp parallel for private(sampleNo_1) schedule(static)
2672  for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2673 // for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2674  int offset_1 = tmp_1->getPointOffset(sampleNo_1,0);
2675  int offset_2 = tmp_2->getPointOffset(sampleNo_1,0);
2676 // const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2677  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2678  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2679  tensor_binary_operation(size, ptr_0, ptr_1, ptr_2, operation);
2680 
2681 // }
2682  }
2683 
2684  }
2685  else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2686 
2687  // Borrow DataTagged input from Data object
2688  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2689 
2690  // Prepare the DataConstant input
2691  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2692 
2693  // Prepare a DataTagged output 2
2694  res = Data(0.0, shape1, arg_0_Z.getFunctionSpace()); // DataTagged output
2695  res.tag();
2696  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2697 
2698  // Prepare offset into DataConstant
2699  int offset_1 = tmp_1->getPointOffset(0,0);
2700  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2701 
2702  // Get the pointers to the actual data
2703  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2704  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2705 
2706 
2707  // Compute a result for the default
2708  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2709  // Compute a result for each tag
2710  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2711  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2712  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2713  tmp_2->addTag(i->first);
2714  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2715  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2716 
2717  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2718  }
2719 
2720  }
2721  else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2722 
2723  // Borrow DataTagged input from Data object
2724  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2725 
2726  // Borrow DataTagged input from Data object
2727  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2728 
2729  // Prepare a DataTagged output 2
2730  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());
2731  res.tag(); // DataTagged output
2732  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2733 
2734  // Get the pointers to the actual data
2735  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2736  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2737  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2738 
2739  // Compute a result for the default
2740  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2741  // Merge the tags
2742  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2743  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2744  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2745  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2746  tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2747  }
2748  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2749  tmp_2->addTag(i->first);
2750  }
2751  // Compute a result for each tag
2752  const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2753  for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2754  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2755  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2756  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2757 
2758  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2759  }
2760 
2761  }
2762  else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2763 
2764  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2765  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2766  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2767  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2768  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2769 
2770  int sampleNo_0,dataPointNo_0;
2771  int numSamples_0 = arg_0_Z.getNumSamples();
2772  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2773  res.requireWrite();
2774  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2775  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2776  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2777  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2778  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2779  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2780  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2781  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2782  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2783  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2784  }
2785  }
2786 
2787  }
2788  else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2789  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2790  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2791  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2792  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2793 
2794  int sampleNo_0,dataPointNo_0;
2795  int numSamples_0 = arg_0_Z.getNumSamples();
2796  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2797  int offset_1 = tmp_1->getPointOffset(0,0);
2798  res.requireWrite();
2799  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2800  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2801  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2802  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2803  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2804  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2805  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2806  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2807  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2808  }
2809  }
2810 
2811 
2812  }
2813  else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2814 
2815  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2816  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2817  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2818  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2819  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2820 
2821  int sampleNo_0,dataPointNo_0;
2822  int numSamples_0 = arg_0_Z.getNumSamples();
2823  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2824  res.requireWrite();
2825  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2826  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2827  int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2828  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2829  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2830  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2831  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2832  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2833  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2834  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2835  }
2836  }
2837 
2838  }
2839  else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2840 
2841  // After finding a common function space above the two inputs have the same numSamples and num DPPS
2842  res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2843  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2844  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2845  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2846 
2847  int sampleNo_0,dataPointNo_0;
2848  int numSamples_0 = arg_0_Z.getNumSamples();
2849  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2850  res.requireWrite();
2851  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2852  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2853  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2854  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2855  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2856  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2857  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2858  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2859  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2860  tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2861  }
2862  }
2863 
2864  }
2865  else {
2866  throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2867  }
2868 
2869  } else if (0 == rank1) {
2870  if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2871  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
2872  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2873  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2874  double *ptr_2 = &(res.getDataAtOffsetRW(0));
2875  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2876  }
2877  else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2878 
2879  // Prepare the DataConstant input
2880  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2881 
2882  // Borrow DataTagged input from Data object
2883  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2884 
2885  // Prepare a DataTagged output 2
2886  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
2887  res.tag();
2888  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2889 
2890  // Prepare offset into DataConstant
2891  int offset_0 = tmp_0->getPointOffset(0,0);
2892  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2893 
2894  //Get the pointers to the actual data
2895  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2896  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2897 
2898  // Compute a result for the default
2899  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2900  // Compute a result for each tag
2901  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2902  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2903  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2904  tmp_2->addTag(i->first);
2905  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2906  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2907  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2908  }
2909  }
2910  else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2911 
2912  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2913  DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2914  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2915  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2916 
2917  int sampleNo_1,dataPointNo_1;
2918  int numSamples_1 = arg_1_Z.getNumSamples();
2919  int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2920  int offset_0 = tmp_0->getPointOffset(0,0);
2921  res.requireWrite();
2922  #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2923  for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2924  for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2925  int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2926  int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2927  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2928  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2929  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2930  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2931  }
2932  }
2933 
2934  }
2935  else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2936 
2937  // Borrow DataTagged input from Data object
2938  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2939 
2940  // Prepare the DataConstant input
2941  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2942 
2943  // Prepare a DataTagged output 2
2944  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
2945  res.tag();
2946  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2947 
2948  // Prepare offset into DataConstant
2949  int offset_1 = tmp_1->getPointOffset(0,0);
2950  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2951  // Get the pointers to the actual data
2952  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2953  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2954  // Compute a result for the default
2955  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2956  // Compute a result for each tag
2957  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2958  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2959  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2960  tmp_2->addTag(i->first);
2961  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2962  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2963  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2964  }
2965 
2966  }
2967  else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2968 
2969  // Borrow DataTagged input from Data object
2970  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2971 
2972  // Borrow DataTagged input from Data object
2973  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2974 
2975  // Prepare a DataTagged output 2
2976  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2977  res.tag(); // DataTagged output
2978  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2979 
2980  // Get the pointers to the actual data
2981  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2982  const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2983  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2984 
2985  // Compute a result for the default
2986  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2987  // Merge the tags
2988  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2989  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2990  const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2991  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2992  tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2993  }
2994  for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2995  tmp_2->addTag(i->first);
2996  }
2997  // Compute a result for each tag
2998  const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2999  for (i=lookup_2.begin();i!=lookup_2.end();i++) {
3000  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
3001  const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
3002  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
3003  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3004  }
3005 
3006  }
3007  else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
3008 
3009  // After finding a common function space above the two inputs have the same numSamples and num DPPS
3010  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3011  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
3012  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
3013  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3014 
3015  int sampleNo_0,dataPointNo_0;
3016  int numSamples_0 = arg_0_Z.getNumSamples();
3017  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3018  res.requireWrite();
3019  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3020  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3021  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
3022  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3023  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3024  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
3025  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3026  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3027  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3028  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3029  }
3030  }
3031 
3032  }
3033  else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
3034  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3035  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3036  DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
3037  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3038 
3039  int sampleNo_0;
3040  int numSamples_0 = arg_0_Z.getNumSamples();
3041  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3042  int offset_1 = tmp_1->getPointOffset(0,0);
3043  const double *ptr_src = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3044  double ptr_1 = ptr_src[0];
3045  int size = size0 * numDataPointsPerSample_0;
3046  res.requireWrite();
3047  #pragma omp parallel for private(sampleNo_0) schedule(static)
3048  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3049 // for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3050  int offset_0 = tmp_0->getPointOffset(sampleNo_0,0);
3051  int offset_2 = tmp_2->getPointOffset(sampleNo_0,0);
3052  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3053 // const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3054  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3055  tensor_binary_operation(size, ptr_0, ptr_1, ptr_2, operation);
3056 // }
3057  }
3058 
3059 
3060  }
3061  else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
3062 
3063  // After finding a common function space above the two inputs have the same numSamples and num DPPS
3064  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3065  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3066  DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
3067  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3068 
3069  int sampleNo_0,dataPointNo_0;
3070  int numSamples_0 = arg_0_Z.getNumSamples();
3071  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3072  res.requireWrite();
3073  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3074  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3075  int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
3076  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3077  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3078  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
3079  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3080  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3081  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3082  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3083  }
3084  }
3085 
3086  }
3087  else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
3088 
3089  // After finding a common function space above the two inputs have the same numSamples and num DPPS
3090  res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
3091  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3092  DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
3093  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3094 
3095  int sampleNo_0,dataPointNo_0;
3096  int numSamples_0 = arg_0_Z.getNumSamples();
3097  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3098  res.requireWrite();
3099  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3100  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3101  for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3102  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
3103  int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
3104  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3105  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3106  const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
3107  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3108  tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
3109  }
3110  }
3111 
3112  }
3113  else {
3114  throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
3115  }
3116 
3117  } else {
3118  throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible shapes");
3119  }
3120 
3121  return res;
3122 }
3123 
3124 template <typename UnaryFunction>
3125 Data
3127  UnaryFunction operation)
3128 {
3129  if (arg_0.isEmpty()) // do this before we attempt to interpolate
3130  {
3131  throw DataException("Error - Operations not permitted on instances of DataEmpty.");
3132  }
3133  if (arg_0.isLazy())
3134  {
3135  throw DataException("Error - Operations not permitted on lazy data.");
3136  }
3137  // Interpolate if necessary and find an appropriate function space
3138  Data arg_0_Z = Data(arg_0);
3139 
3140  // Get rank and shape of inputs
3141  const DataTypes::ShapeType& shape0 = arg_0_Z.getDataPointShape();
3142  int size0 = arg_0_Z.getDataPointSize();
3143 
3144  // Declare output Data object
3145  Data res;
3146 
3147  if (arg_0_Z.isConstant()) {
3148  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataConstant output
3149  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
3150  double *ptr_2 = &(res.getDataAtOffsetRW(0));
3151  tensor_unary_operation(size0, ptr_0, ptr_2, operation);
3152  }
3153  else if (arg_0_Z.isTagged()) {
3154 
3155  // Borrow DataTagged input from Data object
3156  DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
3157 
3158  // Prepare a DataTagged output 2
3159  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
3160  res.tag();
3161  DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
3162 
3163  // Get the pointers to the actual data
3164  const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
3165  double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
3166  // Compute a result for the default
3167  tensor_unary_operation(size0, ptr_0, ptr_2, operation);
3168  // Compute a result for each tag
3169  const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
3170  DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
3171  for (i=lookup_0.begin();i!=lookup_0.end();i++) {
3172  tmp_2->addTag(i->first);
3173  const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
3174  double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
3175  tensor_unary_operation(size0, ptr_0, ptr_2, operation);
3176  }
3177 
3178  }
3179  else if (arg_0_Z.isExpanded()) {
3180 
3181  res = Data(0.0, shape0, arg_0_Z.getFunctionSpace(),true); // DataExpanded output
3182  DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
3183  DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
3184 
3185  int sampleNo_0,dataPointNo_0;
3186  int numSamples_0 = arg_0_Z.getNumSamples();
3187  int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
3188  #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
3189  for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
3190  dataPointNo_0=0;
3191 // for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
3192  int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
3193  int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
3194  const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
3195  double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
3196  tensor_unary_operation(size0*numDataPointsPerSample_0, ptr_0, ptr_2, operation);
3197 // }
3198  }
3199  }
3200  else {
3201  throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");
3202  }
3203 
3204  return res;
3205 }
3206 
3207 }
3208 #endif
double sup_const() const
Definition: Data.cpp:1616
Give a short description of what FunctionSpace does.
Definition: FunctionSpace.h:46
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:2746
DataAbstract_ptr m_data
Definition: Data.h:1758
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:1087
int getNumDataPointsPerSample() const
Return the number of data points per sample.
Definition: Data.h:536
DataTypes::ValueType::reference getDataAtOffsetRW(DataTypes::ValueType::size_type i)
Definition: Data.cpp:3164
Data asin() const
Return the asin of each data point of this Data object.
Definition: Data.cpp:1437
const DataTypes::ShapeType & getDataPointShape() const
Return a reference to the data point shape.
Definition: Data.h:681
double LsupWorker() const
Definition: Data.cpp:1740
Data wherePositive() const
Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.
Definition: Data.cpp:983
void addTag(int tagKey)
addTag - does not modify the default value for this object. ** Not unit tested ** ...
Definition: DataTagged.cpp:429
void setItemO(const boost::python::object &key, const boost::python::object &value)
Definition: Data.cpp:2601
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:3988
Data whereNonNegative() const
Return a Data with a 1 for +ive or 0 values and a 0 for -ive values.
Definition: Data.cpp:997
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:3976
unsigned int getDataPointRank() const
Return the rank of the point data.
Definition: Data.h:501
int getNumSamples() const
Return the number of samples.
Definition: Data.h:524
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:339
boost::python::object __add__(const boost::python::object &right)
wrapper for python add operation
Definition: Data.cpp:4411
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:1853
boost::python::object __rsub__(const boost::python::object &right)
wrapper for python reverse subtract operation
Definition: Data.cpp:4459
DataReady_ptr borrowReadyPtr(void) const
Definition: Data.cpp:3127
void setToZero()
set all values to zero
Definition: Data.cpp:595
Data getItem(const boost::python::object &key) const
Returns a slice from this Data object.
Definition: Data.cpp:2579
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:2269
bool isExpanded() const
Return true if this Data is expanded.
Definition: Data.cpp:847
double Lsup()
Return the maximum absolute value of this Data object.
Definition: Data.cpp:1595
Data tanh() const
Return the tanh of each data point of this Data object.
Definition: Data.cpp:1473
Data delay()
produce a delayed evaluation version of this Data.
Definition: Data.cpp:566
Data & operator*=(const Data &right)
Overloaded operator *=.
Definition: Data.cpp:2363
void tag()
If possible convert this Data to DataTagged. This will only allow Constant data to be converted to ta...
Definition: Data.cpp:934
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:431
boost::python::object integrateToTuple()
Calculate the integral over the function space domain as a python tuple.
Definition: Data.cpp:1368
int getDataPointSize() const
Return the size of the data point. It is the product of the data point shape dimensions.
Definition: Data.cpp:1069
boost::shared_ptr< AbstractDomain > Domain_ptr
Definition: AbstractDomain.h:57
double lazyAlgWorker(double init)
Definition: Data.cpp:1682
Definition: DataReady.h:35
double inf_const() const
Definition: Data.cpp:1647
bool m_shared
Definition: Data.h:1752
std::vector< std::pair< int, int > > RegionType
Definition: DataTypes.h:39
bool probeInterpolation(const FunctionSpace &functionspace) const
Returns true if this can be interpolated to functionspace.
Definition: Data.cpp:1035
escriptDataC getDataC()
Return the C wrapper for the Data object.
Definition: Data.cpp:509
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:3172
double infWorker() const
Definition: Data.cpp:1814
void setItemD(const boost::python::object &key, const Data &value)
Copies slice from value into this Data object.
Definition: Data.cpp:2609
void tensor_unary_operation(const int size, const double *arg1, double *argRes, UnaryFunction operation)
Definition: LocalOps.h:517
virtual ValueType::size_type getPointOffset(int sampleNo, int dataPointNo) const
getPointOffset
Definition: DataTagged.cpp:505
Data atan() const
Return the atan of each data point of this Data object.
Definition: Data.cpp:1452
Data operator-(const Data &left, const Data &right)
Operator- Takes two Data objects.
Definition: Data.cpp:2471
void transpose(const DataTypes::ValueType &in, const DataTypes::ShapeType &inShape, DataTypes::ValueType::size_type inOffset, DataTypes::ValueType &ev, const DataTypes::ShapeType &evShape, DataTypes::ValueType::size_type evOffset, int axis_offset)
Transpose each data point of this Data object around the given axis.
Definition: DataMaths.h:394
Data rpowO(const boost::python::object &left) const
Return the given power of each data point of this boost python object.
Definition: Data.cpp:2439
void print(void)
print the data values to stdout. Used for debugging
Definition: Data.cpp:3928
void typeMatchRight(const Data &right)
Convert the data type of this to match the RHS.
Definition: Data.cpp:2656
void setSlice(const Data &value, const DataTypes::RegionType &region)
Copy the specified slice from the given value into this Data object.
Definition: Data.cpp:2625
Data neg() const
Return the negation of each data point of this Data object.
Definition: Data.cpp:1553
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:3135
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:1195
Data log10() const
Return the log to base 10 of each data point of this Data object.
Definition: Data.cpp:1525
Data interpolateFromTable2DP(boost::python::object table, double Amin, double Astep, Data &B, double Bmin, double Bstep, double undef, bool check_boundaries)
Definition: Data.cpp:3225
Data getSlice(const DataTypes::RegionType &region) const
Return a Data object containing the specified slice of this Data object.
Definition: Data.cpp:2593
bool isLazy() const
Return true if this Data is lazy.
Definition: Data.cpp:882
Data swapaxes(const int axis0, const int axis1) const
swaps the components axis0 and axis1
Definition: Data.cpp:1888
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:2680
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:548
double Lsup_const() const
Definition: Data.cpp:1585
const DataMapType & getTagLookup() const
getTagLookup
Definition: DataTagged.h:660
DataReady_ptr getReadyPtr()
Definition: Data.h:1937
double * getSampleDataRW(ValueType::size_type sampleNo)
Return the sample data for the given sample number.
Definition: DataReady.h:113
DataTypes::ValueType::const_reference getDataByTagRO(int tag, DataTypes::ValueType::size_type i) const
Definition: DataTagged.cpp:542
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:117
Data trace(int axis_offset) const
Return the trace of a matrix.
Definition: Data.cpp:2003
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:1862
Data symmetric() const
Return the symmetric part of a matrix which is half the matrix plus its transpose.
Definition: Data.cpp:1945
(Testing use only) Provides a domain to wrap a collection of values.
Definition: TestDomain.h:45
Data nonuniformslope(boost::python::object in, boost::python::object out, bool check_boundaries)
Definition: Data.cpp:3849
Domain_ptr getDomainPython() const
Return the function space domain. Internal use only! This gets around some python difficulties by cas...
Definition: FunctionSpace.cpp:109
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:4435
DataAbstract_ptr borrowDataPtr(void) const
Definition: Data.cpp:3120
Simulates a full dataset accessible via sampleNo and dataPointNo.
Definition: DataTagged.h:43
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't.
Definition: Data.h:1668
Data interpolateFromTable1D(const WrappedArray &table, double Amin, double Astep, double undef, bool check_boundaries)
Definition: Data.cpp:3243
Data sqrt() const
Return the square root of each data point of this Data object.
Definition: Data.cpp:1578
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't tagged an e...
Definition: Data.h:634
bool isConstant() const
Return true if this Data is constant.
Definition: Data.cpp:875
DataTypes::ValueType & getExpandedVectorReference()
Ensures that the Data is expanded and returns its underlying vector Does not check for exclusive writ...
Definition: Data.cpp:4345
void binaryOp(const Data &right, BinaryFunction operation)
Perform the given binary operation on all of the data's elements. The underlying type of the right ha...
Definition: Data.h:2148
Data cos() const
Return the cos of each data point of this Data object.
Definition: Data.cpp:1423
Data minval() const
Return the minimum value of each data point of this Data object.
Definition: Data.cpp:1880
Data matrixInverse() const
return inverse of matricies.
Definition: Data.cpp:2413
int MPI_Comm
Definition: Esys_MPI.h:29
Data cosh() const
Return the cosh of each data point of this Data object.
Definition: Data.cpp:1466
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:2753
bool m_protected
Definition: Data.h:1751
bool isTagged() const
Return true if this Data is tagged.
Definition: Data.cpp:861
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:1041
Data transpose(int axis_offset) const
Transpose each data point of this Data object around the given axis.
Definition: Data.cpp:2058
Data & operator/=(const Data &right)
Overloaded operator /=.
Definition: Data.cpp:2386
Data pos() const
Return the identity of each data point of this Data object. Simply returns this object unmodified...
Definition: Data.cpp:1560
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:960
Data powD(const Data &right) const
Return the given power of each data point of this Data object.
Definition: Data.cpp:2453
bool isReady() const
Return true if this data is ready.
Definition: Data.cpp:889
DataTypes::ValueType::reference getDataByTagRW(int tag, DataTypes::ValueType::size_type i)
getDataByTag
Definition: DataTagged.cpp:553
Data sin() const
Return the sin of each data point of this Data object.
Definition: Data.cpp:1416
void calc_minGlobalDataPoint(int &ProcNo, int &DataPointNo) const
Definition: Data.cpp:2143
Data minval_nonlazy() const
Definition: Data.cpp:1852
const boost::python::tuple getShapeTuple() const
Return the data point shape as a tuple of integers.
Definition: Data.cpp:525
Data C_TensorUnaryOperation(Data const &arg_0, UnaryFunction operation)
Definition: Data.h:3126
bool isShared() const
Definition: Data.h:1803
Data acos() const
Return the acos of each data point of this Data object.
Definition: Data.cpp:1444
bool isProtected() const
Returns true, if the data object is protected against update.
Definition: Data.cpp:902
double(* UnaryDFunPtr)(double)
Definition: Data.h:77
Data rtruedivO(const boost::python::object &left)
Newer style division operator for python.
Definition: Data.h:2134
Data interpolateFromTable1DP(boost::python::object table, double Amin, double Astep, double undef, bool check_boundaries)
Definition: Data.cpp:3234
std::map< int, int > DataMapType
Definition: DataTagged.h:56
std::string toString() const
Return a text description of the function space.
Definition: FunctionSpace.cpp:118
Data log() const
Return the natural log of each data point of this Data object.
Definition: Data.cpp:1532
Data nonsymmetric() const
Return the nonsymmetric part of a matrix which is half the matrix minus its transpose.
Definition: Data.cpp:1968
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:3215
Provide a wrapper around a Data object so Data may be accessed from C.
Definition: DataC.h:30
void delaySelf()
convert the current data into lazy data.
Definition: Data.cpp:577
bool hasNoSamples() const
Return true if this object contains no samples. This is not the same as isEmpty() ...
Definition: Data.h:717
double rpow(double x, double y)
Definition: Data.h:1982
DataTypes::ValueType::size_type getDataOffset(int sampleNo, int dataPointNo)
Return the offset for the given sample and point within the sample.
Definition: Data.h:668
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:2717
void expand()
Whatever the current Data type make this into a DataExpanded.
Definition: Data.cpp:910
const FunctionSpace & getFunctionSpace() const
Return the function space.
Definition: Data.h:461
Data operator+(const Data &left, const Data &right)
Operator+ Takes two Data objects.
Definition: Data.cpp:2462
boost::python::object integrateToTuple_const() const
Calculate the integral over the function space domain as a python tuple.
Definition: Data.cpp:1358
DataTypes::ValueType::size_type getLength() const
Return the number of doubles stored for this Data.
Definition: Data.cpp:1076
Data represents a collection of datapoints.
Definition: Data.h:71
void updateShareStatus(bool nowshared) const
Update the Data's shared flag This indicates that the DataAbstract used by this object is now shared ...
Definition: Data.h:1784
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:2235
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:1967
Data sinh() const
Return the sinh of each data point of this Data object.
Definition: Data.cpp:1459
double supWorker() const
Definition: Data.cpp:1778
Data maxval() const
Return the maximum value of each data point of this Data object.
Definition: Data.cpp:1872
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:3530
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:1247
#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:3964
const ElementType & const_reference
Definition: DataVector.h:62
double(* BinaryDFunPtr)(double, double)
Definition: Data.h:78
const double * getSampleDataRO(ValueType::size_type sampleNo) const
Definition: DataReady.h:119
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:1956
boost::python::object integrateWorker() const
Definition: Data.cpp:1379
Data & operator=(const Data &other)
Definition: Data.cpp:2332
void forceResolve()
Definition: Data.h:1815
Data C_TensorBinaryOperation(Data const &arg_0, Data const &arg_1, BinaryFunction operation)
Compute a tensor operation with two Data objects.
Definition: Data.h:2317
int getTypeCode() const
Return the function space type code.
Definition: FunctionSpace.cpp:95
int getNumDataPoints() const
Return the number of data points.
Definition: Data.h:513
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:969
Data whereNonPositive() const
Return a Data with a 1 for -ive or 0 values and a 0 for +ive values.
Definition: Data.cpp:1004
Data atanh() const
Return the atanh of each data point of this Data object.
Definition: Data.cpp:1514
Data copySelf()
Return a pointer to a deep copy of this object.
Definition: Data.cpp:550
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:3191
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:646
bool m_lazy
Definition: Data.h:1753
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:868
void setTupleForGlobalDataPoint(int id, int proc, boost::python::object)
Set the value of a global data point.
Definition: Data.cpp:1203
DataTypes::ValueType::const_reference getDefaultValueRO(DataTypes::ValueType::size_type i) const
Definition: DataTagged.h:653
boost::python::object __div__(const boost::python::object &right)
wrapper for python divide operation
Definition: Data.cpp:4508
DataAbstract * borrowData(void) const
return the object produced by the factory, which is a DataConstant or DataExpanded TODO Ownership of ...
Definition: Data.cpp:3113
Data grad() const
Definition: Data.cpp:1059
void setValueOfDataPoint(int dataPointNo, const double)
sets the values of a data-point on this process
Definition: Data.cpp:1284
friend Data condEval(escript::Data &mask, escript::Data &trueval, escript::Data &falseval)
void calc_maxGlobalDataPoint(int &ProcNo, int &DataPointNo) const
Definition: Data.cpp:2229
void copy(const Data &other)
Make this object a deep copy of "other".
Definition: Data.cpp:557
Wraps an expression tree of other DataObjects. The data will be evaluated when required.
Definition: DataLazy.h:102
int getTagNumber(int dpno)
Return the tag number associated with the given data-point.
Definition: Data.cpp:2736
double inf()
Return the minimum value of this Data object.
Definition: Data.cpp:1657
DataException exception class.
Definition: DataException.h:35
#define ESCRIPT_DLL_API
Definition: escriptcore/src/system_dep.h:64
Data nonuniforminterp(boost::python::object in, boost::python::object out, bool check_boundaries)
Definition: Data.cpp:3773
Data erf() const
Return the error function erf of each data point of this Data object.
Definition: Data.cpp:1481
Data acosh() const
Return the acosh of each data point of this Data object.
Definition: Data.cpp:1503
Data condEval(escript::Data &mask, escript::Data &trueval, escript::Data &falseval)
Definition: Data.cpp:4173
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:56
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 sign() const
Return the sign of each data point of this Data object. -1 for negative values, zero for zero values...
Definition: Data.cpp:1539
void exclusiveWrite()
if another object is sharing out member data make a copy to work with instead. This code should only ...
Definition: Data.h:1833
Data oneOver() const
Returns 1./ Data object.
Definition: Data.cpp:976
Data abs() const
Return the absolute value of each data point of this Data object.
Definition: Data.cpp:1546
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:2130
const ValueType * resolveSample(int sampleNo, size_t &roffset) const
Compute the value of the expression for the given sample.
Definition: DataLazy.cpp:1642
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:854
void setProtection()
switches on update protection
Definition: Data.cpp:896
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:4532
Data powO(const boost::python::object &right) const
Return the given power of each data point of this boost python object.
Definition: Data.cpp:2446
Data interpolateFromTable2D(const WrappedArray &table, double Amin, double Astep, double undef, Data &B, double Bmin, double Bstep, bool check_boundaries)
Definition: Data.cpp:3371
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:2220
Data eigenvalues() const
Return the eigenvalues of the symmetric part at each data point of this Data object in increasing val...
Definition: Data.cpp:2082
Data tan() const
Return the tan of each data point of this Data object.
Definition: Data.cpp:1430
Data operator/(const Data &left, const Data &right)
Operator/ Takes two Data objects.
Definition: Data.cpp:2489
Data()
Default constructor. Creates a DataEmpty object.
Definition: Data.cpp:239
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:466
void copyWithMask(const Data &other, const Data &mask)
Copy other Data object into this Data object where mask is positive.
Definition: Data.cpp:617
boost::python::object __mul__(const boost::python::object &right)
wrapper for python multiply operation
Definition: Data.cpp:4484
void typeMatchLeft(Data &right) const
Convert the data type of the RHS to match this.
Definition: Data.cpp:2640
long size_type
Definition: DataVector.h:60
int getNoValues() const
Return the number of values in the shape for this object.
Definition: Data.h:574
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:1011
Data asinh() const
Return the asinh of each data point of this Data object.
Definition: Data.cpp:1492
DataTypes::ValueType::reference getDataPointRW(int sampleNo, int dataPointNo)
Return a reference into the DataVector which points to the specified data point.
Definition: Data.cpp:3207
Data interpolate(const FunctionSpace &functionspace) const
Interpolates this onto the given functionspace and returns the result as a Data object.
Definition: Data.cpp:1029
Data truedivO(const boost::python::object &right)
Newer style division operator for python.
Definition: Data.h:2121
Definition: WrappedArray.h:29
double sup()
Return the maximum value of this Data object.
Definition: Data.cpp:1626
Data exp() const
Return the exponential function of each data point of this Data object.
Definition: Data.cpp:1571
~Data()
Destructor.
Definition: Data.cpp:423
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:529
const boost::python::object getValueOfDataPointAsTuple(int dataPointNo)
Return the value of a data point as a python tuple.
Definition: Data.cpp:1164
const_Domain_ptr getDomain() const
Return the function space domain.
Definition: FunctionSpace.cpp:103
Data truedivD(const Data &right)
Newer style division operator for python.
Definition: Data.h:2109
boost::shared_ptr< const AbstractDomain > const_Domain_ptr
Definition: AbstractDomain.h:60
Domain_ptr getDomainPython() const
Return the domain. TODO: For internal use only. This should be removed.
Definition: Data.h:489
Data & operator-=(const Data &right)
Overloaded operator -=.
Definition: Data.cpp:2340
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:1020
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:2696
const boost::python::object getValueOfGlobalDataPointAsTuple(int procNo, int dataPointNo)
Return a data point across all processors as a python tuple.
Definition: Data.cpp:1306
void dump(const std::string fileName) const
dumps the object into a netCDF file
Definition: Data.cpp:3942
const DataReady * getReady() const
Definition: Data.h:1916
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:2105
Data whereNegative() const
Return a Data with a 1 for -ive values and a 0 for +ive or 0 values.
Definition: Data.cpp:990
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:560
const_Domain_ptr getDomain() const
Return the domain.
Definition: Data.h:474
Data & operator+=(const Data &right)
Overloaded operator +=.
Definition: Data.cpp:2308
void binaryOp(DataTagged &left, const DataConstant &right, BinaryFunction operation)
Perform the given binary operation.
Definition: BinaryOp.h:45