ESYS13  Revision_
DataAlgorithm.h
Go to the documentation of this file.
00001 
00002 /*******************************************************
00003 *
00004 * Copyright (c) 2003-2012 by University of Queensland
00005 * Earth Systems Science Computational Center (ESSCC)
00006 * http://www.uq.edu.au/esscc
00007 *
00008 * Primary Business: Queensland, Australia
00009 * Licensed under the Open Software License version 3.0
00010 * http://www.opensource.org/licenses/osl-3.0.php
00011 *
00012 *******************************************************/
00013 
00014 
00015 #if !defined escript_DataAlgorithm_20040714_H
00016 #define escript_DataAlgorithm_20040714_H
00017 #include "system_dep.h"
00018 
00019 #include "DataExpanded.h"
00020 #include "DataTagged.h"
00021 #include "DataConstant.h"
00022 
00023 #include "DataMaths.h"
00024 
00025 #include <iostream>
00026 #include <algorithm>
00027 #include <list>
00028 
00029 namespace escript {
00030 
00041 template <class BinaryFunction>
00042 class DataAlgorithmAdapter {
00043   public:
00044     DataAlgorithmAdapter(double initialValue):
00045       m_initialValue(initialValue),
00046       m_currentValue(initialValue)
00047     {
00048     }
00049     DataAlgorithmAdapter(const DataAlgorithmAdapter& other):
00050       m_initialValue(other.m_initialValue),
00051       m_currentValue(other.m_initialValue),
00052       operation(other.operation)
00053     {
00054     }
00055     inline void operator()(double value)
00056     {
00057       m_currentValue=operation(m_currentValue,value);
00058       return;
00059     }
00060     inline void resetResult()
00061     {
00062       m_currentValue=m_initialValue;
00063     }
00064     inline double getResult() const
00065     {
00066       return m_currentValue;
00067     }
00068   private:
00069     //
00070     // the initial operation value
00071     double m_initialValue;
00072     //
00073     // the current operation value
00074     double m_currentValue;
00075     //
00076     // The operation to perform
00077     BinaryFunction operation;
00078 };
00079 
00084 struct FMax : public std::binary_function<double,double,double>
00085 {
00086   inline double operator()(double x, double y) const
00087   {
00088     return std::max(x,y);
00089   }
00090 };
00091 
00096 struct FMin : public std::binary_function<double,double,double>
00097 {
00098   inline double operator()(double x, double y) const
00099   {
00100     return std::min(x,y);
00101   }
00102 };
00103 
00108 struct AbsMax : public std::binary_function<double,double,double>
00109 {
00110   inline double operator()(double x, double y) const
00111   {
00112     return std::max(fabs(x),fabs(y));
00113   }
00114 };
00115 
00120 struct AbsMin : public std::binary_function<double,double,double>
00121 {
00122   inline double operator()(double x, double y) const
00123   {
00124     return std::min(fabs(x),fabs(y));
00125   }
00126 };
00127 
00132 struct Length : public std::binary_function<double,double,double>
00133 {
00134   inline double operator()(double x, double y) const
00135   {
00136     return std::sqrt(std::pow(x,2)+std::pow(y,2));
00137   }
00138 };
00139 
00144 struct Trace : public std::binary_function<double,double,double>
00145 {
00146   inline double operator()(double x, double y) const
00147   {
00148     return x+y;
00149   }
00150 };
00151 
00155 struct AbsGT : public std::binary_function<double,double,double>
00156 {
00157   inline double operator()(double x, double y) const
00158   {
00159     return fabs(x)>y;
00160   }
00161 };
00162 
00166 struct AbsLTE : public std::binary_function<double,double,double>
00167 {
00168   inline double operator()(double x, double y) const
00169   {
00170     return fabs(x)<=y;
00171   }
00172 };
00173 
00174 
00180 template <class BinaryFunction>
00181 inline
00182 double
00183 algorithm(const DataExpanded& data,
00184           BinaryFunction operation,
00185       double initial_value)
00186 {
00187   int i,j;
00188   int numDPPSample=data.getNumDPPSample();
00189   int numSamples=data.getNumSamples();
00190   double global_current_value=initial_value;
00191   double local_current_value;
00192 //  DataArrayView dataView=data.getPointDataView();
00193   const DataTypes::ValueType& vec=data.getVectorRO();
00194   const DataTypes::ShapeType& shape=data.getShape();
00195   // calculate the reduction operation value for each data point
00196   // reducing the result for each data-point into the current_value variables
00197   #pragma omp parallel private(local_current_value)
00198   {
00199       local_current_value=initial_value;
00200       #pragma omp for private(i,j) schedule(static)
00201       for (i=0;i<numSamples;i++) {
00202         for (j=0;j<numDPPSample;j++) {
00203 /*          local_current_value=operation(local_current_value,dataView.reductionOp(data.getPointOffset(i,j),operation,initial_value));*/
00204           local_current_value=operation(local_current_value,DataMaths::reductionOp(vec,shape,data.getPointOffset(i,j),operation,initial_value));
00205 
00206         }
00207       }
00208       #pragma omp critical
00209       global_current_value=operation(global_current_value,local_current_value);
00210   }
00211   return global_current_value;
00212 }
00213 
00214 // It is important that the algorithm only be applied to tags which are actually in use.
00215 template <class BinaryFunction>
00216 inline
00217 double
00218 algorithm(DataTagged& data,
00219           BinaryFunction operation,
00220       double initial_value)
00221 {
00222   double current_value=initial_value;
00223 
00224   const DataTypes::ValueType& vec=data.getVectorRO();
00225   const DataTypes::ShapeType& shape=data.getShape();
00226   const DataTagged::DataMapType& lookup=data.getTagLookup();
00227   const std::list<int> used=data.getFunctionSpace().getListOfTagsSTL();
00228   for (std::list<int>::const_iterator i=used.begin();i!=used.end();++i)
00229   {
00230      int tag=*i;
00231      if (tag==0)    // check for the default tag
00232      {
00233     current_value=operation(current_value,DataMaths::reductionOp(vec,shape,data.getDefaultOffset(),operation,initial_value));
00234      }
00235      else
00236      {
00237     DataTagged::DataMapType::const_iterator it=lookup.find(tag);
00238     if (it!=lookup.end())
00239     {
00240         current_value=operation(current_value,DataMaths::reductionOp(vec,shape,it->second,operation,initial_value));
00241     }
00242      }
00243   }
00244   return current_value;
00245 }
00246 
00247 template <class BinaryFunction>
00248 inline
00249 double
00250 algorithm(DataConstant& data,
00251           BinaryFunction operation,
00252       double initial_value)
00253 {
00254   return DataMaths::reductionOp(data.getVectorRO(),data.getShape(),0,operation,initial_value);
00255 }
00256 
00268 template <class BinaryFunction>
00269 inline
00270 void
00271 dp_algorithm(const DataExpanded& data,
00272              DataExpanded& result,
00273              BinaryFunction operation,
00274          double initial_value)
00275 {
00276   int i,j;
00277   int numSamples=data.getNumSamples();
00278   int numDPPSample=data.getNumDPPSample();
00279 //  DataArrayView dataView=data.getPointDataView();
00280 //  DataArrayView resultView=result.getPointDataView();
00281   const DataTypes::ValueType& dataVec=data.getVectorRO();
00282   const DataTypes::ShapeType& shape=data.getShape();
00283   DataTypes::ValueType& resultVec=result.getVectorRW();
00284   // perform the operation on each data-point and assign
00285   // this to the corresponding element in result
00286   #pragma omp parallel for private(i,j) schedule(static)
00287   for (i=0;i<numSamples;i++) {
00288     for (j=0;j<numDPPSample;j++) {
00289 /*      resultView.getData(result.getPointOffset(i,j)) =
00290         dataView.reductionOp(data.getPointOffset(i,j),operation,initial_value);*/
00291       resultVec[result.getPointOffset(i,j)] =
00292         DataMaths::reductionOp(dataVec, shape, data.getPointOffset(i,j),operation,initial_value);
00293 
00294     }
00295   }
00296 }
00297 
00298 template <class BinaryFunction>
00299 inline
00300 void
00301 dp_algorithm(const DataTagged& data,
00302              DataTagged& result,
00303              BinaryFunction operation,
00304          double initial_value)
00305 {
00306   // perform the operation on each tagged value in data
00307   // and assign this to the corresponding element in result
00308   const DataTypes::ShapeType& shape=data.getShape();
00309   const DataTypes::ValueType& vec=data.getVectorRO();
00310   const DataTagged::DataMapType& lookup=data.getTagLookup();
00311   for (DataTagged::DataMapType::const_iterator i=lookup.begin(); i!=lookup.end(); i++) {
00312     result.getDataByTagRW(i->first,0) =
00313     DataMaths::reductionOp(vec,shape,data.getOffsetForTag(i->first),operation,initial_value);
00314   }
00315   // perform the operation on the default data value
00316   // and assign this to the default element in result
00317   result.getVectorRW()[result.getDefaultOffset()] = DataMaths::reductionOp(data.getVectorRO(),data.getShape(),data.getDefaultOffset(),operation,initial_value);
00318 }
00319 
00320 template <class BinaryFunction>
00321 inline
00322 void
00323 dp_algorithm(DataConstant& data,
00324              DataConstant& result,
00325              BinaryFunction operation,
00326          double initial_value)
00327 {
00328   // perform the operation on the data value
00329   // and assign this to the element in result
00330   result.getVectorRW()[0] =
00331     DataMaths::reductionOp(data.getVectorRO(),data.getShape(),0,operation,initial_value);
00332 }
00333 
00334 } // end of namespace
00335 
00336 #endif