00001 # ifndef CPPAD_DIV_OP_INCLUDED 00002 # define CPPAD_DIV_OP_INCLUDED 00003 00004 /* -------------------------------------------------------------------------- 00005 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-07 Bradley M. Bell 00006 00007 CppAD is distributed under multiple licenses. This distribution is under 00008 the terms of the 00009 Common Public License Version 1.0. 00010 00011 A copy of this license is included in the COPYING file of this distribution. 00012 Please visit http://www.coin-or.org/CppAD/ for information on other licenses. 00013 -------------------------------------------------------------------------- */ 00014 00015 /* 00016 $begin ForDivvvOp$$ $comment CppAD Developer Documentation$$ 00017 $spell 00018 Div 00019 Divpv 00020 Divvp 00021 Divvv 00022 Taylor 00023 const 00024 inline 00025 Op 00026 $$ 00027 00028 $index divide, forward operator$$ 00029 $index forward, divide operator$$ 00030 $index operator, divide forward$$ 00031 $index ForDiv$$ 00032 00033 $section Forward Mode Division Operator$$ 00034 00035 $head Syntax$$ 00036 00037 $syntax%inline void ForDivvvOp(size_t %d%, 00038 %Base% *%z%, const %Base% *%x%, const %Base% *%y%)%$$ 00039 $pre 00040 $$ 00041 $syntax%inline void ForDivpvOp(size_t %d%, 00042 %Base% *%z%, const %Base% *%p%, const %Base% *%y%)%$$ 00043 $pre 00044 $$ 00045 $syntax%inline void ForDivvpOp(size_t %d%, 00046 %Base% *%z%, const %Base% *%x%, const %Base% *%p%)%$$ 00047 00048 00049 $head Description$$ 00050 Computes the $italic d$$ order Taylor coefficient for $latex Z$$ where 00051 $table 00052 Operation $cnext Value $rnext 00053 Divvv $cnext $latex Z = X / Y$$ $rnext 00054 Divpv $cnext $latex Z = P / Y$$ $rnext 00055 Divvp $cnext $latex Z = X / P$$ 00056 $tend 00057 00058 $head x$$ 00059 The vector $italic x$$ has length $latex d+1$$ and contains the 00060 $th d$$ order Taylor coefficient row vector for $italic X$$. 00061 00062 $head y$$ 00063 The vector $italic y$$ has length $latex d+1$$ and contains the 00064 $th d$$ order Taylor coefficient row vector for $italic Y$$. 00065 00066 $head p$$ 00067 The scalar $syntax%*%p%$$ contains the value of the parameter $italic P$$. 00068 00069 00070 $head z$$ 00071 The vector $italic z$$ has length $latex d+1$$. 00072 On input it contains the 00073 $th d-1$$ order Taylor coefficient row vector for $italic Z$$. 00074 On output it contains the 00075 $th d$$ order Taylor coefficient row vector for $italic Z$$; i.e., 00076 $syntax%%z%[%d%]%$$ is set equal to the $th d$$ Taylor coefficient for 00077 the function $italic Z$$. 00078 00079 $end 00080 ------------------------------------------------------------------------------ 00081 $begin RevDivvvOp$$ $comment CppAD Developer Documentation$$ 00082 $spell 00083 Div 00084 Divpv 00085 Divvp 00086 Divvv 00087 Taylor 00088 const 00089 inline 00090 Op 00091 px 00092 py 00093 pz 00094 $$ 00095 00096 00097 $index divide, reverse operator$$ 00098 $index reverse, divide operator$$ 00099 $index operator, divide reverse$$ 00100 $index ForDiv$$ 00101 00102 $section Reverse Mode Division Operator$$ 00103 00104 $head Syntax$$ 00105 00106 $syntax%inline void RevDivvvOp(size_t %d%, 00107 const %Base% *%z%, const %Base% *%x%, const %Base% *%y%, 00108 %Base% *%pz%, %Base% *%px%, %Base% *%py%)%$$ 00109 $pre 00110 $$ 00111 $syntax%inline void RevDivpvOp(size_t %d%, 00112 const %Base% *%z%, const %Base% *%p%, const %Base% *%y%, 00113 %Base% *%pz%, %Base% *%py%)%$$ 00114 $pre 00115 $$ 00116 $syntax%inline void RevDivvpOp(size_t %d%, 00117 const %Base% *%z%, const %Base% *%x%, const %Base% *%p%, 00118 const %Base% *%pz%, %Base% *%px%)%$$ 00119 00120 $head Description$$ 00121 We are given the partial derivatives for a function 00122 $latex G(z, x, y)$$ and we wish to compute the partial derivatives for 00123 the function 00124 $latex \[ 00125 H(x, y) = G [ Z(x, y) , x , y ] 00126 \]$$ 00127 where $latex Z(x, y)$$ is defined as the 00128 $th d$$ order Taylor coefficient row vector for $italic Z$$ as 00129 a function of the corresponding vectors for 00130 $italic X$$ and $italic Y$$ where 00131 00132 $table 00133 Operation $cnext Value $rnext 00134 Divvv $cnext $latex Z = X / Y$$ $rnext 00135 Divvp $cnext $latex Z = P / Y$$ $rnext 00136 Divpv $cnext $latex Z = X / P$$ 00137 $tend 00138 00139 Note that $italic Z$$ has been used both the original division 00140 function and for the corresponding mapping of Taylor coefficients. 00141 00142 $head z$$ 00143 The vector $italic z$$ has length $latex d+1$$ and contains the 00144 $th d$$ order Taylor coefficient row vector for $italic Z$$. 00145 00146 $head x$$ 00147 The vector $italic x$$ has length $latex d+1$$ and contains the 00148 $th d$$ order Taylor coefficient row vector for $italic X$$. 00149 00150 $head y$$ 00151 The vector $italic y$$ has length $latex d+1$$ and contains the 00152 $th d$$ order Taylor coefficient row vector for $italic Y$$. 00153 00154 $head p$$ 00155 The scalar $syntax%*%p%$$ contains the value of the parameter $italic P$$. 00156 00157 00158 $head On Input$$ 00159 00160 $subhead px$$ 00161 The vector $italic px$$ has length $latex d+1$$ and 00162 $syntax%%px%[%j%]%$$ contains the partial for $italic G$$ 00163 with respect to the $th j$$ order Taylor coefficient for $italic X$$. 00164 00165 $subhead py$$ 00166 The vector $italic py$$ has length $latex d+1$$ and 00167 $syntax%%py%[%j%]%$$ contains the partial for $italic G$$ 00168 with respect to the $th j$$ order Taylor coefficient for $italic Y$$. 00169 00170 $head pz$$ 00171 The vector $italic pz$$ has length $latex d+1$$ and 00172 $syntax%%pz%[%j%]%$$ contains the partial for $italic G$$ 00173 with respect to the $th j$$ order Taylor coefficient for $italic Z$$. 00174 00175 $head On Output$$ 00176 00177 $subhead px$$ 00178 If present, 00179 the vector $italic px$$ has length $latex d+1$$ and 00180 $syntax%%px%[%j%]%$$ contains the partial for $italic H$$ 00181 with respect to the $th j$$ order Taylor coefficient for $italic X$$. 00182 00183 $subhead py$$ 00184 If present, 00185 the vector $italic py$$ has length $latex d+1$$ and 00186 $syntax%%py%[%j%]%$$ contains the partial for $italic H$$ 00187 with respect to the $th j$$ order Taylor coefficient for $italic Y$$. 00188 00189 00190 $subhead pz$$ 00191 The vector $italic pz$$ has length $latex d+1$$ and 00192 its contents are no longer specified value; i.e., it has been used 00193 for work space. 00194 00195 $end 00196 ------------------------------------------------------------------------------ 00197 */ 00198 00199 // BEGIN CppAD namespace 00200 namespace CppAD { 00201 00202 // --------------------------- Divvv ----------------------------------------- 00203 00204 template <class Base> 00205 inline void ForDivvvOp(size_t j, 00206 Base *z, const Base *x, const Base *y) 00207 { 00208 size_t k; 00209 00210 # if 0 00211 // 05-01-08 removed divide by zero check for better use of CondExp 00212 CPPAD_ASSERT_KNOWN( 00213 y[0] != Base(0), 00214 "Attempt to divide by zero" 00215 ); 00216 # endif 00217 00218 z[j] = x[j]; 00219 for(k = 1; k <= j; k++) 00220 z[j] -= z[j-k] * y[k]; 00221 z[j] /= y[0]; 00222 } 00223 00224 template <class Base> 00225 inline void RevDivvvOp(size_t d, 00226 const Base *z, const Base *x, const Base *y, 00227 Base *pz, Base *px, Base *py) 00228 { size_t k; 00229 00230 // number of indices to access 00231 size_t j = d + 1; 00232 00233 # if 0 00234 // 05-01-08 removed divide by zero check for better use of CondExp 00235 CPPAD_ASSERT_KNOWN( 00236 y[0] != Base(0), 00237 "Attempt to divide by zero" 00238 ); 00239 # endif 00240 00241 while(j) 00242 { --j; 00243 // scale partial w.r.t. z[j] 00244 pz[j] /= y[0]; 00245 00246 px[j] += pz[j]; 00247 for(k = 1; k <= j; k++) 00248 { pz[j-k] -= pz[j] * y[k]; 00249 py[k] -= pz[j] * z[j-k]; 00250 } 00251 py[0] -= pz[j] * z[j]; 00252 } 00253 } 00254 00255 // --------------------------- Divpv ---------------------------------------- 00256 00257 template <class Base> 00258 inline void ForDivpvOp(size_t j, 00259 Base *z, const Base *p, const Base *y) 00260 { size_t k; 00261 00262 # if 0 00263 // 05-01-08 removed divide by zero check for better use of CondExp 00264 CPPAD_ASSERT_KNOWN( 00265 y[0] != Base(0), 00266 "Attempt to divide by zero" 00267 ); 00268 # endif 00269 00270 if( j == 0 ) 00271 z[j] = (*p); 00272 else z[j] = Base(0); 00273 for(k = 1; k <= j; k++) 00274 z[j] -= z[j-k] * y[k]; 00275 z[j] /= y[0]; 00276 00277 } 00278 00279 template <class Base> 00280 inline void RevDivpvOp(size_t d, 00281 const Base *z, const Base *p, const Base *y, 00282 Base *pz, Base *py) 00283 { size_t k; 00284 00285 // number of indices to access 00286 size_t j = d + 1; 00287 00288 # if 0 00289 // 05-01-08 removed divide by zero check for better use of CondExp 00290 CPPAD_ASSERT_KNOWN( 00291 y[0] != Base(0), 00292 "Attempt to divide by zero" 00293 ); 00294 # endif 00295 00296 while(j) 00297 { --j; 00298 // scale partial w.r.t z[j] 00299 pz[j] /= y[0]; 00300 00301 for(k = 1; k <= j; k++) 00302 { pz[j-k] -= pz[j] * y[k]; 00303 py[k] -= pz[j] * z[j-k]; 00304 } 00305 py[0] -= pz[j] * z[j]; 00306 } 00307 } 00308 00309 // --------------------------- Divvp ------------------------------------------ 00310 00311 template <class Base> 00312 inline void ForDivvpOp(size_t j, 00313 Base *z, const Base *x, const Base *p) 00314 { 00315 z[j] = x[j] / (*p); 00316 } 00317 00318 template <class Base> 00319 inline void RevDivvpOp(size_t d, 00320 const Base *z, const Base *x, const Base *p, 00321 const Base *pz, Base *px) 00322 { 00323 // number of indices to access 00324 size_t j = d + 1; 00325 00326 # if 0 00327 // 05-01-08 removed divide by zero check for better use of CondExp 00328 // should catch this case in zero order forward mode 00329 CPPAD_ASSERT_UNKNOWN( (*p) != Base(0) ); 00330 # endif 00331 00332 while(j) 00333 { --j; 00334 px[j] += pz[j] / (*p); 00335 } 00336 } 00337 00338 } // END CppAD namespace 00339 00340 # endif