/[escript]/branches/clazy/escriptcore/src/DataLazy.cpp
ViewVC logotype

Diff of /branches/clazy/escriptcore/src/DataLazy.cpp

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

trunk/escriptcore/src/DataLazy.cpp revision 5985 by jfenwick, Fri Feb 26 02:21:44 2016 UTC branches/clazy/escriptcore/src/DataLazy.cpp revision 6511 by jfenwick, Fri Mar 3 01:41:39 2017 UTC
# Line 5  Line 5 
5  * http://www.uq.edu.au  * http://www.uq.edu.au
6  *  *
7  * Primary Business: Queensland, Australia  * Primary Business: Queensland, Australia
8  * Licensed under the Open Software License version 3.0  * Licensed under the Apache License, version 2.0
9  * http://www.opensource.org/licenses/osl-3.0.php  * http://www.apache.org/licenses/LICENSE-2.0
10  *  *
11  * Development until 2012 by Earth Systems Science Computational Center (ESSCC)  * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12  * Development 2012-2013 by School of Earth Sciences  * Development 2012-2013 by School of Earth Sciences
# Line 14  Line 14 
14  *  *
15  *****************************************************************************/  *****************************************************************************/
16    
 #define ESNEEDPYTHON  
 #include "esysUtils/first.h"  
   
17  #include "DataLazy.h"  #include "DataLazy.h"
18  #include "Data.h"  #include "Data.h"
19  #include "DataTypes.h"  #include "DataTypes.h"
20  #include "EscriptParams.h"  #include "EscriptParams.h"
21  #include "FunctionSpace.h"  #include "FunctionSpace.h"
 #include "UnaryFuncs.h"    // for escript::fsign  
22  #include "Utils.h"  #include "Utils.h"
23    #include "DataVectorOps.h"
24    
25  #include "esysUtils/Esys_MPI.h"  #include <iomanip> // for some fancy formatting in debug
   
 #ifdef _OPENMP  
 #include <omp.h>  
 #endif  
 #ifdef USE_NETCDF  
 #include <netcdfcpp.h>  
 #endif  
   
 #include <iomanip>              // for some fancy formatting in debug  
26    
27  using namespace escript::DataTypes;  using namespace escript::DataTypes;
28    
# Line 50  bool privdebug=false; Line 38  bool privdebug=false;
38  #define DISABLEDEBUG privdebug=false;  #define DISABLEDEBUG privdebug=false;
39  }  }
40    
41  // #define SIZELIMIT if ((m_height>escript::escriptParams.getTOO_MANY_LEVELS()) || (m_children>escript::escriptParams.getTOO_MANY_NODES())) {cerr << "\n!!!!!!! SIZE LIMIT EXCEEDED " << m_children << ";" << m_height << endl << toString() << endl;resolveToIdentity();}  //#define SIZELIMIT if ((m_height>escript::escriptParams.getInt("TOO_MANY_LEVELS")) || (m_children>escript::escriptParams.getInt("TOO_MANY_NODES"))) {cerr << "\n!!!!!!! SIZE LIMIT EXCEEDED " << m_children << ";" << m_height << endl << toString() << endl;resolveToIdentity();}
42    
43  // #define SIZELIMIT if ((m_height>escript::escriptParams.getTOO_MANY_LEVELS()) || (m_children>escript::escriptParams.getTOO_MANY_NODES())) {cerr << "SIZE LIMIT EXCEEDED " << m_height << endl;resolveToIdentity();}  //#define SIZELIMIT if ((m_height>escript::escriptParams.getInt("TOO_MANY_LEVELS")) || (m_children>escript::escriptParams.getInt("TOO_MANY_NODES"))) {cerr << "SIZE LIMIT EXCEEDED " << m_height << endl;resolveToIdentity();}
44    
45    
46  #define SIZELIMIT if (m_height>escript::escriptParams.getTOO_MANY_LEVELS())  {if (escript::escriptParams.getLAZY_VERBOSE()){cerr << "SIZE LIMIT EXCEEDED height=" << m_height << endl;}resolveToIdentity();}  #define SIZELIMIT \
47        if (m_height > escript::escriptParams.getTooManyLevels()) {\
48            if (escript::escriptParams.getLazyVerbose()) {\
49                cerr << "SIZE LIMIT EXCEEDED height=" << m_height << endl;\
50            }\
51            resolveToIdentity();\
52        }
53    
54  /*  /*
55  How does DataLazy work?  How does DataLazy work?
# Line 139  std::vector<void*> stackend(getNumberOfT Line 133  std::vector<void*> stackend(getNumberOfT
133  size_t maxstackuse=0;  size_t maxstackuse=0;
134  #endif  #endif
135    
 enum ES_opgroup  
 {  
    G_UNKNOWN,  
    G_IDENTITY,  
    G_BINARY,            // pointwise operations with two arguments  
    G_UNARY,             // pointwise operations with one argument  
    G_UNARY_P,           // pointwise operations with one argument, requiring a parameter  
    G_NP1OUT,            // non-pointwise op with one output  
    G_NP1OUT_P,          // non-pointwise op with one output requiring a parameter  
    G_TENSORPROD,        // general tensor product  
    G_NP1OUT_2P,         // non-pointwise op with one output requiring two params  
    G_REDUCTION,         // non-pointwise unary op with a scalar output  
    G_CONDEVAL  
 };  
   
   
   
   
 string ES_opstrings[]={"UNKNOWN","IDENTITY","+","-","*","/","^",  
                         "sin","cos","tan",  
                         "asin","acos","atan","sinh","cosh","tanh","erf",  
                         "asinh","acosh","atanh",  
                         "log10","log","sign","abs","neg","pos","exp","sqrt",  
                         "1/","where>0","where<0","where>=0","where<=0", "where<>0","where=0",  
                         "symmetric","nonsymmetric",  
                         "prod",  
                         "transpose", "trace",  
                         "swapaxes",  
                         "minval", "maxval",  
                         "condEval"};  
 int ES_opcount=44;  
 ES_opgroup opgroups[]={G_UNKNOWN,G_IDENTITY,G_BINARY,G_BINARY,G_BINARY,G_BINARY, G_BINARY,  
                         G_UNARY,G_UNARY,G_UNARY, //10  
                         G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,        // 17  
                         G_UNARY,G_UNARY,G_UNARY,                                        // 20  
                         G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,        // 28  
                         G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY, G_UNARY_P, G_UNARY_P,          // 35  
                         G_NP1OUT,G_NP1OUT,  
                         G_TENSORPROD,  
                         G_NP1OUT_P, G_NP1OUT_P,  
                         G_NP1OUT_2P,  
                         G_REDUCTION, G_REDUCTION,  
                         G_CONDEVAL};  
 inline  
 ES_opgroup  
 getOpgroup(ES_optype op)  
 {  
   return opgroups[op];  
 }  
136    
137  // return the FunctionSpace of the result of "left op right"  // return the FunctionSpace of the result of "left op right"
138  FunctionSpace  FunctionSpace
# Line 440  GTPShape(DataAbstract_ptr left, DataAbst Line 385  GTPShape(DataAbstract_ptr left, DataAbst
385    
386  }       // end anonymous namespace  }       // end anonymous namespace
387    
   
   
 // Return a string representing the operation  
 const std::string&  
 opToString(ES_optype op)  
 {  
   if (op<0 || op>=ES_opcount)  
   {  
     op=UNKNOWNOP;  
   }  
   return ES_opstrings[op];  
 }  
   
388  void DataLazy::LazyNodeSetup()  void DataLazy::LazyNodeSetup()
389  {  {
390  #ifdef _OPENMP  #ifdef _OPENMP
# Line 486  DataLazy::DataLazy(DataAbstract_ptr p) Line 418  DataLazy::DataLazy(DataAbstract_ptr p)
418     }     }
419     else     else
420     {     {
         p->makeLazyShared();  
421          DataReady_ptr dr=dynamic_pointer_cast<DataReady>(p);          DataReady_ptr dr=dynamic_pointer_cast<DataReady>(p);
422          makeIdentity(dr);          makeIdentity(dr);
423  LAZYDEBUG(cout << "Wrapping " << dr.get() << " id=" << m_id.get() << endl;)  LAZYDEBUG(cout << "Wrapping " << dr.get() << " id=" << m_id.get() << endl;)
# Line 951  DataLazy::collapseToReady() const Line 882  DataLazy::collapseToReady() const
882          result=left.symmetric();          result=left.symmetric();
883          break;          break;
884      case NSYM:      case NSYM:
885          result=left.nonsymmetric();          result=left.antisymmetric();
886          break;          break;
887      case PROD:      case PROD:
888          result=C_GeneralTensorProduct(left,right,m_axis_offset, m_transpose);          result=C_GeneralTensorProduct(left,right,m_axis_offset, m_transpose);
# Line 971  DataLazy::collapseToReady() const Line 902  DataLazy::collapseToReady() const
902      case MAXVAL:      case MAXVAL:
903          result=left.minval();          result=left.minval();
904          break;          break;
905        case HER:
906        result=left.hermitian();
907        break;
908      default:      default:
909          throw DataException("Programmer error - collapseToReady does not know how to resolve operator "+opToString(m_op)+".");          throw DataException("Programmer error - collapseToReady does not know how to resolve operator "+opToString(m_op)+".");
910    }    }
# Line 998  DataLazy::collapse() const Line 932  DataLazy::collapse() const
932    m_op=IDENTITY;    m_op=IDENTITY;
933  }  }
934    
   
   
   
   
   
 #define PROC_OP(TYPE,X)                               \  
         for (int j=0;j<onumsteps;++j)\  
         {\  
           for (int i=0;i<numsteps;++i,resultp+=resultStep) \  
           { \  
 LAZYDEBUG(cout << "[left,right]=[" << lroffset << "," << rroffset << "]" << endl;)\  
 LAZYDEBUG(cout << "{left,right}={" << (*left)[lroffset] << "," << (*right)[rroffset] << "}\n";)\  
              tensor_binary_operation< TYPE >(chunksize, &((*left)[lroffset]), &((*right)[rroffset]), resultp, X); \  
 LAZYDEBUG(cout << " result=      " << resultp[0] << endl;) \  
              lroffset+=leftstep; \  
              rroffset+=rightstep; \  
           }\  
           lroffset+=oleftstep;\  
           rroffset+=orightstep;\  
         }  
   
   
935  // The result will be stored in m_samples  // The result will be stored in m_samples
936  // The return value is a pointer to the DataVector, offset is the offset within the return value  // The return value is a pointer to the DataVector, offset is the offset within the return value
937  const DataTypes::RealVectorType*  const DataTypes::RealVectorType*
# Line 1033  LAZYDEBUG(cout << "Resolve sample " << t Line 945  LAZYDEBUG(cout << "Resolve sample " << t
945    }    }
946    if (m_op==IDENTITY)      if (m_op==IDENTITY)  
947    {    {
948      const ValueType& vec=m_id->getVectorRO();      const RealVectorType& vec=m_id->getVectorRO();
949      roffset=m_id->getPointOffset(sampleNo, 0);      roffset=m_id->getPointOffset(sampleNo, 0);
950  #ifdef LAZY_STACK_PROF  #ifdef LAZY_STACK_PROF
951  int x;  int x;
# Line 1090  DataLazy::resolveNodeUnary(int tid, int Line 1002  DataLazy::resolveNodeUnary(int tid, int
1002    const double* left=&((*leftres)[roffset]);    const double* left=&((*leftres)[roffset]);
1003    roffset=m_samplesize*tid;    roffset=m_samplesize*tid;
1004    double* result=&(m_samples[roffset]);    double* result=&(m_samples[roffset]);
1005    switch (m_op)    if (m_op==POS)
1006    {    {
1007      case SIN:        // this should be prevented earlier
1008          tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::sin);      // operation is meaningless for lazy
1009          break;          throw DataException("Programmer error - POS not supported for lazy data.");    
1010      case COS:    }
1011          tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::cos);    tensor_unary_array_operation(m_samplesize,
1012          break;                               left,
1013      case TAN:                               result,
1014          tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::tan);                               m_op,
1015          break;                               m_tol);  
     case ASIN:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::asin);  
         break;  
     case ACOS:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::acos);  
         break;  
     case ATAN:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::atan);  
         break;  
     case SINH:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::sinh);  
         break;  
     case COSH:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::cosh);  
         break;  
     case TANH:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::tanh);  
         break;  
     case ERF:  
 #if defined (_WIN32) && !defined(__INTEL_COMPILER)  
         throw DataException("Error - Data:: erf function is not supported on _WIN32 platforms.");  
 #else  
         tensor_unary_operation(m_samplesize, left, result, ::erf);  
         break;  
 #endif  
    case ASINH:  
 #if defined (_WIN32) && !defined(__INTEL_COMPILER)  
         tensor_unary_operation(m_samplesize, left, result, escript::asinh_substitute);  
 #else  
         tensor_unary_operation(m_samplesize, left, result, ::asinh);  
 #endif    
         break;  
    case ACOSH:  
 #if defined (_WIN32) && !defined(__INTEL_COMPILER)  
         tensor_unary_operation(m_samplesize, left, result, escript::acosh_substitute);  
 #else  
         tensor_unary_operation(m_samplesize, left, result, ::acosh);  
 #endif    
         break;  
    case ATANH:  
 #if defined (_WIN32) && !defined(__INTEL_COMPILER)  
         tensor_unary_operation(m_samplesize, left, result, escript::atanh_substitute);  
 #else  
         tensor_unary_operation(m_samplesize, left, result, ::atanh);  
 #endif    
         break;  
     case LOG10:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::log10);  
         break;  
     case LOG:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::log);  
         break;  
     case SIGN:  
         tensor_unary_operation(m_samplesize, left, result, escript::fsign);  
         break;  
     case ABS:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::fabs);  
         break;  
     case NEG:  
         tensor_unary_operation(m_samplesize, left, result, negate<double>());  
         break;  
     case POS:  
         // it doesn't mean anything for delayed.  
         // it will just trigger a deep copy of the lazy object  
         throw DataException("Programmer error - POS not supported for lazy data.");  
         break;  
     case EXP:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::exp);  
         break;  
     case SQRT:  
         tensor_unary_operation<double (*)(double)>(m_samplesize, left, result, ::sqrt);  
         break;  
     case RECIP:  
         tensor_unary_operation(m_samplesize, left, result, bind1st(divides<double>(),1.));  
         break;  
     case GZ:  
         tensor_unary_operation(m_samplesize, left, result, bind2nd(greater<double>(),0.0));  
         break;  
     case LZ:  
         tensor_unary_operation(m_samplesize, left, result, bind2nd(less<double>(),0.0));  
         break;  
     case GEZ:  
         tensor_unary_operation(m_samplesize, left, result, bind2nd(greater_equal<double>(),0.0));  
         break;  
     case LEZ:  
         tensor_unary_operation(m_samplesize, left, result, bind2nd(less_equal<double>(),0.0));  
         break;  
 // There are actually G_UNARY_P but I don't see a compelling reason to treat them differently  
     case NEZ:  
         tensor_unary_operation(m_samplesize, left, result, bind2nd(AbsGT(),m_tol));  
         break;  
     case EZ:  
         tensor_unary_operation(m_samplesize, left, result, bind2nd(AbsLTE(),m_tol));  
         break;  
   
     default:  
         throw DataException("Programmer error - resolveUnary can not resolve operator "+opToString(m_op)+".");  
   }  
1016    return &(m_samples);    return &(m_samples);
1017  }  }
1018    
# Line 1232  DataLazy::resolveNodeReduction(int tid, Line 1046  DataLazy::resolveNodeReduction(int tid,
1046            for (unsigned int z=0;z<ndpps;++z)            for (unsigned int z=0;z<ndpps;++z)
1047            {            {
1048              FMin op;              FMin op;
1049              *result=DataMaths::reductionOp(*leftres, m_left->getShape(), loffset, op, numeric_limits<double>::max());              *result=escript::reductionOpVector(*leftres, m_left->getShape(), loffset, op, numeric_limits<double>::max());
1050              loffset+=psize;              loffset+=psize;
1051              result++;              result++;
1052            }            }
# Line 1243  DataLazy::resolveNodeReduction(int tid, Line 1057  DataLazy::resolveNodeReduction(int tid,
1057            for (unsigned int z=0;z<ndpps;++z)            for (unsigned int z=0;z<ndpps;++z)
1058            {            {
1059            FMax op;            FMax op;
1060            *result=DataMaths::reductionOp(*leftres, m_left->getShape(), loffset, op, numeric_limits<double>::max()*-1);            *result=escript::reductionOpVector(*leftres, m_left->getShape(), loffset, op, numeric_limits<double>::max()*-1);
1061            loffset+=psize;            loffset+=psize;
1062            result++;            result++;
1063            }            }
# Line 1270  DataLazy::resolveNodeNP1OUT(int tid, int Line 1084  DataLazy::resolveNodeNP1OUT(int tid, int
1084      throw DataException("Programmer error - resolveNodeNP1OUT should not be called on identity nodes.");      throw DataException("Programmer error - resolveNodeNP1OUT should not be called on identity nodes.");
1085    }    }
1086    size_t subroffset;    size_t subroffset;
1087    const ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);    const RealVectorType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1088    roffset=m_samplesize*tid;    roffset=m_samplesize*tid;
1089    size_t loop=0;    size_t loop=0;
1090    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;
# Line 1281  DataLazy::resolveNodeNP1OUT(int tid, int Line 1095  DataLazy::resolveNodeNP1OUT(int tid, int
1095      case SYM:      case SYM:
1096          for (loop=0;loop<numsteps;++loop)          for (loop=0;loop<numsteps;++loop)
1097          {          {
1098              DataMaths::symmetric(*leftres,m_left->getShape(),subroffset, m_samples, getShape(), offset);              escript::symmetric(*leftres,m_left->getShape(),subroffset, m_samples, getShape(), offset);
1099              subroffset+=step;              subroffset+=step;
1100              offset+=step;              offset+=step;
1101          }          }
# Line 1289  DataLazy::resolveNodeNP1OUT(int tid, int Line 1103  DataLazy::resolveNodeNP1OUT(int tid, int
1103      case NSYM:      case NSYM:
1104          for (loop=0;loop<numsteps;++loop)          for (loop=0;loop<numsteps;++loop)
1105          {          {
1106              DataMaths::nonsymmetric(*leftres,m_left->getShape(),subroffset, m_samples, getShape(), offset);              escript::antisymmetric(*leftres,m_left->getShape(),subroffset, m_samples, getShape(), offset);
1107              subroffset+=step;              subroffset+=step;
1108              offset+=step;              offset+=step;
1109          }          }
# Line 1316  DataLazy::resolveNodeNP1OUT_P(int tid, i Line 1130  DataLazy::resolveNodeNP1OUT_P(int tid, i
1130    }    }
1131    size_t subroffset;    size_t subroffset;
1132    size_t offset;    size_t offset;
1133    const ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);    const RealVectorType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1134    roffset=m_samplesize*tid;    roffset=m_samplesize*tid;
1135    offset=roffset;    offset=roffset;
1136    size_t loop=0;    size_t loop=0;
# Line 1328  DataLazy::resolveNodeNP1OUT_P(int tid, i Line 1142  DataLazy::resolveNodeNP1OUT_P(int tid, i
1142      case TRACE:      case TRACE:
1143          for (loop=0;loop<numsteps;++loop)          for (loop=0;loop<numsteps;++loop)
1144          {          {
1145              DataMaths::trace(*leftres,m_left->getShape(),subroffset, m_samples ,getShape(),offset,m_axis_offset);              escript::trace(*leftres,m_left->getShape(),subroffset, m_samples ,getShape(),offset,m_axis_offset);
1146              subroffset+=instep;              subroffset+=instep;
1147              offset+=outstep;              offset+=outstep;
1148          }          }
# Line 1336  DataLazy::resolveNodeNP1OUT_P(int tid, i Line 1150  DataLazy::resolveNodeNP1OUT_P(int tid, i
1150      case TRANS:      case TRANS:
1151          for (loop=0;loop<numsteps;++loop)          for (loop=0;loop<numsteps;++loop)
1152          {          {
1153              DataMaths::transpose(*leftres,m_left->getShape(),subroffset, m_samples, getShape(),offset,m_axis_offset);              escript::transpose(*leftres,m_left->getShape(),subroffset, m_samples, getShape(),offset,m_axis_offset);
1154              subroffset+=instep;              subroffset+=instep;
1155              offset+=outstep;              offset+=outstep;
1156          }          }
# Line 1361  DataLazy::resolveNodeNP1OUT_2P(int tid, Line 1175  DataLazy::resolveNodeNP1OUT_2P(int tid,
1175    }    }
1176    size_t subroffset;    size_t subroffset;
1177    size_t offset;    size_t offset;
1178    const ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);    const RealVectorType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1179    roffset=m_samplesize*tid;    roffset=m_samplesize*tid;
1180    offset=roffset;    offset=roffset;
1181    size_t loop=0;    size_t loop=0;
# Line 1373  DataLazy::resolveNodeNP1OUT_2P(int tid, Line 1187  DataLazy::resolveNodeNP1OUT_2P(int tid,
1187      case SWAP:      case SWAP:
1188          for (loop=0;loop<numsteps;++loop)          for (loop=0;loop<numsteps;++loop)
1189          {          {
1190              DataMaths::swapaxes(*leftres,m_left->getShape(),subroffset, m_samples, getShape(),offset, m_axis_offset, m_transpose);              escript::swapaxes(*leftres,m_left->getShape(),subroffset, m_samples, getShape(),offset, m_axis_offset, m_transpose);
1191              subroffset+=instep;              subroffset+=instep;
1192              offset+=outstep;              offset+=outstep;
1193          }          }
# Line 1397  DataLazy::resolveNodeCondEval(int tid, i Line 1211  DataLazy::resolveNodeCondEval(int tid, i
1211    }    }
1212    size_t subroffset;    size_t subroffset;
1213    
1214    const ValueType* maskres=m_mask->resolveNodeSample(tid, sampleNo, subroffset);    const RealVectorType* maskres=m_mask->resolveNodeSample(tid, sampleNo, subroffset);
1215    const ValueType* srcres=0;    const RealVectorType* srcres=0;
1216    if ((*maskres)[subroffset]>0)    if ((*maskres)[subroffset]>0)
1217    {    {
1218          srcres=m_left->resolveNodeSample(tid, sampleNo, subroffset);          srcres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
# Line 1541  LAZYDEBUG(cout << "Resolve binary: " << Line 1355  LAZYDEBUG(cout << "Resolve binary: " <<
1355    
1356    int resultStep=max(leftstep,rightstep);       // only one (at most) should be !=0    int resultStep=max(leftstep,rightstep);       // only one (at most) should be !=0
1357          // Get the values of sub-expressions          // Get the values of sub-expressions
1358    const ValueType* left=m_left->resolveNodeSample(tid,sampleNo,lroffset);          const RealVectorType* left=m_left->resolveNodeSample(tid,sampleNo,lroffset);      
1359    const ValueType* right=m_right->resolveNodeSample(tid,sampleNo,rroffset);    const RealVectorType* right=m_right->resolveNodeSample(tid,sampleNo,rroffset);
1360  LAZYDEBUG(cout << "Post sub calls in " << toString() << endl;)  LAZYDEBUG(cout << "Post sub calls in " << toString() << endl;)
1361  LAZYDEBUG(cout << "shapes=" << DataTypes::shapeToString(m_left->getShape()) << "," << DataTypes::shapeToString(m_right->getShape()) << endl;)  LAZYDEBUG(cout << "shapes=" << DataTypes::shapeToString(m_left->getShape()) << "," << DataTypes::shapeToString(m_right->getShape()) << endl;)
1362  LAZYDEBUG(cout << "chunksize=" << chunksize << endl << "leftstep=" << leftstep << " rightstep=" << rightstep;)  LAZYDEBUG(cout << "chunksize=" << chunksize << endl << "leftstep=" << leftstep << " rightstep=" << rightstep;)
# Line 1560  LAZYDEBUG(cout << "Right res["<< rroffse Line 1374  LAZYDEBUG(cout << "Right res["<< rroffse
1374    switch(m_op)    switch(m_op)
1375    {    {
1376      case ADD:      case ADD:
1377          PROC_OP(NO_ARG,plus<double>());          //PROC_OP(NO_ARG,plus<double>());
1378          escript::binaryOpVectorLazyHelper<real_t, real_t, real_t>(resultp,
1379                 &(*left)[0],
1380                 &(*right)[0],
1381                 chunksize,
1382                 onumsteps,
1383                 numsteps,
1384                 resultStep,
1385                 leftstep,
1386                 rightstep,
1387                 oleftstep,
1388                 orightstep,
1389                 lroffset,
1390                 rroffset,
1391                 escript::ES_optype::ADD);  
1392          break;          break;
1393      case SUB:      case SUB:
1394          PROC_OP(NO_ARG,minus<double>());        escript::binaryOpVectorLazyHelper<real_t, real_t, real_t>(resultp,
1395                 &(*left)[0],
1396                 &(*right)[0],
1397                 chunksize,
1398                 onumsteps,
1399                 numsteps,
1400                 resultStep,
1401                 leftstep,
1402                 rightstep,
1403                 oleftstep,
1404                 orightstep,
1405                 lroffset,
1406                 rroffset,
1407                 escript::ES_optype::SUB);        
1408            //PROC_OP(NO_ARG,minus<double>());
1409          break;          break;
1410      case MUL:      case MUL:
1411          PROC_OP(NO_ARG,multiplies<double>());          //PROC_OP(NO_ARG,multiplies<double>());
1412          escript::binaryOpVectorLazyHelper<real_t, real_t, real_t>(resultp,
1413                 &(*left)[0],
1414                 &(*right)[0],
1415                 chunksize,
1416                 onumsteps,
1417                 numsteps,
1418                 resultStep,
1419                 leftstep,
1420                 rightstep,
1421                 oleftstep,
1422                 orightstep,
1423                 lroffset,
1424                 rroffset,
1425                 escript::ES_optype::MUL);        
1426          break;          break;
1427      case DIV:      case DIV:
1428          PROC_OP(NO_ARG,divides<double>());          //PROC_OP(NO_ARG,divides<double>());
1429          escript::binaryOpVectorLazyHelper<real_t, real_t, real_t>(resultp,
1430                 &(*left)[0],
1431                 &(*right)[0],
1432                 chunksize,
1433                 onumsteps,
1434                 numsteps,
1435                 resultStep,
1436                 leftstep,
1437                 rightstep,
1438                 oleftstep,
1439                 orightstep,
1440                 lroffset,
1441                 rroffset,
1442                 escript::ES_optype::DIV);        
1443          break;          break;
1444      case POW:      case POW:
1445         PROC_OP(double (double,double),::pow);         //PROC_OP(double (double,double),::pow);
1446          escript::binaryOpVectorLazyHelper<real_t, real_t, real_t>(resultp,
1447                 &(*left)[0],
1448                 &(*right)[0],
1449                 chunksize,
1450                 onumsteps,
1451                 numsteps,
1452                 resultStep,
1453                 leftstep,
1454                 rightstep,
1455                 oleftstep,
1456                 orightstep,
1457                 lroffset,
1458                 rroffset,
1459                 escript::ES_optype::POW);        
1460          break;          break;
1461      default:      default:
1462          throw DataException("Programmer error - resolveBinary can not resolve operator "+opToString(m_op)+".");          throw DataException("Programmer error - resolveBinary can not resolve operator "+opToString(m_op)+".");
# Line 1602  LAZYDEBUG(cout << "Resolve TensorProduct Line 1486  LAZYDEBUG(cout << "Resolve TensorProduct
1486    roffset=m_samplesize*tid;    roffset=m_samplesize*tid;
1487    size_t offset=roffset;    size_t offset=roffset;
1488    
1489    const ValueType* left=m_left->resolveNodeSample(tid, sampleNo, lroffset);    const RealVectorType* left=m_left->resolveNodeSample(tid, sampleNo, lroffset);
1490    
1491    const ValueType* right=m_right->resolveNodeSample(tid, sampleNo, rroffset);    const RealVectorType* right=m_right->resolveNodeSample(tid, sampleNo, rroffset);
1492    
1493  LAZYDEBUG(cerr << "[Left shape]=" << DataTypes::shapeToString(m_left->getShape()) << "\n[Right shape]=" << DataTypes::shapeToString(m_right->getShape()) << " result=" <<DataTypes::shapeToString(getShape()) <<  endl;  LAZYDEBUG(cerr << "[Left shape]=" << DataTypes::shapeToString(m_left->getShape()) << "\n[Right shape]=" << DataTypes::shapeToString(m_right->getShape()) << " result=" <<DataTypes::shapeToString(getShape()) <<  endl;
1494  cout << getNoValues() << endl;)  cout << getNoValues() << endl;)
# Line 1741  DataLazy::resolveGroupWorker(std::vector Line 1625  DataLazy::resolveGroupWorker(std::vector
1625    {             // it is possible that dats[0] is one of the objects which we discarded and    {             // it is possible that dats[0] is one of the objects which we discarded and
1626                  // all the other functionspaces match.                  // all the other functionspaces match.
1627          vector<DataExpanded*> dep;          vector<DataExpanded*> dep;
1628          vector<ValueType*> vecs;          vector<RealVectorType*> vecs;
1629          for (int i=0;i<work.size();++i)          for (int i=0;i<work.size();++i)
1630          {          {
1631                  dep.push_back(new DataExpanded(fs,work[i]->getShape(), ValueType(work[i]->getNoValues())));                  dep.push_back(new DataExpanded(fs,work[i]->getShape(), RealVectorType(work[i]->getNoValues())));
1632                  vecs.push_back(&(dep[i]->getVectorRW()));                  vecs.push_back(&(dep[i]->getVectorRW()));
1633          }          }
1634          int totalsamples=work[0]->getNumSamples();          int totalsamples=work[0]->getNumSamples();
1635          const ValueType* res=0; // Storage for answer          const RealVectorType* res=0; // Storage for answer
1636          int sample;          int sample;
1637          #pragma omp parallel private(sample, res)          #pragma omp parallel private(sample, res)
1638          {          {
# Line 1800  DataLazy::resolveNodeWorker() Line 1684  DataLazy::resolveNodeWorker()
1684      return m_id;      return m_id;
1685    }    }
1686          // from this point on we must have m_op!=IDENTITY and m_readytype=='E'          // from this point on we must have m_op!=IDENTITY and m_readytype=='E'
1687    DataExpanded* result=new DataExpanded(getFunctionSpace(),getShape(),  ValueType(getNoValues()));    DataExpanded* result=new DataExpanded(getFunctionSpace(),getShape(),  RealVectorType(getNoValues()));
1688    ValueType& resvec=result->getVectorRW();    RealVectorType& resvec=result->getVectorRW();
1689    DataReady_ptr resptr=DataReady_ptr(result);    DataReady_ptr resptr=DataReady_ptr(result);
1690    
1691    int sample;    int sample;
1692    int totalsamples=getNumSamples();    int totalsamples=getNumSamples();
1693    const ValueType* res=0;       // Storage for answer    const RealVectorType* res=0;       // Storage for answer
1694  LAZYDEBUG(cout << "Total number of samples=" <<totalsamples << endl;)  LAZYDEBUG(cout << "Total number of samples=" <<totalsamples << endl;)
1695    #pragma omp parallel private(sample,res)    #pragma omp parallel private(sample,res)
1696    {    {
# Line 1850  DataLazy::toString() const Line 1734  DataLazy::toString() const
1734  {  {
1735    ostringstream oss;    ostringstream oss;
1736    oss << "Lazy Data: [depth=" << m_height<< "] ";    oss << "Lazy Data: [depth=" << m_height<< "] ";
1737    switch (escriptParams.getLAZY_STR_FMT())    switch (escriptParams.getLazyStrFmt())
1738    {    {
1739    case 1:       // tree format    case 1:       // tree format
1740          oss << endl;          oss << endl;
# Line 2011  DataLazy::deepCopy() const Line 1895  DataLazy::deepCopy() const
1895    }    }
1896  }  }
1897    
1898    // For this, we don't care what op you were doing because the answer is now zero
1899    DataAbstract*
1900    DataLazy::zeroedCopy() const
1901    {
1902      return new DataLazy(m_id->zeroedCopy()->getPtr());
1903    }
1904    
1905  // There is no single, natural interpretation of getLength on DataLazy.  // There is no single, natural interpretation of getLength on DataLazy.
1906  // Instances of DataReady can look at the size of their vectors.  // Instances of DataReady can look at the size of their vectors.
# Line 2111  DataLazy::actsExpanded() const Line 2000  DataLazy::actsExpanded() const
2000          return (m_readytype=='E');          return (m_readytype=='E');
2001  }  }
2002    
2003  }       // end namespace  } // end namespace
2004    

Legend:
Removed from v.5985  
changed lines
  Added in v.6511

  ViewVC Help
Powered by ViewVC 1.1.26