00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #ifndef SUNDANCE_ADREAL_H
00032 #define SUNDANCE_ADREAL_H
00033
00034
00035 #include "SundanceDefs.hpp"
00036 #include "SundancePoint.hpp"
00037 #include "SundanceStdMathFunctors.hpp"
00038 #include "Teuchos_RefCountPtr.hpp"
00039
00040 namespace Sundance
00041 {
00042
00043
00044
00045
00046 class ADReal
00047 {
00048 public:
00049
00050 ADReal() : value_(0), gradient_(0.0, 0.0, 0.0){;}
00051
00052
00053 ADReal(double value, int direction, int spatialDimension)
00054 : value_(value), gradient_()
00055 {
00056 if (spatialDimension==1) gradient_ = Point(0.0);
00057 if (spatialDimension==2) gradient_ = Point(0.0, 0.0);
00058 if (spatialDimension==3) gradient_ = Point(0.0, 0.0, 0.0);
00059 gradient_[direction] = 1.0;
00060 }
00061
00062 ADReal(double value, int spatialDimension)
00063 : value_(value), gradient_()
00064 {
00065 if (spatialDimension==1) gradient_ = Point(0.0);
00066 if (spatialDimension==2) gradient_ = Point(0.0, 0.0);
00067 if (spatialDimension==3) gradient_ = Point(0.0, 0.0, 0.0);
00068 }
00069
00070
00071 ADReal(double value, const Point& gradient)
00072 : value_(value), gradient_(gradient) {;}
00073
00074
00075 ADReal operator-() const ;
00076
00077 ADReal& operator+=(const ADReal& other) ;
00078
00079 ADReal& operator-=(const ADReal& other) ;
00080
00081 ADReal& operator*=(const ADReal& other) ;
00082
00083 ADReal& operator/=(const ADReal& other) ;
00084
00085
00086 ADReal& operator+=(const double& scalar) ;
00087
00088 ADReal& operator-=(const double& scalar) ;
00089
00090 ADReal& operator*=(const double& scalar) ;
00091
00092 ADReal& operator/=(const double& scalar) ;
00093
00094
00095 ADReal operator+(const ADReal& other) const ;
00096
00097 ADReal operator-(const ADReal& other) const ;
00098
00099 ADReal operator*(const ADReal& other) const ;
00100
00101 ADReal operator/(const ADReal& other) const ;
00102
00103
00104 ADReal operator+(const double& scalar) const ;
00105
00106 ADReal operator-(const double& scalar) const ;
00107
00108 ADReal operator*(const double& scalar) const ;
00109
00110 ADReal operator/(const double& scalar) const ;
00111
00112
00113 const double& value() const {return value_;}
00114
00115 const Point& gradient() const {return gradient_;}
00116
00117 void reciprocate() ;
00118
00119 static double& totalFlops() {static double rtn = 0; return rtn;}
00120
00121
00122 static void addFlops(const double& flops) {totalFlops() += flops;}
00123 private:
00124 double value_;
00125 Point gradient_;
00126 };
00127
00128
00129
00130
00131 ADReal operator+(const double& scalar,
00132 const ADReal& a);
00133
00134 ADReal operator-(const double& scalar,
00135 const ADReal& a);
00136
00137 ADReal operator*(const double& scalar,
00138 const ADReal& a);
00139
00140 ADReal operator/(const double& scalar,
00141 const ADReal& a);
00142
00143
00144 inline ADReal pow(const ADReal& x, const double& y)
00145 {
00146 Teuchos::RCP<UnaryFunctor> func = Teuchos::rcp(new PowerFunctor(y));
00147 double f;
00148 double df;
00149 double xVal = x.value();
00150 func->eval1(&xVal, 1, &f, &df);
00151 return ADReal(f, df*x.gradient());
00152 }
00153
00154 #define SUNDANCE_AD_FUNCTOR(opName, functorName) \
00155 inline ADReal opName(const ADReal& x)\
00156 {\
00157 Teuchos::RCP<UnaryFunctor> func = Teuchos::rcp(new functorName());\
00158 double f;\
00159 double df;\
00160 double xVal = x.value();\
00161 func->eval1(&xVal, 1, &f, &df);\
00162 return ADReal(f, df*x.gradient());}
00163
00164
00165 SUNDANCE_AD_FUNCTOR(exp, StdExp)
00166
00167 SUNDANCE_AD_FUNCTOR(log, StdLog)
00168
00169 SUNDANCE_AD_FUNCTOR(sqrt, StdSqrt)
00170
00171 SUNDANCE_AD_FUNCTOR(sin, StdSin)
00172
00173 SUNDANCE_AD_FUNCTOR(cos, StdCos)
00174
00175
00176
00177 }
00178
00179
00180 namespace std
00181 {
00182 inline ostream& operator<<(std::ostream& os, const Sundance::ADReal& x)
00183 {
00184 std::cerr << "AD[" << x.value() << ", grad="<< x.gradient() << "]";
00185 return os;
00186 }
00187 }
00188
00189
00190
00191 #endif
00192
00193
00194