/[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

revision 2496 by jfenwick, Fri Jun 26 06:09:47 2009 UTC revision 4154 by jfenwick, Tue Jan 22 09:30:23 2013 UTC
# Line 1  Line 1 
1    
2  /*******************************************************  /*****************************************************************************
3  *  *
4  * Copyright (c) 2003-2008 by University of Queensland  * Copyright (c) 2003-2013 by University of Queensland
5  * Earth Systems Science Computational Center (ESSCC)  * http://www.uq.edu.au
 * http://www.uq.edu.au/esscc  
6  *  *
7  * Primary Business: Queensland, Australia  * Primary Business: Queensland, Australia
8  * Licensed under the Open Software License version 3.0  * Licensed under the Open Software License version 3.0
9  * http://www.opensource.org/licenses/osl-3.0.php  * http://www.opensource.org/licenses/osl-3.0.php
10  *  *
11  *******************************************************/  * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12    * Development since 2012 by School of Earth Sciences
13    *
14    *****************************************************************************/
15    
16    
17  #include "DataLazy.h"  #include "DataLazy.h"
18  #ifdef USE_NETCDF  #include "esysUtils/Esys_MPI.h"
 #include <netcdfcpp.h>  
 #endif  
 #ifdef PASO_MPI  
 #include <mpi.h>  
 #endif  
19  #ifdef _OPENMP  #ifdef _OPENMP
20  #include <omp.h>  #include <omp.h>
21  #endif  #endif
# Line 30  Line 27 
27    
28  #include "EscriptParams.h"  #include "EscriptParams.h"
29    
30    #ifdef USE_NETCDF
31    #include <netcdfcpp.h>
32    #endif
33    
34  #include <iomanip>      // for some fancy formatting in debug  #include <iomanip>      // for some fancy formatting in debug
35    
36  // #define LAZYDEBUG(X) if (privdebug){X;}  // #define LAZYDEBUG(X) if (privdebug){X;}
# Line 42  bool privdebug=false; Line 43  bool privdebug=false;
43  #define DISABLEDEBUG privdebug=false;  #define DISABLEDEBUG privdebug=false;
44  }  }
45    
46  // #define SIZELIMIT  // #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();}
47  // #define SIZELIMIT if ((m_height>7) || (m_children>127)) {cerr << "\n!!!!!!! SIZE LIMIT EXCEEDED " << m_children << ";" << m_height << endl << toString() << endl; resolveToIdentity();}  
48  //#define SIZELIMIT if ((m_height>7) || (m_children>127)) {resolveToIdentity();}  // #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();}
49    
50  #define SIZELIMIT if ((m_height>escript::escriptParams.getTOO_MANY_LEVELS()) || (m_children>escript::escriptParams.getTOO_MANY_NODES())) {resolveToIdentity();}  
51    #define SIZELIMIT if (m_height>escript::escriptParams.getTOO_MANY_LEVELS())  {if (escript::escriptParams.getLAZY_VERBOSE()){cerr << "SIZE LIMIT EXCEEDED height=" << m_height << endl;}resolveToIdentity();}
52    
53  /*  /*
54  How does DataLazy work?  How does DataLazy work?
# Line 110  namespace escript Line 112  namespace escript
112  namespace  namespace
113  {  {
114    
115    
116    // enabling this will print out when ever the maximum stacksize used by resolve increases
117    // it assumes _OPENMP is also in use
118    //#define LAZY_STACK_PROF
119    
120    
121    
122    #ifndef _OPENMP
123      #ifdef LAZY_STACK_PROF
124      #undef LAZY_STACK_PROF
125      #endif
126    #endif
127    
128    
129    #ifdef LAZY_STACK_PROF
130    std::vector<void*> stackstart(getNumberOfThreads());
131    std::vector<void*> stackend(getNumberOfThreads());
132    size_t maxstackuse=0;
133    #endif
134    
135  enum ES_opgroup  enum ES_opgroup
136  {  {
137     G_UNKNOWN,     G_UNKNOWN,
# Line 120  enum ES_opgroup Line 142  enum ES_opgroup
142     G_NP1OUT,        // non-pointwise op with one output     G_NP1OUT,        // non-pointwise op with one output
143     G_NP1OUT_P,      // non-pointwise op with one output requiring a parameter     G_NP1OUT_P,      // non-pointwise op with one output requiring a parameter
144     G_TENSORPROD,    // general tensor product     G_TENSORPROD,    // general tensor product
145     G_NP1OUT_2P      // non-pointwise op with one output requiring two params     G_NP1OUT_2P,     // non-pointwise op with one output requiring two params
146       G_REDUCTION,     // non-pointwise unary op with a scalar output
147       G_CONDEVAL
148  };  };
149    
150    
# Line 135  string ES_opstrings[]={"UNKNOWN","IDENTI Line 159  string ES_opstrings[]={"UNKNOWN","IDENTI
159              "symmetric","nonsymmetric",              "symmetric","nonsymmetric",
160              "prod",              "prod",
161              "transpose", "trace",              "transpose", "trace",
162              "swapaxes"};              "swapaxes",
163  int ES_opcount=41;              "minval", "maxval",
164                "condEval"};
165    int ES_opcount=44;
166  ES_opgroup opgroups[]={G_UNKNOWN,G_IDENTITY,G_BINARY,G_BINARY,G_BINARY,G_BINARY, G_BINARY,  ES_opgroup opgroups[]={G_UNKNOWN,G_IDENTITY,G_BINARY,G_BINARY,G_BINARY,G_BINARY, G_BINARY,
167              G_UNARY,G_UNARY,G_UNARY, //10              G_UNARY,G_UNARY,G_UNARY, //10
168              G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,    // 17              G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,    // 17
# Line 146  ES_opgroup opgroups[]={G_UNKNOWN,G_IDENT Line 172  ES_opgroup opgroups[]={G_UNKNOWN,G_IDENT
172              G_NP1OUT,G_NP1OUT,              G_NP1OUT,G_NP1OUT,
173              G_TENSORPROD,              G_TENSORPROD,
174              G_NP1OUT_P, G_NP1OUT_P,              G_NP1OUT_P, G_NP1OUT_P,
175              G_NP1OUT_2P};              G_NP1OUT_2P,
176                G_REDUCTION, G_REDUCTION,
177                G_CONDEVAL};
178  inline  inline
179  ES_opgroup  ES_opgroup
180  getOpgroup(ES_optype op)  getOpgroup(ES_optype op)
# Line 191  resultShape(DataAbstract_ptr left, DataA Line 219  resultShape(DataAbstract_ptr left, DataA
219        {        {
220          throw DataException("Shapes not the name - shapes must match for (point)binary operations.");          throw DataException("Shapes not the name - shapes must match for (point)binary operations.");
221        }        }
222    
223        if (left->getRank()==0)   // we need to allow scalar * anything        if (left->getRank()==0)   // we need to allow scalar * anything
224        {        {
225          return right->getShape();          return right->getShape();
# Line 218  resultShape(DataAbstract_ptr left, ES_op Line 247  resultShape(DataAbstract_ptr left, ES_op
247          int rank=left->getRank();          int rank=left->getRank();
248          if (axis_offset<0 || axis_offset>rank)          if (axis_offset<0 || axis_offset>rank)
249          {          {
250                 throw DataException("Error - Data::transpose must have 0 <= axis_offset <= rank=" + rank);              stringstream e;
251              }              e << "Error - Data::transpose must have 0 <= axis_offset <= rank=" << rank;
252              for (int i=0; i<rank; i++)              throw DataException(e.str());
253            }
254            for (int i=0; i<rank; i++)
255          {          {
256             int index = (axis_offset+i)%rank;             int index = (axis_offset+i)%rank;
257                 sh.push_back(s[index]); // Append to new shape             sh.push_back(s[index]); // Append to new shape
258              }          }
259          return sh;          return sh;
260         }         }
261      break;      break;
# Line 303  SwapShape(DataAbstract_ptr left, const i Line 334  SwapShape(DataAbstract_ptr left, const i
334          throw DataException("Error - Data::swapaxes argument must have at least rank 2.");          throw DataException("Error - Data::swapaxes argument must have at least rank 2.");
335       }       }
336       if (axis0<0 || axis0>rank-1) {       if (axis0<0 || axis0>rank-1) {
337          throw DataException("Error - Data::swapaxes: axis0 must be between 0 and rank-1=" + rank-1);          stringstream e;
338            e << "Error - Data::swapaxes: axis0 must be between 0 and rank-1=" << (rank-1);
339            throw DataException(e.str());
340       }       }
341       if (axis1<0 || axis1>rank-1) {       if (axis1<0 || axis1>rank-1) {
342           throw DataException("Error - Data::swapaxes: axis1 must be between 0 and rank-1=" + rank-1);          stringstream e;
343            e << "Error - Data::swapaxes: axis1 must be between 0 and rank-1=" << (rank-1);
344            throw DataException(e.str());
345       }       }
346       if (axis0 == axis1) {       if (axis0 == axis1) {
347           throw DataException("Error - Data::swapaxes: axis indices must be different.");           throw DataException("Error - Data::swapaxes: axis indices must be different.");
# Line 395  GTPShape(DataAbstract_ptr left, DataAbst Line 430  GTPShape(DataAbstract_ptr left, DataAbst
430    return shape2;    return shape2;
431  }  }
432    
 // determine the number of samples requires to evaluate an expression combining left and right  
 // NP1OUT needs an extra buffer because we can't write the answers over the top of the input.  
 // The same goes for G_TENSORPROD  
 // It might seem that pointwise binary ops (G_BINARY) could be written over the top of the lefts.  
 // This would be true were it not for the possibility that the LHS could be a scalar which needs to be examined  
 // multiple times  
 int  
 calcBuffs(const DataLazy_ptr& left, const DataLazy_ptr& right, ES_optype op)  
 {  
    switch(getOpgroup(op))  
    {  
    case G_IDENTITY: return 1;  
    case G_BINARY: return 1+max(left->getBuffsRequired(),right->getBuffsRequired()+1);  
    case G_UNARY:  
    case G_UNARY_P: return max(left->getBuffsRequired(),1);  
    case G_NP1OUT: return 1+max(left->getBuffsRequired(),1);  
    case G_NP1OUT_P: return 1+max(left->getBuffsRequired(),1);  
    case G_TENSORPROD: return 1+max(left->getBuffsRequired(),right->getBuffsRequired()+1);  
    case G_NP1OUT_2P: return 1+max(left->getBuffsRequired(),1);  
    default:  
     throw DataException("Programmer Error - attempt to calcBuffs() for operator "+opToString(op)+".");  
    }  
 }  
   
   
433  }   // end anonymous namespace  }   // end anonymous namespace
434    
435    
# Line 435  opToString(ES_optype op) Line 445  opToString(ES_optype op)
445    return ES_opstrings[op];    return ES_opstrings[op];
446  }  }
447    
448    void DataLazy::LazyNodeSetup()
449    {
450    #ifdef _OPENMP
451        int numthreads=omp_get_max_threads();
452        m_samples.resize(numthreads*m_samplesize);
453        m_sampleids=new int[numthreads];
454        for (int i=0;i<numthreads;++i)
455        {
456            m_sampleids[i]=-1;  
457        }
458    #else
459        m_samples.resize(m_samplesize);
460        m_sampleids=new int[1];
461        m_sampleids[0]=-1;
462    #endif  // _OPENMP
463    }
464    
465    
466    // Creates an identity node
467  DataLazy::DataLazy(DataAbstract_ptr p)  DataLazy::DataLazy(DataAbstract_ptr p)
468      : parent(p->getFunctionSpace(),p->getShape())      : parent(p->getFunctionSpace(),p->getShape())
469        ,m_sampleids(0),
470        m_samples(1)
471  {  {
472     if (p->isLazy())     if (p->isLazy())
473     {     {
# Line 456  LAZYDEBUG(cout << "Wrapping " << dr.get( Line 486  LAZYDEBUG(cout << "Wrapping " << dr.get(
486  LAZYDEBUG(cout << "(1)Lazy created with " << m_samplesize << endl;)  LAZYDEBUG(cout << "(1)Lazy created with " << m_samplesize << endl;)
487  }  }
488    
   
   
   
489  DataLazy::DataLazy(DataAbstract_ptr left, ES_optype op)  DataLazy::DataLazy(DataAbstract_ptr left, ES_optype op)
490      : parent(left->getFunctionSpace(),left->getShape()),      : parent(left->getFunctionSpace(),(getOpgroup(op)!=G_REDUCTION)?left->getShape():DataTypes::scalarShape),
491      m_op(op),      m_op(op),
492      m_axis_offset(0),      m_axis_offset(0),
493      m_transpose(0),      m_transpose(0),
494      m_SL(0), m_SM(0), m_SR(0)      m_SL(0), m_SM(0), m_SR(0)
495  {  {
496     if ((getOpgroup(op)!=G_UNARY) && (getOpgroup(op)!=G_NP1OUT))     if ((getOpgroup(op)!=G_UNARY) && (getOpgroup(op)!=G_NP1OUT) && (getOpgroup(op)!=G_REDUCTION))
497     {     {
498      throw DataException("Programmer error - constructor DataLazy(left, op) will only process UNARY operations.");      throw DataException("Programmer error - constructor DataLazy(left, op) will only process UNARY operations.");
499     }     }
# Line 482  DataLazy::DataLazy(DataAbstract_ptr left Line 509  DataLazy::DataLazy(DataAbstract_ptr left
509     }     }
510     m_readytype=lleft->m_readytype;     m_readytype=lleft->m_readytype;
511     m_left=lleft;     m_left=lleft;
    m_buffsRequired=calcBuffs(m_left, m_right,m_op); // yeah m_right will be null at this point  
512     m_samplesize=getNumDPPSample()*getNoValues();     m_samplesize=getNumDPPSample()*getNoValues();
    m_maxsamplesize=max(m_samplesize,m_left->getMaxSampleSize());  
513     m_children=m_left->m_children+1;     m_children=m_left->m_children+1;
514     m_height=m_left->m_height+1;     m_height=m_left->m_height+1;
515       LazyNodeSetup();
516     SIZELIMIT     SIZELIMIT
517  }  }
518    
# Line 553  LAZYDEBUG(cout << "Right " << right.get( Line 579  LAZYDEBUG(cout << "Right " << right.get(
579      m_readytype='C';      m_readytype='C';
580     }     }
581     m_samplesize=getNumDPPSample()*getNoValues();     m_samplesize=getNumDPPSample()*getNoValues();
    m_maxsamplesize=max(max(m_samplesize,m_right->getMaxSampleSize()),m_left->getMaxSampleSize());    
    m_buffsRequired=calcBuffs(m_left, m_right,m_op);  
582     m_children=m_left->m_children+m_right->m_children+2;     m_children=m_left->m_children+m_right->m_children+2;
583     m_height=max(m_left->m_height,m_right->m_height)+1;     m_height=max(m_left->m_height,m_right->m_height)+1;
584       LazyNodeSetup();
585     SIZELIMIT     SIZELIMIT
586  LAZYDEBUG(cout << "(3)Lazy created with " << m_samplesize << endl;)  LAZYDEBUG(cout << "(3)Lazy created with " << m_samplesize << endl;)
587  }  }
# Line 620  DataLazy::DataLazy(DataAbstract_ptr left Line 645  DataLazy::DataLazy(DataAbstract_ptr left
645      m_readytype='C';      m_readytype='C';
646     }     }
647     m_samplesize=getNumDPPSample()*getNoValues();     m_samplesize=getNumDPPSample()*getNoValues();
    m_maxsamplesize=max(max(m_samplesize,m_right->getMaxSampleSize()),m_left->getMaxSampleSize());    
    m_buffsRequired=calcBuffs(m_left, m_right,m_op);  
648     m_children=m_left->m_children+m_right->m_children+2;     m_children=m_left->m_children+m_right->m_children+2;
649     m_height=max(m_left->m_height,m_right->m_height)+1;     m_height=max(m_left->m_height,m_right->m_height)+1;
650       LazyNodeSetup();
651     SIZELIMIT     SIZELIMIT
652  LAZYDEBUG(cout << "(4)Lazy created with " << m_samplesize << endl;)  LAZYDEBUG(cout << "(4)Lazy created with " << m_samplesize << endl;)
653  }  }
# Line 651  DataLazy::DataLazy(DataAbstract_ptr left Line 675  DataLazy::DataLazy(DataAbstract_ptr left
675     }     }
676     m_readytype=lleft->m_readytype;     m_readytype=lleft->m_readytype;
677     m_left=lleft;     m_left=lleft;
    m_buffsRequired=calcBuffs(m_left, m_right,m_op); // yeah m_right will be null at this point  
678     m_samplesize=getNumDPPSample()*getNoValues();     m_samplesize=getNumDPPSample()*getNoValues();
    m_maxsamplesize=max(m_samplesize,m_left->getMaxSampleSize());  
679     m_children=m_left->m_children+1;     m_children=m_left->m_children+1;
680     m_height=m_left->m_height+1;     m_height=m_left->m_height+1;
681       LazyNodeSetup();
682     SIZELIMIT     SIZELIMIT
683  LAZYDEBUG(cout << "(5)Lazy created with " << m_samplesize << endl;)  LAZYDEBUG(cout << "(5)Lazy created with " << m_samplesize << endl;)
684  }  }
# Line 682  DataLazy::DataLazy(DataAbstract_ptr left Line 705  DataLazy::DataLazy(DataAbstract_ptr left
705     }     }
706     m_readytype=lleft->m_readytype;     m_readytype=lleft->m_readytype;
707     m_left=lleft;     m_left=lleft;
    m_buffsRequired=calcBuffs(m_left, m_right,m_op); // yeah m_right will be null at this point  
708     m_samplesize=getNumDPPSample()*getNoValues();     m_samplesize=getNumDPPSample()*getNoValues();
    m_maxsamplesize=max(m_samplesize,m_left->getMaxSampleSize());  
709     m_children=m_left->m_children+1;     m_children=m_left->m_children+1;
710     m_height=m_left->m_height+1;     m_height=m_left->m_height+1;
711       LazyNodeSetup();
712     SIZELIMIT     SIZELIMIT
713  LAZYDEBUG(cout << "(6)Lazy created with " << m_samplesize << endl;)  LAZYDEBUG(cout << "(6)Lazy created with " << m_samplesize << endl;)
714  }  }
# Line 714  DataLazy::DataLazy(DataAbstract_ptr left Line 736  DataLazy::DataLazy(DataAbstract_ptr left
736     }     }
737     m_readytype=lleft->m_readytype;     m_readytype=lleft->m_readytype;
738     m_left=lleft;     m_left=lleft;
    m_buffsRequired=calcBuffs(m_left, m_right,m_op); // yeah m_right will be null at this point  
739     m_samplesize=getNumDPPSample()*getNoValues();     m_samplesize=getNumDPPSample()*getNoValues();
    m_maxsamplesize=max(m_samplesize,m_left->getMaxSampleSize());  
740     m_children=m_left->m_children+1;     m_children=m_left->m_children+1;
741     m_height=m_left->m_height+1;     m_height=m_left->m_height+1;
742       LazyNodeSetup();
743     SIZELIMIT     SIZELIMIT
744  LAZYDEBUG(cout << "(7)Lazy created with " << m_samplesize << endl;)  LAZYDEBUG(cout << "(7)Lazy created with " << m_samplesize << endl;)
745  }  }
746    
747  DataLazy::~DataLazy()  
748    namespace
749  {  {
 }  
750    
751        inline int max3(int a, int b, int c)
752        {
753        int t=(a>b?a:b);
754        return (t>c?t:c);
755    
756  int      }
 DataLazy::getBuffsRequired() const  
 {  
     return m_buffsRequired;  
757  }  }
758    
759    DataLazy::DataLazy(DataAbstract_ptr mask, DataAbstract_ptr left, DataAbstract_ptr right/*, double tol*/)
760  size_t      : parent(left->getFunctionSpace(), left->getShape()),
761  DataLazy::getMaxSampleSize() const      m_op(CONDEVAL),
762        m_axis_offset(0),
763        m_transpose(0),
764        m_tol(0)
765  {  {
766      return m_maxsamplesize;  
767       DataLazy_ptr lmask;
768       DataLazy_ptr lleft;
769       DataLazy_ptr lright;
770       if (!mask->isLazy())
771       {
772        lmask=DataLazy_ptr(new DataLazy(mask));
773       }
774       else
775       {
776        lmask=dynamic_pointer_cast<DataLazy>(mask);
777       }
778       if (!left->isLazy())
779       {
780        lleft=DataLazy_ptr(new DataLazy(left));
781       }
782       else
783       {
784        lleft=dynamic_pointer_cast<DataLazy>(left);
785       }
786       if (!right->isLazy())
787       {
788        lright=DataLazy_ptr(new DataLazy(right));
789       }
790       else
791       {
792        lright=dynamic_pointer_cast<DataLazy>(right);
793       }
794       m_readytype=lmask->m_readytype;
795       if ((lleft->m_readytype!=lright->m_readytype) || (lmask->m_readytype!=lleft->m_readytype))
796       {
797        throw DataException("Programmer Error - condEval arguments must have the same readytype");
798       }
799       m_left=lleft;
800       m_right=lright;
801       m_mask=lmask;
802       m_samplesize=getNumDPPSample()*getNoValues();
803       m_children=m_left->m_children+m_right->m_children+m_mask->m_children+1;
804       m_height=max3(m_left->m_height,m_right->m_height,m_mask->m_height)+1;
805       LazyNodeSetup();
806       SIZELIMIT
807    LAZYDEBUG(cout << "(8)Lazy created with " << m_samplesize << endl;)
808  }  }
809    
810    
811    
812  size_t  DataLazy::~DataLazy()
 DataLazy::getSampleBufferSize() const  
813  {  {
814      return m_maxsamplesize*(max(1,m_buffsRequired));     delete[] m_sampleids;
815  }  }
816    
817    
818  /*  /*
819    \brief Evaluates the expression using methods on Data.    \brief Evaluates the expression using methods on Data.
820    This does the work for the collapse method.    This does the work for the collapse method.
# Line 891  DataLazy::collapseToReady() Line 957  DataLazy::collapseToReady()
957      case SWAP:      case SWAP:
958      result=left.swapaxes(m_axis_offset, m_transpose);      result=left.swapaxes(m_axis_offset, m_transpose);
959      break;      break;
960        case MINVAL:
961        result=left.minval();
962        break;
963        case MAXVAL:
964        result=left.minval();
965        break;
966      default:      default:
967      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)+".");
968    }    }
# Line 918  DataLazy::collapse() Line 990  DataLazy::collapse()
990    m_op=IDENTITY;    m_op=IDENTITY;
991  }  }
992    
993  /*  
994    \brief Compute the value of the expression (unary operation) for the given sample.  
995    \return Vector which stores the value of the subexpression for the given sample.  
996    \param v A vector to store intermediate results.  
997    \param offset Index in v to begin storing results.  
998    \param sampleNo Sample number to evaluate.  #define PROC_OP(TYPE,X)                               \
999    \param roffset (output parameter) the offset in the return vector where the result begins.      for (int j=0;j<onumsteps;++j)\
1000        {\
1001    The return value will be an existing vector so do not deallocate it.        for (int i=0;i<numsteps;++i,resultp+=resultStep) \
1002    If the result is stored in v it should be stored at the offset given.        { \
1003    Everything from offset to the end of v should be considered available for this method to use.  LAZYDEBUG(cout << "[left,right]=[" << lroffset << "," << rroffset << "]" << endl;)\
1004  */  LAZYDEBUG(cout << "{left,right}={" << (*left)[lroffset] << "," << (*right)[rroffset] << "}\n";)\
1005  DataTypes::ValueType*           tensor_binary_operation< TYPE >(chunksize, &((*left)[lroffset]), &((*right)[rroffset]), resultp, X); \
1006  DataLazy::resolveUnary(ValueType& v, size_t offset, int sampleNo, size_t& roffset) const  LAZYDEBUG(cout << " result=      " << resultp[0] << endl;) \
1007             lroffset+=leftstep; \
1008             rroffset+=rightstep; \
1009          }\
1010          lroffset+=oleftstep;\
1011          rroffset+=orightstep;\
1012        }
1013    
1014    
1015    // The result will be stored in m_samples
1016    // The return value is a pointer to the DataVector, offset is the offset within the return value
1017    const DataTypes::ValueType*
1018    DataLazy::resolveNodeSample(int tid, int sampleNo, size_t& roffset)
1019    {
1020    LAZYDEBUG(cout << "Resolve sample " << toString() << endl;)
1021        // collapse so we have a 'E' node or an IDENTITY for some other type
1022      if (m_readytype!='E' && m_op!=IDENTITY)
1023      {
1024        collapse();
1025      }
1026      if (m_op==IDENTITY)  
1027      {
1028        const ValueType& vec=m_id->getVectorRO();
1029        roffset=m_id->getPointOffset(sampleNo, 0);
1030    #ifdef LAZY_STACK_PROF
1031    int x;
1032    if (&x<stackend[omp_get_thread_num()])
1033    {
1034           stackend[omp_get_thread_num()]=&x;
1035    }
1036    #endif
1037        return &(vec);
1038      }
1039      if (m_readytype!='E')
1040      {
1041        throw DataException("Programmer Error - Collapse did not produce an expanded node.");
1042      }
1043      if (m_sampleids[tid]==sampleNo)
1044      {
1045        roffset=tid*m_samplesize;
1046        return &(m_samples);        // sample is already resolved
1047      }
1048      m_sampleids[tid]=sampleNo;
1049    
1050      switch (getOpgroup(m_op))
1051      {
1052      case G_UNARY:
1053      case G_UNARY_P: return resolveNodeUnary(tid, sampleNo, roffset);
1054      case G_BINARY: return resolveNodeBinary(tid, sampleNo, roffset);
1055      case G_NP1OUT: return resolveNodeNP1OUT(tid, sampleNo, roffset);
1056      case G_NP1OUT_P: return resolveNodeNP1OUT_P(tid, sampleNo, roffset);
1057      case G_TENSORPROD: return resolveNodeTProd(tid, sampleNo, roffset);
1058      case G_NP1OUT_2P: return resolveNodeNP1OUT_2P(tid, sampleNo, roffset);
1059      case G_REDUCTION: return resolveNodeReduction(tid, sampleNo, roffset);
1060      case G_CONDEVAL: return resolveNodeCondEval(tid, sampleNo, roffset);
1061      default:
1062        throw DataException("Programmer Error - resolveSample does not know how to process "+opToString(m_op)+".");
1063      }
1064    }
1065    
1066    const DataTypes::ValueType*
1067    DataLazy::resolveNodeUnary(int tid, int sampleNo, size_t& roffset)
1068  {  {
1069      // we assume that any collapsing has been done before we get here      // we assume that any collapsing has been done before we get here
1070      // since we only have one argument we don't need to think about only      // since we only have one argument we don't need to think about only
1071      // processing single points.      // processing single points.
1072        // we will also know we won't get identity nodes
1073    if (m_readytype!='E')    if (m_readytype!='E')
1074    {    {
1075      throw DataException("Programmer error - resolveUnary should only be called on expanded Data.");      throw DataException("Programmer error - resolveUnary should only be called on expanded Data.");
1076    }    }
1077    const ValueType* vleft=m_left->resolveSample(v,offset,sampleNo,roffset);    if (m_op==IDENTITY)
1078    const double* left=&((*vleft)[roffset]);    {
1079    double* result=&(v[offset]);      throw DataException("Programmer error - resolveNodeUnary should not be called on identity nodes.");
1080    roffset=offset;    }
1081      const DataTypes::ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, roffset);
1082      const double* left=&((*leftres)[roffset]);
1083      roffset=m_samplesize*tid;
1084      double* result=&(m_samples[roffset]);
1085    switch (m_op)    switch (m_op)
1086    {    {
1087      case SIN:        case SIN:  
# Line 1053  DataLazy::resolveUnary(ValueType& v, siz Line 1191  DataLazy::resolveUnary(ValueType& v, siz
1191      default:      default:
1192      throw DataException("Programmer error - resolveUnary can not resolve operator "+opToString(m_op)+".");      throw DataException("Programmer error - resolveUnary can not resolve operator "+opToString(m_op)+".");
1193    }    }
1194    return &v;    return &(m_samples);
1195  }  }
1196    
1197    
1198    const DataTypes::ValueType*
1199    DataLazy::resolveNodeReduction(int tid, int sampleNo, size_t& roffset)
1200    {
1201        // we assume that any collapsing has been done before we get here
1202        // since we only have one argument we don't need to think about only
1203        // processing single points.
1204        // we will also know we won't get identity nodes
1205      if (m_readytype!='E')
1206      {
1207        throw DataException("Programmer error - resolveUnary should only be called on expanded Data.");
1208      }
1209      if (m_op==IDENTITY)
1210      {
1211        throw DataException("Programmer error - resolveNodeUnary should not be called on identity nodes.");
1212      }
1213      size_t loffset=0;
1214      const DataTypes::ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, loffset);
1215    
1216      roffset=m_samplesize*tid;
1217      unsigned int ndpps=getNumDPPSample();
1218      unsigned int psize=DataTypes::noValues(m_left->getShape());
1219      double* result=&(m_samples[roffset]);
1220      switch (m_op)
1221      {
1222        case MINVAL:
1223        {
1224          for (unsigned int z=0;z<ndpps;++z)
1225          {
1226            FMin op;
1227            *result=DataMaths::reductionOp(*leftres, m_left->getShape(), loffset, op, numeric_limits<double>::max());
1228            loffset+=psize;
1229            result++;
1230          }
1231        }
1232        break;
1233        case MAXVAL:
1234        {
1235          for (unsigned int z=0;z<ndpps;++z)
1236          {
1237          FMax op;
1238          *result=DataMaths::reductionOp(*leftres, m_left->getShape(), loffset, op, numeric_limits<double>::max()*-1);
1239          loffset+=psize;
1240          result++;
1241          }
1242        }
1243        break;
1244        default:
1245        throw DataException("Programmer error - resolveUnary can not resolve operator "+opToString(m_op)+".");
1246      }
1247      return &(m_samples);
1248    }
1249    
1250    const DataTypes::ValueType*
1251    DataLazy::resolveNodeNP1OUT(int tid, int sampleNo, size_t& roffset)
 /*  
   \brief Compute the value of the expression (unary operation) for the given sample.  
   \return Vector which stores the value of the subexpression for the given sample.  
   \param v A vector to store intermediate results.  
   \param offset Index in v to begin storing results.  
   \param sampleNo Sample number to evaluate.  
   \param roffset (output parameter) the offset in the return vector where the result begins.  
   
   The return value will be an existing vector so do not deallocate it.  
   If the result is stored in v it should be stored at the offset given.  
   Everything from offset to the end of v should be considered available for this method to use.  
 */  
 DataTypes::ValueType*  
 DataLazy::resolveNP1OUT(ValueType& v, size_t offset, int sampleNo, size_t& roffset) const  
1252  {  {
1253      // we assume that any collapsing has been done before we get here      // we assume that any collapsing has been done before we get here
1254      // since we only have one argument we don't need to think about only      // since we only have one argument we don't need to think about only
1255      // processing single points.      // processing single points.
1256    if (m_readytype!='E')    if (m_readytype!='E')
1257    {    {
1258      throw DataException("Programmer error - resolveNP1OUT should only be called on expanded Data.");      throw DataException("Programmer error - resolveNodeNP1OUT should only be called on expanded Data.");
1259    }    }
1260      // since we can't write the result over the input, we need a result offset further along    if (m_op==IDENTITY)
1261    size_t subroffset=roffset+m_samplesize;    {
1262  LAZYDEBUG(cerr << "subroffset=" << subroffset << endl;)      throw DataException("Programmer error - resolveNodeNP1OUT should not be called on identity nodes.");
1263    const ValueType* vleft=m_left->resolveSample(v,offset+m_left->m_samplesize,sampleNo,subroffset);    }
1264    roffset=offset;    size_t subroffset;
1265      const ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1266      roffset=m_samplesize*tid;
1267    size_t loop=0;    size_t loop=0;
1268    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;
1269    size_t step=getNoValues();    size_t step=getNoValues();
1270      size_t offset=roffset;
1271    switch (m_op)    switch (m_op)
1272    {    {
1273      case SYM:      case SYM:
1274      for (loop=0;loop<numsteps;++loop)      for (loop=0;loop<numsteps;++loop)
1275      {      {
1276          DataMaths::symmetric(*vleft,m_left->getShape(),subroffset, v, getShape(), offset);          DataMaths::symmetric(*leftres,m_left->getShape(),subroffset, m_samples, getShape(), offset);
1277          subroffset+=step;          subroffset+=step;
1278          offset+=step;          offset+=step;
1279      }      }
# Line 1104  LAZYDEBUG(cerr << "subroffset=" << subro Line 1281  LAZYDEBUG(cerr << "subroffset=" << subro
1281      case NSYM:      case NSYM:
1282      for (loop=0;loop<numsteps;++loop)      for (loop=0;loop<numsteps;++loop)
1283      {      {
1284          DataMaths::nonsymmetric(*vleft,m_left->getShape(),subroffset, v, getShape(), offset);          DataMaths::nonsymmetric(*leftres,m_left->getShape(),subroffset, m_samples, getShape(), offset);
1285          subroffset+=step;          subroffset+=step;
1286          offset+=step;          offset+=step;
1287      }      }
# Line 1112  LAZYDEBUG(cerr << "subroffset=" << subro Line 1289  LAZYDEBUG(cerr << "subroffset=" << subro
1289      default:      default:
1290      throw DataException("Programmer error - resolveNP1OUT can not resolve operator "+opToString(m_op)+".");      throw DataException("Programmer error - resolveNP1OUT can not resolve operator "+opToString(m_op)+".");
1291    }    }
1292    return &v;    return &m_samples;
1293  }  }
1294    
1295  /*  const DataTypes::ValueType*
1296    \brief Compute the value of the expression (unary operation) for the given sample.  DataLazy::resolveNodeNP1OUT_P(int tid, int sampleNo, size_t& roffset)
   \return Vector which stores the value of the subexpression for the given sample.  
   \param v A vector to store intermediate results.  
   \param offset Index in v to begin storing results.  
   \param sampleNo Sample number to evaluate.  
   \param roffset (output parameter) the offset in the return vector where the result begins.  
   
   The return value will be an existing vector so do not deallocate it.  
   If the result is stored in v it should be stored at the offset given.  
   Everything from offset to the end of v should be considered available for this method to use.  
 */  
 DataTypes::ValueType*  
 DataLazy::resolveNP1OUT_P(ValueType& v, size_t offset, int sampleNo, size_t& roffset) const  
1297  {  {
1298      // we assume that any collapsing has been done before we get here      // we assume that any collapsing has been done before we get here
1299      // since we only have one argument we don't need to think about only      // since we only have one argument we don't need to think about only
1300      // processing single points.      // processing single points.
1301    if (m_readytype!='E')    if (m_readytype!='E')
1302    {    {
1303      throw DataException("Programmer error - resolveNP1OUT_P should only be called on expanded Data.");      throw DataException("Programmer error - resolveNodeNP1OUT_P should only be called on expanded Data.");
1304      }
1305      if (m_op==IDENTITY)
1306      {
1307        throw DataException("Programmer error - resolveNodeNP1OUT_P should not be called on identity nodes.");
1308    }    }
     // since we can't write the result over the input, we need a result offset further along  
1309    size_t subroffset;    size_t subroffset;
1310    const ValueType* vleft=m_left->resolveSample(v,offset+m_left->m_samplesize,sampleNo,subroffset);    size_t offset;
1311  LAZYDEBUG(cerr << "srcsamplesize=" << offset+m_left->m_samplesize << " beg=" << subroffset << endl;)    const ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1312  LAZYDEBUG(cerr << "Offset for 5800=" << getPointOffset(5800/getNumDPPSample(),5800%getNumDPPSample()) << endl;)    roffset=m_samplesize*tid;
1313    roffset=offset;    offset=roffset;
1314    size_t loop=0;    size_t loop=0;
1315    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;
1316    size_t outstep=getNoValues();    size_t outstep=getNoValues();
1317    size_t instep=m_left->getNoValues();    size_t instep=m_left->getNoValues();
 LAZYDEBUG(cerr << "instep=" << instep << " outstep=" << outstep<< " numsteps=" << numsteps << endl;)  
1318    switch (m_op)    switch (m_op)
1319    {    {
1320      case TRACE:      case TRACE:
1321      for (loop=0;loop<numsteps;++loop)      for (loop=0;loop<numsteps;++loop)
1322      {      {
1323  size_t zz=sampleNo*getNumDPPSample()+loop;              DataMaths::trace(*leftres,m_left->getShape(),subroffset, m_samples ,getShape(),offset,m_axis_offset);
 if (zz==5800)  
 {  
 LAZYDEBUG(cerr << "point=" <<  zz<< endl;)  
 LAZYDEBUG(cerr << "Input to  trace=" << DataTypes::pointToString(*vleft,m_left->getShape(),subroffset,"") << endl;)  
 LAZYDEBUG(cerr << "Offset for point=" << getPointOffset(5800/getNumDPPSample(),5800%getNumDPPSample()) << " vs ";)  
 LAZYDEBUG(cerr << subroffset << endl;)  
 LAZYDEBUG(cerr << "output=" << offset << endl;)  
 }  
             DataMaths::trace(*vleft,m_left->getShape(),subroffset, v ,getShape(),offset,m_axis_offset);  
 if (zz==5800)  
 {  
 LAZYDEBUG(cerr << "Result of trace=" << DataTypes::pointToString(v,getShape(),offset,"") << endl;)  
 }  
1324          subroffset+=instep;          subroffset+=instep;
1325          offset+=outstep;          offset+=outstep;
1326      }      }
# Line 1174  LAZYDEBUG(cerr << "Result of trace=" << Line 1328  LAZYDEBUG(cerr << "Result of trace=" <<
1328      case TRANS:      case TRANS:
1329      for (loop=0;loop<numsteps;++loop)      for (loop=0;loop<numsteps;++loop)
1330      {      {
1331              DataMaths::transpose(*vleft,m_left->getShape(),subroffset, v,getShape(),offset,m_axis_offset);              DataMaths::transpose(*leftres,m_left->getShape(),subroffset, m_samples, getShape(),offset,m_axis_offset);
1332          subroffset+=instep;          subroffset+=instep;
1333          offset+=outstep;          offset+=outstep;
1334      }      }
# Line 1182  LAZYDEBUG(cerr << "Result of trace=" << Line 1336  LAZYDEBUG(cerr << "Result of trace=" <<
1336      default:      default:
1337      throw DataException("Programmer error - resolveNP1OUTP can not resolve operator "+opToString(m_op)+".");      throw DataException("Programmer error - resolveNP1OUTP can not resolve operator "+opToString(m_op)+".");
1338    }    }
1339    return &v;    return &m_samples;
1340  }  }
1341    
1342    
1343  /*  const DataTypes::ValueType*
1344    \brief Compute the value of the expression (unary operation with int params) for the given sample.  DataLazy::resolveNodeNP1OUT_2P(int tid, int sampleNo, size_t& roffset)
   \return Vector which stores the value of the subexpression for the given sample.  
   \param v A vector to store intermediate results.  
   \param offset Index in v to begin storing results.  
   \param sampleNo Sample number to evaluate.  
   \param roffset (output parameter) the offset in the return vector where the result begins.  
   
   The return value will be an existing vector so do not deallocate it.  
   If the result is stored in v it should be stored at the offset given.  
   Everything from offset to the end of v should be considered available for this method to use.  
 */  
 DataTypes::ValueType*  
 DataLazy::resolveNP1OUT_2P(ValueType& v, size_t offset, int sampleNo, size_t& roffset) const  
1345  {  {
     // we assume that any collapsing has been done before we get here  
     // since we only have one argument we don't need to think about only  
     // processing single points.  
1346    if (m_readytype!='E')    if (m_readytype!='E')
1347    {    {
1348      throw DataException("Programmer error - resolveNP1OUT_2P should only be called on expanded Data.");      throw DataException("Programmer error - resolveNodeNP1OUT_2P should only be called on expanded Data.");
1349      }
1350      if (m_op==IDENTITY)
1351      {
1352        throw DataException("Programmer error - resolveNodeNP1OUT_2P should not be called on identity nodes.");
1353    }    }
     // since we can't write the result over the input, we need a result offset further along  
1354    size_t subroffset;    size_t subroffset;
1355    const ValueType* vleft=m_left->resolveSample(v,offset+m_left->m_samplesize,sampleNo,subroffset);    size_t offset;
1356  LAZYDEBUG(cerr << "srcsamplesize=" << offset+m_left->m_samplesize << " beg=" << subroffset << endl;)    const ValueType* leftres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1357  LAZYDEBUG(cerr << "Offset for 5800=" << getPointOffset(5800/getNumDPPSample(),5800%getNumDPPSample()) << endl;)    roffset=m_samplesize*tid;
1358    roffset=offset;    offset=roffset;
1359    size_t loop=0;    size_t loop=0;
1360    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;    size_t numsteps=(m_readytype=='E')?getNumDPPSample():1;
1361    size_t outstep=getNoValues();    size_t outstep=getNoValues();
1362    size_t instep=m_left->getNoValues();    size_t instep=m_left->getNoValues();
 LAZYDEBUG(cerr << "instep=" << instep << " outstep=" << outstep<< " numsteps=" << numsteps << endl;)  
1363    switch (m_op)    switch (m_op)
1364    {    {
1365      case SWAP:      case SWAP:
1366      for (loop=0;loop<numsteps;++loop)      for (loop=0;loop<numsteps;++loop)
1367      {      {
1368              DataMaths::swapaxes(*vleft,m_left->getShape(),subroffset, v,getShape(),offset,m_axis_offset, m_transpose);              DataMaths::swapaxes(*leftres,m_left->getShape(),subroffset, m_samples, getShape(),offset, m_axis_offset, m_transpose);
1369          subroffset+=instep;          subroffset+=instep;
1370          offset+=outstep;          offset+=outstep;
1371      }      }
1372      break;      break;
1373      default:      default:
1374      throw DataException("Programmer error - resolveNP1OUT2P can not resolve operator "+opToString(m_op)+".");      throw DataException("Programmer error - resolveNodeNP1OUT2P can not resolve operator "+opToString(m_op)+".");
1375    }    }
1376    return &v;    return &m_samples;
1377  }  }
1378    
1379    const DataTypes::ValueType*
1380    DataLazy::resolveNodeCondEval(int tid, int sampleNo, size_t& roffset)
1381    {
1382      if (m_readytype!='E')
1383      {
1384        throw DataException("Programmer error - resolveNodeCondEval should only be called on expanded Data.");
1385      }
1386      if (m_op!=CONDEVAL)
1387      {
1388        throw DataException("Programmer error - resolveNodeCondEval should only be called on CONDEVAL nodes.");
1389      }
1390      size_t subroffset;
1391    
1392      const ValueType* maskres=m_mask->resolveNodeSample(tid, sampleNo, subroffset);
1393      const ValueType* srcres=0;
1394      if ((*maskres)[subroffset]>0)
1395      {
1396        srcres=m_left->resolveNodeSample(tid, sampleNo, subroffset);
1397      }
1398      else
1399      {
1400        srcres=m_right->resolveNodeSample(tid, sampleNo, subroffset);
1401      }
1402    
1403  #define PROC_OP(TYPE,X)                               \    // Now we need to copy the result
1404      for (int j=0;j<onumsteps;++j)\  
1405      {\    roffset=m_samplesize*tid;
1406        for (int i=0;i<numsteps;++i,resultp+=resultStep) \    for (int i=0;i<m_samplesize;++i)
1407        { \    {
1408  LAZYDEBUG(cout << "[left,right]=[" << lroffset << "," << rroffset << "]" << endl;)\      m_samples[roffset+i]=(*srcres)[subroffset+i];  
1409  LAZYDEBUG(cout << "{left,right}={" << (*left)[lroffset] << "," << (*right)[rroffset] << "}\n";)\    }
1410           tensor_binary_operation< TYPE >(chunksize, &((*left)[lroffset]), &((*right)[rroffset]), resultp, X); \  
1411  LAZYDEBUG(cout << " result=      " << resultp[0] << endl;) \    return &m_samples;
1412           lroffset+=leftstep; \  }
          rroffset+=rightstep; \  
       }\  
       lroffset+=oleftstep;\  
       rroffset+=orightstep;\  
     }  
1413    
 /*  
   \brief Compute the value of the expression (binary operation) for the given sample.  
   \return Vector which stores the value of the subexpression for the given sample.  
   \param v A vector to store intermediate results.  
   \param offset Index in v to begin storing results.  
   \param sampleNo Sample number to evaluate.  
   \param roffset (output parameter) the offset in the return vector where the result begins.  
   
   The return value will be an existing vector so do not deallocate it.  
   If the result is stored in v it should be stored at the offset given.  
   Everything from offset to the end of v should be considered available for this method to use.  
 */  
1414  // This method assumes that any subexpressions which evaluate to Constant or Tagged Data  // This method assumes that any subexpressions which evaluate to Constant or Tagged Data
1415  // have already been collapsed to IDENTITY. So we must have at least one expanded child.  // have already been collapsed to IDENTITY. So we must have at least one expanded child.
1416  // If both children are expanded, then we can process them in a single operation (we treat  // If both children are expanded, then we can process them in a single operation (we treat
# Line 1274  LAZYDEBUG(cout << " result=      " << re Line 1420  LAZYDEBUG(cout << " result=      " << re
1420  // There is an additional complication when scalar operations are considered.  // There is an additional complication when scalar operations are considered.
1421  // For example, 2+Vector.  // For example, 2+Vector.
1422  // In this case each double within the point is treated individually  // In this case each double within the point is treated individually
1423  DataTypes::ValueType*  const DataTypes::ValueType*
1424  DataLazy::resolveBinary(ValueType& v,  size_t offset, int sampleNo, size_t& roffset) const  DataLazy::resolveNodeBinary(int tid, int sampleNo, size_t& roffset)
1425  {  {
1426  LAZYDEBUG(cout << "Resolve binary: " << toString() << endl;)  LAZYDEBUG(cout << "Resolve binary: " << toString() << endl;)
1427    
# Line 1387  LAZYDEBUG(cout << "Resolve binary: " << Line 1533  LAZYDEBUG(cout << "Resolve binary: " <<
1533    
1534    int resultStep=max(leftstep,rightstep);   // only one (at most) should be !=0    int resultStep=max(leftstep,rightstep);   // only one (at most) should be !=0
1535      // Get the values of sub-expressions      // Get the values of sub-expressions
1536    const ValueType* left=m_left->resolveSample(v,offset+getMaxSampleSize(),sampleNo,lroffset);   // see note on    const ValueType* left=m_left->resolveNodeSample(tid,sampleNo,lroffset);  
1537      // calcBufss for why we can't put offset as the 2nd param above    const ValueType* right=m_right->resolveNodeSample(tid,sampleNo,rroffset);
   const ValueType* right=m_right->resolveSample(v,offset+2*getMaxSampleSize(),sampleNo,rroffset); // Note  
     // the right child starts further along.  
1538  LAZYDEBUG(cout << "Post sub calls in " << toString() << endl;)  LAZYDEBUG(cout << "Post sub calls in " << toString() << endl;)
1539  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;)
1540  LAZYDEBUG(cout << "chunksize=" << chunksize << endl << "leftstep=" << leftstep << " rightstep=" << rightstep;)  LAZYDEBUG(cout << "chunksize=" << chunksize << endl << "leftstep=" << leftstep << " rightstep=" << rightstep;)
# Line 1399  LAZYDEBUG(cout << "onumsteps=" << onumst Line 1543  LAZYDEBUG(cout << "onumsteps=" << onumst
1543  LAZYDEBUG(cout << " DPPS=" << m_left->getNumDPPSample() << "," <<m_right->getNumDPPSample() << endl;)  LAZYDEBUG(cout << " DPPS=" << m_left->getNumDPPSample() << "," <<m_right->getNumDPPSample() << endl;)
1544  LAZYDEBUG(cout << "" << LS << RS << LN << RN << LES << RES <<LEN << REN <<   endl;)  LAZYDEBUG(cout << "" << LS << RS << LN << RN << LES << RES <<LEN << REN <<   endl;)
1545    
1546    LAZYDEBUG(cout << "Left res["<< lroffset<< "]=" << (*left)[lroffset] << endl;)
1547    LAZYDEBUG(cout << "Right res["<< rroffset<< "]=" << (*right)[rroffset] << endl;)
1548    
1549    double* resultp=&(v[offset]);     // results are stored at the vector offset we recieved  
1550      roffset=m_samplesize*tid;
1551      double* resultp=&(m_samples[roffset]);        // results are stored at the vector offset we recieved
1552    switch(m_op)    switch(m_op)
1553    {    {
1554      case ADD:      case ADD:
# Line 1421  LAZYDEBUG(cout << "" << LS << RS << LN < Line 1569  LAZYDEBUG(cout << "" << LS << RS << LN <
1569      default:      default:
1570      throw DataException("Programmer error - resolveBinary can not resolve operator "+opToString(m_op)+".");      throw DataException("Programmer error - resolveBinary can not resolve operator "+opToString(m_op)+".");
1571    }    }
1572    roffset=offset;  LAZYDEBUG(cout << "Result res[" << roffset<< "]" << m_samples[roffset] << endl;)
1573    return &v;    return &m_samples;
1574  }  }
1575    
1576    
   
 /*  
   \brief Compute the value of the expression (tensor product) for the given sample.  
   \return Vector which stores the value of the subexpression for the given sample.  
   \param v A vector to store intermediate results.  
   \param offset Index in v to begin storing results.  
   \param sampleNo Sample number to evaluate.  
   \param roffset (output parameter) the offset in the return vector where the result begins.  
   
   The return value will be an existing vector so do not deallocate it.  
   If the result is stored in v it should be stored at the offset given.  
   Everything from offset to the end of v should be considered available for this method to use.  
 */  
1577  // This method assumes that any subexpressions which evaluate to Constant or Tagged Data  // This method assumes that any subexpressions which evaluate to Constant or Tagged Data
1578  // have already been collapsed to IDENTITY. So we must have at least one expanded child.  // have already been collapsed to IDENTITY. So we must have at least one expanded child.
1579  // unlike the other resolve helpers, we must treat these datapoints separately.  // unlike the other resolve helpers, we must treat these datapoints separately.
1580  DataTypes::ValueType*  const DataTypes::ValueType*
1581  DataLazy::resolveTProd(ValueType& v,  size_t offset, int sampleNo, size_t& roffset) const  DataLazy::resolveNodeTProd(int tid, int sampleNo, size_t& roffset)
1582  {  {
1583  LAZYDEBUG(cout << "Resolve TensorProduct: " << toString()  << " to offset " << offset<< endl;)  LAZYDEBUG(cout << "Resolve TensorProduct: " << toString() << endl;)
1584    
1585    size_t lroffset=0, rroffset=0;    // offsets in the left and right result vectors    size_t lroffset=0, rroffset=0;    // offsets in the left and right result vectors
1586      // first work out which of the children are expanded      // first work out which of the children are expanded
1587    bool leftExp=(m_left->m_readytype=='E');    bool leftExp=(m_left->m_readytype=='E');
1588    bool rightExp=(m_right->m_readytype=='E');    bool rightExp=(m_right->m_readytype=='E');
1589    int steps=getNumDPPSample();    int steps=getNumDPPSample();
 /*  int leftStep=((leftExp && !rightExp)? m_right->getNoValues() : 0);  
   int rightStep=((rightExp && !leftExp)? m_left->getNoValues() : 0);*/  
1590    int leftStep=(leftExp? m_left->getNoValues() : 0);        // do not have scalars as input to this method    int leftStep=(leftExp? m_left->getNoValues() : 0);        // do not have scalars as input to this method
1591    int rightStep=(rightExp?m_right->getNoValues() : 0);    int rightStep=(rightExp?m_right->getNoValues() : 0);
1592    
1593    int resultStep=getNoValues();    int resultStep=getNoValues();
1594      // Get the values of sub-expressions (leave a gap of one sample for the result).    roffset=m_samplesize*tid;
1595    int gap=offset+m_samplesize;      size_t offset=roffset;
   
 LAZYDEBUG(cout << "Query left with offset=" << gap << endl;)  
   
   const ValueType* left=m_left->resolveSample(v,gap,sampleNo,lroffset);  
   gap+=m_left->getMaxSampleSize();  
   
1596    
1597  LAZYDEBUG(cout << "Query right with offset=" << gap << endl;)    const ValueType* left=m_left->resolveNodeSample(tid, sampleNo, lroffset);
1598    
1599      const ValueType* right=m_right->resolveNodeSample(tid, sampleNo, rroffset);
   const ValueType* right=m_right->resolveSample(v,gap,sampleNo,rroffset);  
1600    
1601  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;
1602  cout << getNoValues() << endl;)  cout << getNoValues() << endl;)
 LAZYDEBUG(cerr << "Result of left=";)  
 LAZYDEBUG(cerr << "[" << lroffset << " .. " << lroffset+m_left->getNoValues() << "]" << endl;  
1603    
1604  for (int i=lroffset, limit=lroffset+(leftExp?m_left->getNoValues()*m_left->getNumDPPSample():m_left->getNoValues());i<limit;++i)  
 {  
 cout << "[" << setw(2) << i-lroffset << "] " << setw(10) << (*left)[i] << " ";  
 if (i%4==0) cout << endl;  
 })  
 LAZYDEBUG(cerr << "\nResult of right=" << endl;)  
 LAZYDEBUG(  
 for (int i=rroffset, limit=rroffset+(rightExp?m_right->getNoValues()*m_right->getNumDPPSample():m_right->getNoValues());i<limit;++i)  
 {  
 cerr << "[" <<  setw(2)<< i-rroffset << "] " << setw(10) << (*right)[i] << " ";  
 if (i%4==0) cout << endl;  
 }  
 cerr << endl;  
 )  
1605  LAZYDEBUG(cerr << "Post sub calls: " << toString() << endl;)  LAZYDEBUG(cerr << "Post sub calls: " << toString() << endl;)
1606  LAZYDEBUG(cout << "LeftExp=" << leftExp << " rightExp=" << rightExp << endl;)  LAZYDEBUG(cout << "LeftExp=" << leftExp << " rightExp=" << rightExp << endl;)
1607  LAZYDEBUG(cout << "LeftR=" << m_left->getRank() << " rightExp=" << m_right->getRank() << endl;)  LAZYDEBUG(cout << "LeftR=" << m_left->getRank() << " rightExp=" << m_right->getRank() << endl;)
# Line 1499  LAZYDEBUG(cout << "m_samplesize=" << m_s Line 1610  LAZYDEBUG(cout << "m_samplesize=" << m_s
1610  LAZYDEBUG(cout << "outputshape=" << DataTypes::shapeToString(getShape()) << endl;)  LAZYDEBUG(cout << "outputshape=" << DataTypes::shapeToString(getShape()) << endl;)
1611  LAZYDEBUG(cout << "DPPS=" << m_right->getNumDPPSample() <<"."<<endl;)  LAZYDEBUG(cout << "DPPS=" << m_right->getNumDPPSample() <<"."<<endl;)
1612    
1613    double* resultp=&(v[offset]);     // results are stored at the vector offset we recieved    double* resultp=&(m_samples[offset]);     // results are stored at the vector offset we recieved
1614    switch(m_op)    switch(m_op)
1615    {    {
1616      case PROD:      case PROD:
1617      for (int i=0;i<steps;++i,resultp+=resultStep)      for (int i=0;i<steps;++i,resultp+=resultStep)
1618      {      {
   
 LAZYDEBUG(cout << "lroffset=" << lroffset << "rroffset=" << rroffset << endl;)  
 LAZYDEBUG(cout << "l*=" << left << " r*=" << right << endl;)  
 LAZYDEBUG(cout << "m_SL=" << m_SL << " m_SM=" << m_SM << " m_SR=" << m_SR << endl;)  
   
1619            const double *ptr_0 = &((*left)[lroffset]);            const double *ptr_0 = &((*left)[lroffset]);
1620            const double *ptr_1 = &((*right)[rroffset]);            const double *ptr_1 = &((*right)[rroffset]);
1621    
# Line 1518  LAZYDEBUG(cout << DataTypes::pointToStri Line 1624  LAZYDEBUG(cout << DataTypes::pointToStri
1624    
1625            matrix_matrix_product(m_SL, m_SM, m_SR, ptr_0, ptr_1, resultp, m_transpose);            matrix_matrix_product(m_SL, m_SM, m_SR, ptr_0, ptr_1, resultp, m_transpose);
1626    
 LAZYDEBUG(cout << "Results--\n";  
 {  
   DataVector dv(getNoValues());  
 for (int z=0;z<getNoValues();++z)  
 {  
   cout << "[" << setw(2) << z<< "] " << setw(10) << resultp[z] << " ";  
   if (z%4==0) cout << endl;  
   dv[z]=resultp[z];  
 }  
 cout << endl << DataTypes::pointToString(dv,getShape(),0,"RESLT");  
 cout << "\nWritten to: " << resultp << " resultStep=" << resultStep << endl;  
 }  
 )  
1627        lroffset+=leftStep;        lroffset+=leftStep;
1628        rroffset+=rightStep;        rroffset+=rightStep;
1629      }      }
# Line 1539  cout << "\nWritten to: " << resultp << " Line 1632  cout << "\nWritten to: " << resultp << "
1632      throw DataException("Programmer error - resolveTProduct can not resolve operator "+opToString(m_op)+".");      throw DataException("Programmer error - resolveTProduct can not resolve operator "+opToString(m_op)+".");
1633    }    }
1634    roffset=offset;    roffset=offset;
1635    return &v;    return &m_samples;
1636  }  }
1637    
1638    
   
 /*  
   \brief Compute the value of the expression for the given sample.  
   \return Vector which stores the value of the subexpression for the given sample.  
   \param v A vector to store intermediate results.  
   \param offset Index in v to begin storing results.  
   \param sampleNo Sample number to evaluate.  
   \param roffset (output parameter) the offset in the return vector where the result begins.  
   
   The return value will be an existing vector so do not deallocate it.  
 */  
 // the vector and the offset are a place where the method could write its data if it wishes  
 // it is not obligated to do so. For example, if it has its own storage already, it can use that.  
 // Hence the return value to indicate where the data is actually stored.  
 // Regardless, the storage should be assumed to be used, even if it isn't.  
   
 // the roffset is the offset within the returned vector where the data begins  
1639  const DataTypes::ValueType*  const DataTypes::ValueType*
1640  DataLazy::resolveSample(ValueType& v, size_t offset, int sampleNo, size_t& roffset)  DataLazy::resolveSample(int sampleNo, size_t& roffset)
 {  
 LAZYDEBUG(cout << "Resolve sample " << toString() << endl;)  
     // collapse so we have a 'E' node or an IDENTITY for some other type  
   if (m_readytype!='E' && m_op!=IDENTITY)  
   {  
     collapse();  
   }  
   if (m_op==IDENTITY)    
   {  
     const ValueType& vec=m_id->getVectorRO();  
     if (m_readytype=='C')  
     {  
     roffset=0;  
 LAZYDEBUG(cout << "Finish  sample " << toString() << endl;)  
     return &(vec);  
     }  
     roffset=m_id->getPointOffset(sampleNo, 0);  
 LAZYDEBUG(cout << "Finish  sample " << toString() << endl;)  
     return &(vec);  
   }  
   if (m_readytype!='E')  
   {  
     throw DataException("Programmer Error - Collapse did not produce an expanded node.");  
   }  
   switch (getOpgroup(m_op))  
   {  
   case G_UNARY:  
   case G_UNARY_P: return resolveUnary(v, offset,sampleNo,roffset);  
   case G_BINARY: return resolveBinary(v, offset,sampleNo,roffset);  
   case G_NP1OUT: return resolveNP1OUT(v, offset, sampleNo,roffset);  
   case G_NP1OUT_P: return resolveNP1OUT_P(v, offset, sampleNo,roffset);  
   case G_TENSORPROD: return resolveTProd(v,offset, sampleNo,roffset);  
   case G_NP1OUT_2P: return resolveNP1OUT_2P(v, offset, sampleNo, roffset);  
   default:  
     throw DataException("Programmer Error - resolveSample does not know how to process "+opToString(m_op)+".");  
   }  
   
 }  
   
 const DataTypes::ValueType*  
 DataLazy::resolveSample(BufferGroup& bg, int sampleNo, size_t& roffset)  
1641  {  {
1642  #ifdef _OPENMP  #ifdef _OPENMP
1643      int tid=omp_get_thread_num();      int tid=omp_get_thread_num();
1644  #else  #else
1645      int tid=0;      int tid=0;
1646  #endif  #endif
1647      return resolveSample(bg.getBuffer(tid),bg.getOffset(tid),sampleNo,roffset);  
1648    #ifdef LAZY_STACK_PROF
1649        stackstart[tid]=&tid;
1650        stackend[tid]=&tid;
1651        const DataTypes::ValueType* r=resolveNodeSample(tid, sampleNo, roffset);
1652        size_t d=(size_t)stackstart[tid]-(size_t)stackend[tid];
1653        #pragma omp critical
1654        if (d>maxstackuse)
1655        {
1656    cout << "Max resolve Stack use " << d << endl;
1657            maxstackuse=d;
1658        }
1659        return r;
1660    #else
1661        return resolveNodeSample(tid, sampleNo, roffset);
1662    #endif
1663  }  }
1664    
1665    
1666  // This needs to do the work of the idenity constructor  // This needs to do the work of the identity constructor
1667  void  void
1668  DataLazy::resolveToIdentity()  DataLazy::resolveToIdentity()
1669  {  {
1670     if (m_op==IDENTITY)     if (m_op==IDENTITY)
1671      return;      return;
1672     DataReady_ptr p=resolve();     DataReady_ptr p=resolveNodeWorker();
1673     makeIdentity(p);     makeIdentity(p);
1674  }  }
1675    
# Line 1635  void DataLazy::makeIdentity(const DataRe Line 1685  void DataLazy::makeIdentity(const DataRe
1685     else if(p->isExpanded()) {m_readytype='E';}     else if(p->isExpanded()) {m_readytype='E';}
1686     else if (p->isTagged()) {m_readytype='T';}     else if (p->isTagged()) {m_readytype='T';}
1687     else {throw DataException("Unknown DataReady instance in convertToIdentity constructor.");}     else {throw DataException("Unknown DataReady instance in convertToIdentity constructor.");}
    m_buffsRequired=1;  
1688     m_samplesize=p->getNumDPPSample()*p->getNoValues();     m_samplesize=p->getNumDPPSample()*p->getNoValues();
    m_maxsamplesize=m_samplesize;  
1689     m_left.reset();     m_left.reset();
1690     m_right.reset();     m_right.reset();
1691  }  }
1692    
1693  // To simplify the memory management, all threads operate on one large vector, rather than one each.  
 // Each sample is evaluated independently and copied into the result DataExpanded.  
1694  DataReady_ptr  DataReady_ptr
1695  DataLazy::resolve()  DataLazy::resolve()
1696  {  {
1697        resolveToIdentity();
1698        return m_id;
1699    }
1700    
1701    
1702    /* This is really a static method but I think that caused problems in windows */
1703    void
1704    DataLazy::resolveGroupWorker(std::vector<DataLazy*>& dats)
1705    {
1706      if (dats.empty())
1707      {
1708        return;
1709      }
1710      vector<DataLazy*> work;
1711      FunctionSpace fs=dats[0]->getFunctionSpace();
1712      bool match=true;
1713      for (int i=dats.size()-1;i>=0;--i)
1714      {
1715        if (dats[i]->m_readytype!='E')
1716        {
1717            dats[i]->collapse();
1718        }
1719        if (dats[i]->m_op!=IDENTITY)
1720        {
1721            work.push_back(dats[i]);
1722            if (fs!=dats[i]->getFunctionSpace())
1723            {
1724                match=false;
1725            }
1726        }
1727      }
1728      if (work.empty())
1729      {
1730        return;     // no work to do
1731      }
1732      if (match)    // all functionspaces match.  Yes I realise this is overly strict
1733      {     // it is possible that dats[0] is one of the objects which we discarded and
1734            // all the other functionspaces match.
1735        vector<DataExpanded*> dep;
1736        vector<ValueType*> vecs;
1737        for (int i=0;i<work.size();++i)
1738        {
1739            dep.push_back(new DataExpanded(fs,work[i]->getShape(), ValueType(work[i]->getNoValues())));
1740            vecs.push_back(&(dep[i]->getVectorRW()));
1741        }
1742        int totalsamples=work[0]->getNumSamples();
1743        const ValueType* res=0; // Storage for answer
1744        int sample;
1745        #pragma omp parallel private(sample, res)
1746        {
1747            size_t roffset=0;
1748            #pragma omp for schedule(static)
1749            for (sample=0;sample<totalsamples;++sample)
1750            {
1751            roffset=0;
1752            int j;
1753            for (j=work.size()-1;j>=0;--j)
1754            {
1755    #ifdef _OPENMP
1756                    res=work[j]->resolveNodeSample(omp_get_thread_num(),sample,roffset);
1757    #else
1758                    res=work[j]->resolveNodeSample(0,sample,roffset);
1759    #endif
1760                    DataVector::size_type outoffset=dep[j]->getPointOffset(sample,0);
1761                    memcpy(&((*vecs[j])[outoffset]),&((*res)[roffset]),work[j]->m_samplesize*sizeof(DataVector::ElementType));
1762            }
1763            }
1764        }
1765        // Now we need to load the new results as identity ops into the lazy nodes
1766        for (int i=work.size()-1;i>=0;--i)
1767        {
1768            work[i]->makeIdentity(boost::dynamic_pointer_cast<DataReady>(dep[i]->getPtr()));
1769        }
1770      }
1771      else  // functionspaces do not match
1772      {
1773        for (int i=0;i<work.size();++i)
1774        {
1775            work[i]->resolveToIdentity();
1776        }
1777      }
1778    }
1779    
1780    
1781  LAZYDEBUG(cout << "Sample size=" << m_samplesize << endl;)  
1782  LAZYDEBUG(cout << "Buffers=" << m_buffsRequired << endl;)  // This version of resolve uses storage in each node to hold results
1783    DataReady_ptr
1784    DataLazy::resolveNodeWorker()
1785    {
1786    if (m_readytype!='E')     // if the whole sub-expression is Constant or Tagged, then evaluate it normally    if (m_readytype!='E')     // if the whole sub-expression is Constant or Tagged, then evaluate it normally
1787    {    {
1788      collapse();      collapse();
# Line 1659  LAZYDEBUG(cout << "Buffers=" << m_buffsR Line 1792  LAZYDEBUG(cout << "Buffers=" << m_buffsR
1792      return m_id;      return m_id;
1793    }    }
1794      // 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'
   size_t threadbuffersize=m_maxsamplesize*(max(1,m_buffsRequired)); // Each thread needs to have enough  
     // storage to evaluate its expression  
   int numthreads=1;  
 #ifdef _OPENMP  
   numthreads=omp_get_max_threads();  
 #endif  
   ValueType v(numthreads*threadbuffersize);  
 LAZYDEBUG(cout << "Buffer created with size=" << v.size() << endl;)  
1795    DataExpanded* result=new DataExpanded(getFunctionSpace(),getShape(),  ValueType(getNoValues()));    DataExpanded* result=new DataExpanded(getFunctionSpace(),getShape(),  ValueType(getNoValues()));
1796    ValueType& resvec=result->getVectorRW();    ValueType& resvec=result->getVectorRW();
1797    DataReady_ptr resptr=DataReady_ptr(result);    DataReady_ptr resptr=DataReady_ptr(result);
1798    
1799    int sample;    int sample;
   size_t outoffset;     // offset in the output data  
1800    int totalsamples=getNumSamples();    int totalsamples=getNumSamples();
1801    const ValueType* res=0;   // Vector storing the answer    const ValueType* res=0;   // Storage for answer
   size_t resoffset=0;       // where in the vector to find the answer  
1802  LAZYDEBUG(cout << "Total number of samples=" <<totalsamples << endl;)  LAZYDEBUG(cout << "Total number of samples=" <<totalsamples << endl;)
1803    #pragma omp parallel for private(sample,resoffset,outoffset,res) schedule(static)    #pragma omp parallel private(sample,res)
   for (sample=0;sample<totalsamples;++sample)  
1804    {    {
1805        if (sample==0)  {ENABLEDEBUG}      size_t roffset=0;
1806    #ifdef LAZY_STACK_PROF
1807  //       if (sample==5800/getNumDPPSample())  {ENABLEDEBUG}      stackstart[omp_get_thread_num()]=&roffset;
1808  LAZYDEBUG(cout << "################################# " << sample << endl;)      stackend[omp_get_thread_num()]=&roffset;
1809    #endif
1810        #pragma omp for schedule(static)
1811        for (sample=0;sample<totalsamples;++sample)
1812        {
1813            roffset=0;
1814  #ifdef _OPENMP  #ifdef _OPENMP
1815      res=resolveSample(v,threadbuffersize*omp_get_thread_num(),sample,resoffset);              res=resolveNodeSample(omp_get_thread_num(),sample,roffset);
1816  #else  #else
1817      res=resolveSample(v,0,sample,resoffset);   // res would normally be v, but not if its a single IDENTITY op.              res=resolveNodeSample(0,sample,roffset);
1818  #endif  #endif
1819  LAZYDEBUG(cerr << "-------------------------------- " << endl;)  LAZYDEBUG(cout << "Sample #" << sample << endl;)
1820  LAZYDEBUG(cerr<< "Copying sample#" << sample << endl;)  LAZYDEBUG(cout << "Final res[" << roffset<< "]=" << (*res)[roffset] << (*res)[roffset]<< endl; )
1821      outoffset=result->getPointOffset(sample,0);              DataVector::size_type outoffset=result->getPointOffset(sample,0);
1822  LAZYDEBUG(cerr << "offset=" << outoffset << " from offset=" << resoffset << " " << m_samplesize << " doubles" << endl;)              memcpy(&(resvec[outoffset]),&((*res)[roffset]),m_samplesize*sizeof(DataVector::ElementType));
1823      for (unsigned int i=0;i<m_samplesize;++i,++outoffset,++resoffset)   // copy values into the output vector      }
1824      {    }
1825  LAZYDEBUG(cerr << "outoffset=" << outoffset << " resoffset=" << resoffset << " " << (*res)[resoffset]<< endl;)  #ifdef LAZY_STACK_PROF
1826      resvec[outoffset]=(*res)[resoffset];    for (int i=0;i<getNumberOfThreads();++i)
1827      }    {
1828  LAZYDEBUG(cerr << DataTypes::pointToString(resvec,getShape(),outoffset-m_samplesize+DataTypes::noValues(getShape()),"Final result:") << endl;)      size_t r=((size_t)stackstart[i] - (size_t)stackend[i]);
1829  LAZYDEBUG(cerr << "*********************************" << endl;)  //  cout << i << " " << stackstart[i] << " .. " << stackend[i] << " = " <<  r << endl;
1830      DISABLEDEBUG      if (r>maxstackuse)
1831        {
1832            maxstackuse=r;
1833        }
1834    }    }
1835      cout << "Max resolve Stack use=" << maxstackuse << endl;
1836    #endif
1837    return resptr;    return resptr;
1838  }  }
1839    
# Line 1708  std::string Line 1841  std::string
1841  DataLazy::toString() const  DataLazy::toString() const
1842  {  {
1843    ostringstream oss;    ostringstream oss;
1844    oss << "Lazy Data:";    oss << "Lazy Data: [depth=" << m_height<< "] ";
1845    intoString(oss);    switch (escriptParams.getLAZY_STR_FMT())
1846      {
1847      case 1:   // tree format
1848        oss << endl;
1849        intoTreeString(oss,"");
1850        break;
1851      case 2:   // just the depth
1852        break;
1853      default:
1854        intoString(oss);
1855        break;
1856      }
1857    return oss.str();    return oss.str();
1858  }  }
1859    
# Line 1750  DataLazy::intoString(ostringstream& oss) Line 1894  DataLazy::intoString(ostringstream& oss)
1894    case G_UNARY_P:    case G_UNARY_P:
1895    case G_NP1OUT:    case G_NP1OUT:
1896    case G_NP1OUT_P:    case G_NP1OUT_P:
1897      case G_REDUCTION:
1898      oss << opToString(m_op) << '(';      oss << opToString(m_op) << '(';
1899      m_left->intoString(oss);      m_left->intoString(oss);
1900      oss << ')';      oss << ')';
# Line 1767  DataLazy::intoString(ostringstream& oss) Line 1912  DataLazy::intoString(ostringstream& oss)
1912      oss << ", " << m_axis_offset << ", " << m_transpose;      oss << ", " << m_axis_offset << ", " << m_transpose;
1913      oss << ')';      oss << ')';
1914      break;      break;
1915      case G_CONDEVAL:
1916        oss << opToString(m_op)<< '(' ;
1917        m_mask->intoString(oss);
1918        oss << " ? ";
1919        m_left->intoString(oss);
1920        oss << " : ";
1921        m_right->intoString(oss);
1922        oss << ')';
1923        break;
1924    default:    default:
1925      oss << "UNKNOWN";      oss << "UNKNOWN";
1926    }    }
1927  }  }
1928    
1929    
1930    void
1931    DataLazy::intoTreeString(ostringstream& oss, string indent) const
1932    {
1933      oss << '[' << m_rank << ':' << setw(3) << m_samplesize << "] " << indent;
1934      switch (getOpgroup(m_op))
1935      {
1936      case G_IDENTITY:
1937        if (m_id->isExpanded())
1938        {
1939           oss << "E";
1940        }
1941        else if (m_id->isTagged())
1942        {
1943          oss << "T";
1944        }
1945        else if (m_id->isConstant())
1946        {
1947          oss << "C";
1948        }
1949        else
1950        {
1951          oss << "?";
1952        }
1953        oss << '@' << m_id.get() << endl;
1954        break;
1955      case G_BINARY:
1956        oss << opToString(m_op) << endl;
1957        indent+='.';
1958        m_left->intoTreeString(oss, indent);
1959        m_right->intoTreeString(oss, indent);
1960        break;
1961      case G_UNARY:
1962      case G_UNARY_P:
1963      case G_NP1OUT:
1964      case G_NP1OUT_P:
1965      case G_REDUCTION:
1966        oss << opToString(m_op) << endl;
1967        indent+='.';
1968        m_left->intoTreeString(oss, indent);
1969        break;
1970      case G_TENSORPROD:
1971        oss << opToString(m_op) << endl;
1972        indent+='.';
1973        m_left->intoTreeString(oss, indent);
1974        m_right->intoTreeString(oss, indent);
1975        break;
1976      case G_NP1OUT_2P:
1977        oss << opToString(m_op) << ", " << m_axis_offset << ", " << m_transpose<< endl;
1978        indent+='.';
1979        m_left->intoTreeString(oss, indent);
1980        break;
1981      default:
1982        oss << "UNKNOWN";
1983      }
1984    }
1985    
1986    
1987  DataAbstract*  DataAbstract*
1988  DataLazy::deepCopy()  DataLazy::deepCopy()
1989  {  {
1990    switch (getOpgroup(m_op))    switch (getOpgroup(m_op))
1991    {    {
1992    case G_IDENTITY:  return new DataLazy(m_id->deepCopy()->getPtr());    case G_IDENTITY:  return new DataLazy(m_id->deepCopy()->getPtr());
1993    case G_UNARY: return new DataLazy(m_left->deepCopy()->getPtr(),m_op);    case G_UNARY:
1994      case G_REDUCTION:      return new DataLazy(m_left->deepCopy()->getPtr(),m_op);
1995      case G_UNARY_P:   return new DataLazy(m_left->deepCopy()->getPtr(), m_op, m_tol);
1996    case G_BINARY:    return new DataLazy(m_left->deepCopy()->getPtr(),m_right->deepCopy()->getPtr(),m_op);    case G_BINARY:    return new DataLazy(m_left->deepCopy()->getPtr(),m_right->deepCopy()->getPtr(),m_op);
1997    case G_NP1OUT: return new DataLazy(m_left->deepCopy()->getPtr(), m_right->deepCopy()->getPtr(),m_op);    case G_NP1OUT: return new DataLazy(m_left->deepCopy()->getPtr(), m_right->deepCopy()->getPtr(),m_op);
1998    case G_TENSORPROD: return new DataLazy(m_left->deepCopy()->getPtr(), m_right->deepCopy()->getPtr(), m_op, m_axis_offset, m_transpose);    case G_TENSORPROD: return new DataLazy(m_left->deepCopy()->getPtr(), m_right->deepCopy()->getPtr(), m_op, m_axis_offset, m_transpose);
1999      case G_NP1OUT_P:   return new DataLazy(m_left->deepCopy()->getPtr(),m_op,  m_axis_offset);
2000      case G_NP1OUT_2P:  return new DataLazy(m_left->deepCopy()->getPtr(), m_op, m_axis_offset, m_transpose);
2001    default:    default:
2002      throw DataException("Programmer error - do not know how to deepcopy operator "+opToString(m_op)+".");      throw DataException("Programmer error - do not know how to deepcopy operator "+opToString(m_op)+".");
2003    }    }
2004  }  }
2005    
2006    
2007    
2008  // There is no single, natural interpretation of getLength on DataLazy.  // There is no single, natural interpretation of getLength on DataLazy.
2009  // Instances of DataReady can look at the size of their vectors.  // Instances of DataReady can look at the size of their vectors.
2010  // For lazy though, it could be the size the data would be if it were resolved;  // For lazy though, it could be the size the data would be if it were resolved;
# Line 1876  DataLazy::setToZero() Line 2093  DataLazy::setToZero()
2093  //   m_readytype='C';  //   m_readytype='C';
2094  //   m_buffsRequired=1;  //   m_buffsRequired=1;
2095    
2096      privdebug=privdebug;  // to stop the compiler complaining about unused privdebug
2097    throw DataException("Programmer error - setToZero not supported for DataLazy (DataLazy objects should be read only).");    throw DataException("Programmer error - setToZero not supported for DataLazy (DataLazy objects should be read only).");
   
2098  }  }
2099    
2100  bool  bool

Legend:
Removed from v.2496  
changed lines
  Added in v.4154

  ViewVC Help
Powered by ViewVC 1.1.26