CppAD: A C++ Algorithmic Differentiation Package 20110419
|
00001 /* $Id$ */ 00002 # ifndef CPPAD_POW_INCLUDED 00003 # define CPPAD_POW_INCLUDED 00004 00005 /* -------------------------------------------------------------------------- 00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-11 Bradley M. Bell 00007 00008 CppAD is distributed under multiple licenses. This distribution is under 00009 the terms of the 00010 Common Public License Version 1.0. 00011 00012 A copy of this license is included in the COPYING file of this distribution. 00013 Please visit http://www.coin-or.org/CppAD/ for information on other licenses. 00014 -------------------------------------------------------------------------- */ 00015 00016 /* 00017 $begin pow$$ 00018 $spell 00019 Vec 00020 std 00021 namespace 00022 CppAD 00023 const 00024 $$ 00025 00026 $index pow, AD$$ 00027 $index exponent, AD function$$ 00028 00029 $section The AD Power Function$$ 00030 00031 $head Syntax$$ 00032 $syntax%%z% = pow(%x%, %y%)%$$ 00033 00034 00035 $head Purpose$$ 00036 Determines the value of the power function which is defined by 00037 $latex \[ 00038 {\rm pow} (x, y) = x^y 00039 \] $$ 00040 This version of the $code pow$$ function may use 00041 logarithms and exponentiation to compute derivatives. 00042 This will not work if $italic x$$ is less than or equal zero. 00043 If the value of $italic y$$ is an integer, 00044 the $cref/pow_int/$$ function is used to compute this value 00045 using only multiplication (and division if $italic y$$ is negative). 00046 (This will work even if $italic x$$ is less than or equal zero.) 00047 00048 $head x$$ 00049 The argument $italic x$$ has the following prototype 00050 $syntax% 00051 const %Type% &%x% 00052 %$$ 00053 where $italic Type$$ is 00054 $syntax%VecAD<%Base%>::reference%$$, 00055 $syntax%AD<%Base%>%$$, 00056 $syntax%%Base%$$, 00057 $syntax%double%$$, 00058 or 00059 $syntax%int%$$. 00060 00061 $head y$$ 00062 The argument $italic y$$ has the following prototype 00063 $syntax% 00064 const %Type% &%y% 00065 %$$ 00066 where $italic Type$$ is 00067 $syntax%VecAD<%Base%>::reference%$$, 00068 $syntax%AD<%Base%>%$$, 00069 $syntax%%Base%$$, 00070 $syntax%double%$$, 00071 or 00072 $syntax%int%$$. 00073 00074 $head z$$ 00075 The result $italic z$$ has prototype 00076 $syntax% 00077 AD<%Base%> %z% 00078 %$$ 00079 00080 $head Standard Types$$ 00081 A definition for the $code pow$$ function is included 00082 in the CppAD namespace for the case where both $italic x$$ 00083 and $italic y$$ have the same type and that type is 00084 $code float$$ or $code double$$. 00085 00086 $head Operation Sequence$$ 00087 This is an AD of $italic Base$$ 00088 $xref/glossary/Operation/Atomic/atomic operation/1/$$ 00089 and hence is part of the current 00090 AD of $italic Base$$ 00091 $xref/glossary/Operation/Sequence/operation sequence/1/$$. 00092 00093 $head Example$$ 00094 $children% 00095 example/pow.cpp% 00096 example/pow_int.cpp 00097 %$$ 00098 The files 00099 $xref/Pow.cpp/$$, $cref/pow_int.cpp/$$ 00100 contain an examples and tests of this function. 00101 They returns true if they succeed and false otherwise. 00102 00103 $end 00104 ------------------------------------------------------------------------------- 00105 */ 00106 00107 // BEGIN CppAD namespace 00108 namespace CppAD { 00109 00110 // copy of standard functions in CppAD namespace 00111 inline float pow(const float &x, const float &y) 00112 { return std::pow(x, y); } 00113 00114 inline double pow(const double &x, const double &y) 00115 { return std::pow(x, y); } 00116 00117 // case where x and y are AD<Base> ----------------------------------------- 00118 template <class Base> AD<Base> 00119 pow(const AD<Base> &x, const AD<Base> &y) 00120 { ADTape<Base> *tape = AD<Base>::tape_ptr(); 00121 size_t tape_id = 0; 00122 if( tape != CPPAD_NULL ) 00123 tape_id = tape->id_; 00124 00125 // id_ setting for parameters cannot match 0 00126 bool var_x = x.id_ == tape_id; 00127 bool var_y = y.id_ == tape_id; 00128 CPPAD_ASSERT_KNOWN( 00129 Parameter(x) || var_x , 00130 "pow: first argument is a variable for a different thread" 00131 ); 00132 CPPAD_ASSERT_KNOWN( 00133 Parameter(y) || var_y , 00134 "pow: second argument is a variable for a different thread" 00135 ); 00136 00137 AD<Base> result; 00138 result.value_ = pow(x.value_, y.value_); 00139 CPPAD_ASSERT_UNKNOWN( Parameter(result) ); 00140 00141 if( var_x ) 00142 { if( var_y ) 00143 { // result = variable^variable 00144 CPPAD_ASSERT_UNKNOWN( NumRes(PowvvOp) == 3 ); 00145 CPPAD_ASSERT_UNKNOWN( NumArg(PowvvOp) == 2 ); 00146 00147 // put operand addresses in tape 00148 tape->Rec_.PutArg(x.taddr_, y.taddr_); 00149 00150 // put operator in the tape 00151 result.taddr_ = tape->Rec_.PutOp(PowvvOp); 00152 00153 // make result a variable 00154 result.id_ = tape_id; 00155 } 00156 else if( IdenticalZero( y.value_ ) ) 00157 { // result = variable^0 00158 } 00159 else 00160 { // result = variable^parameter 00161 CPPAD_ASSERT_UNKNOWN( NumRes(PowvpOp) == 3 ); 00162 CPPAD_ASSERT_UNKNOWN( NumArg(PowvpOp) == 2 ); 00163 00164 // put operand addresses in tape 00165 size_t p = tape->Rec_.PutPar(y.value_); 00166 tape->Rec_.PutArg(x.taddr_, p); 00167 00168 // put operator in the tape 00169 result.taddr_ = tape->Rec_.PutOp(PowvpOp); 00170 00171 // make result a variable 00172 result.id_ = tape_id; 00173 } 00174 } 00175 else if( var_y ) 00176 { if( IdenticalZero(x.value_) ) 00177 { // result = 0^variable 00178 } 00179 else 00180 { // result = variable^parameter 00181 CPPAD_ASSERT_UNKNOWN( NumRes(PowpvOp) == 3 ); 00182 CPPAD_ASSERT_UNKNOWN( NumArg(PowpvOp) == 2 ); 00183 00184 // put operand addresses in tape 00185 size_t p = tape->Rec_.PutPar(x.value_); 00186 tape->Rec_.PutArg(p, y.taddr_); 00187 00188 // put operator in the tape 00189 result.taddr_ = tape->Rec_.PutOp(PowpvOp); 00190 00191 // make result a variable 00192 result.id_ = tape_id; 00193 } 00194 } 00195 return result; 00196 } 00197 // ========================================================================= 00198 // Fold operations in same way as CPPAD_FOLD_AD_VALUED_BINARY_OPERATION(Op) 00199 // ------------------------------------------------------------------------- 00200 // Operations with VecAD_reference<Base> and AD<Base> only 00201 00202 template <class Base> AD<Base> 00203 pow(const AD<Base> &x, const VecAD_reference<Base> &y) 00204 { return pow(x, y.ADBase()); } 00205 00206 template <class Base> AD<Base> 00207 pow(const VecAD_reference<Base> &x, const VecAD_reference<Base> &y) 00208 { return pow(x.ADBase(), y.ADBase()); } 00209 00210 template <class Base> AD<Base> 00211 pow(const VecAD_reference<Base> &x, const AD<Base> &y) 00212 { return pow(x.ADBase(), y); } 00213 // ------------------------------------------------------------------------- 00214 // Operations with Base 00215 00216 template <class Base> AD<Base> 00217 pow(const Base &x, const AD<Base> &y) 00218 { return pow(AD<Base>(x), y); } 00219 00220 template <class Base> AD<Base> 00221 pow(const Base &x, const VecAD_reference<Base> &y) 00222 { return pow(AD<Base>(x), y.ADBase()); } 00223 00224 template <class Base> AD<Base> 00225 pow(const AD<Base> &x, const Base &y) 00226 { return pow(x, AD<Base>(y)); } 00227 00228 template <class Base> AD<Base> 00229 pow(const VecAD_reference<Base> &x, const Base &y) 00230 { return pow(x.ADBase(), AD<Base>(y)); } 00231 // ------------------------------------------------------------------------- 00232 // Operations with double 00233 00234 template <class Base> AD<Base> 00235 pow(const double &x, const AD<Base> &y) 00236 { return pow(AD<Base>(x), y); } 00237 00238 template <class Base> AD<Base> 00239 pow(const double &x, const VecAD_reference<Base> &y) 00240 { return pow(AD<Base>(x), y.ADBase()); } 00241 00242 template <class Base> AD<Base> 00243 pow(const AD<Base> &x, const double &y) 00244 { return pow(x, AD<Base>(y)); } 00245 00246 template <class Base> AD<Base> 00247 pow(const VecAD_reference<Base> &x, const double &y) 00248 { return pow(x.ADBase(), AD<Base>(y)); } 00249 // ------------------------------------------------------------------------- 00250 // Special case to avoid ambuigity when Base is double 00251 00252 inline AD<double> 00253 pow(const double &x, const AD<double> &y) 00254 { return pow(AD<double>(x), y); } 00255 00256 inline AD<double> 00257 pow(const double &x, const VecAD_reference<double> &y) 00258 { return pow(AD<double>(x), y.ADBase()); } 00259 00260 inline AD<double> 00261 pow(const AD<double> &x, const double &y) 00262 { return pow(x, AD<double>(y)); } 00263 00264 inline AD<double> 00265 pow(const VecAD_reference<double> &x, const double &y) 00266 { return pow(x.ADBase(), AD<double>(y)); } 00267 00268 // ========================================================================= 00269 // Fold operations for the cases where x is an int, 00270 // but let cppad/pow_int.hpp handle the cases where y is an int. 00271 // ------------------------------------------------------------------------- 00272 template <class Base> AD<Base> pow 00273 (const int &x, const VecAD_reference<Base> &y) 00274 { return pow(AD<Base>(x), y.ADBase()); } 00275 00276 template <class Base> AD<Base> pow 00277 (const int &x, const AD<Base> &y) 00278 { return pow(AD<Base>(x), y); } 00279 00280 } // END CppAD namespace 00281 00282 # endif