00001 # ifndef CPPAD_ACOS_OP_INCLUDED 00002 # define CPPAD_ACOS_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 ForAcosOp$$ $comment CppAD Developer Documentation$$ 00017 $spell 00018 sqrt 00019 acos 00020 Taylor 00021 const 00022 inline 00023 Op 00024 $$ 00025 00026 $index acos, forward$$ 00027 $index forward, acos$$ 00028 $index ForAcosOp$$ 00029 00030 $section Forward Mode Acos Function$$ 00031 00032 $head Syntax$$ 00033 00034 $syntax%inline void ForAcosOp(size_t %d%, 00035 %Base% *%z%, %Base% *%b%, const %Base% *%x%)%$$ 00036 00037 $head Description$$ 00038 Computes the $italic d$$ order Taylor coefficient for 00039 $latex Z$$ and $latex B$$ where 00040 $syntax% 00041 %Z% = acos(%X%) 00042 %B% = \sqrt{ 1 - %X% * %X% } 00043 %$$ 00044 00045 $head x$$ 00046 The vector $italic x$$ has length $latex d+1$$ and contains the 00047 $th d$$ order Taylor coefficient row vector for $italic X$$. 00048 00049 $head z$$ 00050 The vector $italic z$$ has length $latex d+1$$. 00051 On input it contains the 00052 $th d-1$$ order Taylor coefficient row vector for $italic Z$$. 00053 On output it contains the 00054 $th d$$ order Taylor coefficient row vector for $italic Z$$; i.e., 00055 $syntax%%z%[%d%]%$$ is set equal to the $th d$$ Taylor coefficient for 00056 the function $italic S$$. 00057 00058 $head b$$ 00059 The vector $italic c$$ has length $latex d+1$$. 00060 On input it contains the 00061 $th d-1$$ order Taylor coefficient row vector for $italic B$$. 00062 On output it contains the 00063 $th d$$ order Taylor coefficient row vector for $italic B$$; i.e., 00064 $syntax%%b%[%d%]%$$ is set equal to the $th d$$ Taylor coefficient for 00065 the function $italic B$$. 00066 00067 $end 00068 ------------------------------------------------------------------------------ 00069 $begin RevAcosOp$$ $comment CppAD Developer Documentation$$ 00070 $spell 00071 ps 00072 Acos 00073 pb 00074 Sin 00075 Taylor 00076 const 00077 inline 00078 Op 00079 px 00080 py 00081 pz 00082 $$ 00083 00084 00085 $index acos, reverse$$ 00086 $index reverse, acos$$ 00087 $index RevAcosOp$$ 00088 00089 $section Reverse Mode Acos Function$$ 00090 00091 $head Syntax$$ 00092 00093 $syntax%inline void RevAcosOp(size_t %d%, 00094 const %Base% *%z%, const %Base% *%b%, const %Base% *%x%, 00095 %Base% *%pz%, %Base% *%pb%, %Base% *%px%)%$$ 00096 00097 $head Description$$ 00098 We are given the partial derivatives for a function 00099 $latex G(z, b, x)$$ and we wish to compute the partial derivatives for 00100 the function 00101 $latex \[ 00102 H(x) = G [ Z(x) , B(x) , x ] 00103 \] $$ 00104 where $latex Z(x)$$ and $latex B(x)$$ are defined as the 00105 $th d$$ order Taylor coefficient row vector for $latex \arccos(x)$$ 00106 and $latex 1 + x * x$$ 00107 as a function of the corresponding row vector for $italic X$$; i.e., 00108 $latex \[ 00109 \begin{array}{rcl} 00110 Z & = & \arccos(X) \\ 00111 B & = & 1 + X * X 00112 \end{array} 00113 \]$$ 00114 Note that $italic Z$$ and $latex B$$ have 00115 been used both the original 00116 functions and for the corresponding mapping of Taylor coefficients. 00117 00118 $head x$$ 00119 The vector $italic x$$ has length $latex d+1$$ and contains the 00120 $th d$$ order Taylor coefficient row vector for $italic X$$. 00121 00122 00123 $head z$$ 00124 The vector $italic z$$ has length $latex d+1$$ and contains 00125 $th d$$ order Taylor coefficient row vector for $italic z$$. 00126 00127 $head b$$ 00128 The vector $italic b$$ has length $latex d+1$$ and contains 00129 $th d$$ order Taylor coefficient row vector for $italic B$$. 00130 00131 00132 $head On Input$$ 00133 00134 $subhead px$$ 00135 The vector $italic px$$ has length $latex d+1$$ and 00136 $syntax%%px%[%j%]%$$ contains the partial for $italic G$$ 00137 with respect to the $th j$$ order Taylor coefficient for $italic X$$. 00138 00139 $subhead pz$$ 00140 The vector $italic pz$$ has length $latex d+1$$ and 00141 $syntax%%pz%[%j%]%$$ contains the partial for $italic G$$ 00142 with respect to the $th j$$ order Taylor coefficient for $italic Z$$. 00143 00144 $subhead pb$$ 00145 The vector $italic pb$$ has length $latex d+1$$ and 00146 $syntax%%pb%[%j%]%$$ contains the partial for $italic G$$ 00147 with respect to the $th j$$ order Taylor coefficient for $italic B$$. 00148 00149 $head On Output$$ 00150 00151 $subhead px$$ 00152 The vector $italic px$$ has length $latex d+1$$ and 00153 $syntax%%px%[%j%]%$$ contains the partial for $italic H$$ 00154 with respect to the $th j$$ order Taylor coefficient for $italic X$$. 00155 00156 $subhead pz$$ 00157 The vector $italic ps$$ has length $latex d+1$$ and 00158 its contents are no longer specified; i.e., it has 00159 been used for work space. 00160 00161 $subhead pb$$ 00162 The vector $italic pb$$ has length $latex d+1$$ and 00163 its contents are no longer specified; i.e., it has 00164 been used for work space. 00165 00166 $end 00167 ------------------------------------------------------------------------------ 00168 */ 00169 00170 // BEGIN CppAD namespace 00171 namespace CppAD { 00172 00173 template <class Base> 00174 inline void ForAcosOp(size_t j, 00175 Base *z, Base *b, const Base *x) 00176 { size_t k; 00177 Base qj; 00178 00179 if( j == 0 ) 00180 { z[j] = acos( x[0] ); 00181 qj = Base(1) - x[0] * x[0]; 00182 b[j] = sqrt( qj ); 00183 } 00184 else 00185 { qj = 0.; 00186 for(k = 0; k <= j; k++) 00187 qj -= x[k] * x[j-k]; 00188 b[j] = Base(0); 00189 z[j] = Base(0); 00190 for(k = 1; k < j; k++) 00191 { b[j] -= Base(k) * b[k] * b[j-k]; 00192 z[j] -= Base(k) * z[k] * b[j-k]; 00193 } 00194 b[j] /= Base(j); 00195 z[j] /= Base(j); 00196 // 00197 b[j] += qj / Base(2); 00198 z[j] -= x[j]; 00199 // 00200 b[j] /= b[0]; 00201 z[j] /= b[0]; 00202 } 00203 } 00204 00205 template <class Base> 00206 inline void RevAcosOp(size_t d, 00207 const Base *z, const Base *b, const Base *x, 00208 Base *pz, Base *pb, Base *px) 00209 { size_t k; 00210 00211 // number of indices to access 00212 size_t j = d; 00213 00214 while(j) 00215 { 00216 // scale partials w.r.t b[j] by 1 / b[0] 00217 pb[j] /= b[0]; 00218 00219 // scale partials w.r.t z[j] by 1 / b[0] 00220 pz[j] /= b[0]; 00221 00222 // update partials w.r.t b^0 00223 pb[0] -= pz[j] * z[j] + pb[j] * b[j]; 00224 00225 // update partial w.r.t. x^0 00226 px[0] -= pb[j] * x[j]; 00227 00228 // update partial w.r.t. x^j 00229 px[j] -= pz[j] + pb[j] * x[0]; 00230 00231 // further scale partial w.r.t. z[j] by 1 / j 00232 pz[j] /= Base(j); 00233 00234 for(k = 1; k < j; k++) 00235 { // update partials w.r.t b^(j-k) 00236 pb[j-k] -= Base(k) * pz[j] * z[k] + pb[j] * b[k]; 00237 00238 // update partials w.r.t. x^k 00239 px[k] -= pb[j] * x[j-k]; 00240 00241 // update partials w.r.t. z^k 00242 pz[k] -= pz[j] * Base(k) * b[j-k]; 00243 } 00244 --j; 00245 } 00246 00247 // j == 0 case 00248 px[0] -= ( pz[0] + pb[0] * x[0]) / b[0]; 00249 } 00250 00251 } // END CppAD namespace 00252 00253 # endif