00001 #include "SundanceADReal.hpp"
00002 #include "PlayaExceptions.hpp"
00003
00004
00005 using namespace Sundance;
00006 using namespace Teuchos;
00007
00008 ADReal ADReal::operator-() const
00009 {
00010 ADReal rtn = *this;
00011 rtn.value_ = -rtn.value_;
00012 rtn.gradient_ = -gradient_;
00013
00014 addFlops(1 + gradient_.dim());
00015 return rtn;
00016 }
00017
00018 ADReal ADReal::operator+(const ADReal& other) const
00019 {
00020 ADReal rtn(*this);
00021 rtn += other;
00022 return rtn;
00023 }
00024
00025 ADReal ADReal::operator-(const ADReal& other) const
00026 {
00027 ADReal rtn(*this);
00028 rtn -= other;
00029 return rtn;
00030 }
00031
00032 ADReal ADReal::operator*(const ADReal& other) const
00033 {
00034 ADReal rtn(*this);
00035 rtn *= other;
00036 return rtn;
00037 }
00038
00039 ADReal ADReal::operator/(const ADReal& other) const
00040 {
00041 ADReal rtn(*this);
00042 rtn /= other;
00043 return rtn;
00044 }
00045
00046 ADReal& ADReal::operator+=(const ADReal& other)
00047 {
00048 value_ += other.value_;
00049 gradient_ += other.gradient_;
00050
00051 addFlops(1 + gradient_.dim());
00052
00053 return *this;
00054 }
00055
00056 ADReal& ADReal::operator-=(const ADReal& other)
00057 {
00058 value_ -= other.value_;
00059 gradient_ -= other.gradient_;
00060
00061 addFlops(1 + gradient_.dim());
00062
00063 return *this;
00064 }
00065
00066 ADReal& ADReal::operator*=(const ADReal& other)
00067 {
00068 gradient_ = (other.value_*gradient_ + value_*other.gradient_);
00069 value_ *= other.value_;
00070
00071 addFlops(1 + 3*gradient_.dim());
00072
00073 return *this;
00074 }
00075
00076 ADReal& ADReal::operator/=(const ADReal& other)
00077 {
00078 TEUCHOS_TEST_FOR_EXCEPTION(other.value_ == 0.0, std::runtime_error,
00079 "ADReal::operator/=() division by zero");
00080
00081 gradient_ = (gradient_/other.value_
00082 - value_/(other.value_*other.value_) * other.gradient_);
00083 value_ /= other.value_;
00084
00085 addFlops(3 + 3*gradient_.dim());
00086
00087 return *this;
00088 }
00089
00090
00091
00092
00093 ADReal& ADReal::operator+=(const double& other)
00094 {
00095 value_ += other;
00096 addFlops(1);
00097 return *this;
00098 }
00099
00100 ADReal& ADReal::operator-=(const double& other)
00101 {
00102 value_ -= other;
00103
00104 addFlops(1);
00105 return *this;
00106 }
00107
00108 ADReal& ADReal::operator*=(const double& other)
00109 {
00110 gradient_ *= other;
00111 value_ *= other;
00112
00113 addFlops(1 + gradient_.dim());
00114
00115 return *this;
00116 }
00117
00118 ADReal& ADReal::operator/=(const double& other)
00119 {
00120 TEUCHOS_TEST_FOR_EXCEPTION(other == 0.0, std::runtime_error,
00121 "ADReal::operator/=() division by zero");
00122
00123 addFlops(2 + gradient_.dim());
00124 gradient_ *= 1.0/other;
00125 value_ /= other;
00126 return *this;
00127 }
00128
00129
00130 ADReal ADReal::operator+(const double& other) const
00131 {
00132 ADReal rtn(*this);
00133 rtn += other;
00134 return rtn;
00135 }
00136
00137 ADReal ADReal::operator-(const double& other) const
00138 {
00139 ADReal rtn(*this);
00140 rtn -= other;
00141 return rtn;
00142 }
00143
00144 ADReal ADReal::operator*(const double& other) const
00145 {
00146 ADReal rtn(*this);
00147 rtn *= other;
00148 return rtn;
00149 }
00150
00151 ADReal ADReal::operator/(const double& other) const
00152 {
00153 ADReal rtn(*this);
00154 rtn /= other;
00155 return rtn;
00156 }
00157
00158
00159 namespace Sundance
00160 {
00161 ADReal operator+(const double& scalar, const ADReal& a)
00162 {
00163 return a+scalar;
00164 }
00165
00166 ADReal operator-(const double& scalar, const ADReal& a)
00167 {
00168 return -a+scalar;
00169 }
00170
00171 ADReal operator*(const double& scalar, const ADReal& a)
00172 {
00173 return a*scalar;
00174 }
00175
00176 ADReal operator/(const double& scalar, const ADReal& a)
00177 {
00178 ADReal rtn(a);
00179 rtn.reciprocate();
00180 return rtn*scalar;
00181 }
00182 }
00183
00184 void ADReal::reciprocate()
00185 {
00186 TEUCHOS_TEST_FOR_EXCEPTION(value_==0.0, std::runtime_error,
00187 "div by 0 in ADReal::reciprocate()");
00188
00189 addFlops(1 + gradient_.dim());
00190
00191 gradient_ /= (value_*value_);
00192 gradient_ *= -1.0;
00193 value_ = 1.0/value_;
00194 }
00195
00196
00197
00198
00199
00200
00201
00202