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