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