CppAD: A C++ Algorithmic Differentiation Package
20130102
|
00001 /* $Id$ */ 00002 # ifndef CPPAD_LOG_OP_INCLUDED 00003 # define CPPAD_LOG_OP_INCLUDED 00004 00005 /* -------------------------------------------------------------------------- 00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-12 Bradley M. Bell 00007 00008 CppAD is distributed under multiple licenses. This distribution is under 00009 the terms of the 00010 Eclipse 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 CPPAD_BEGIN_NAMESPACE 00017 /*! 00018 \defgroup log_op_hpp log_op.hpp 00019 \{ 00020 \file log_op.hpp 00021 Forward and reverse mode calculations for z = log(x). 00022 */ 00023 00024 /*! 00025 Compute forward mode Taylor coefficient for result of op = LogOp. 00026 00027 The C++ source code corresponding to this operation is 00028 \verbatim 00029 z = log(x) 00030 \endverbatim 00031 00032 \copydetails forward_unary1_op 00033 */ 00034 template <class Base> 00035 inline void forward_log_op( 00036 size_t j , 00037 size_t i_z , 00038 size_t i_x , 00039 size_t nc_taylor , 00040 Base* taylor ) 00041 { 00042 size_t k; 00043 00044 // check assumptions 00045 CPPAD_ASSERT_UNKNOWN( NumArg(LogOp) == 1 ); 00046 CPPAD_ASSERT_UNKNOWN( NumRes(LogOp) == 1 ); 00047 CPPAD_ASSERT_UNKNOWN( i_x < i_z ); 00048 CPPAD_ASSERT_UNKNOWN( j < nc_taylor ); 00049 00050 // Taylor coefficients corresponding to argument and result 00051 Base* x = taylor + i_x * nc_taylor; 00052 Base* z = taylor + i_z * nc_taylor; 00053 00054 if( j == 0 ) 00055 z[0] = log( x[0] ); 00056 else if ( j == 1 ) 00057 z[1] = x[1] / x[0]; 00058 else 00059 { 00060 z[j] = -z[1] * x[j-1]; 00061 for(k = 2; k < j; k++) 00062 z[j] -= Base(k) * z[k] * x[j-k]; 00063 z[j] /= Base(j); 00064 z[j] += x[j]; 00065 z[j] /= x[0]; 00066 } 00067 } 00068 00069 /*! 00070 Compute zero order forward mode Taylor coefficient for result of op = LogOp. 00071 00072 The C++ source code corresponding to this operation is 00073 \verbatim 00074 z = log(x) 00075 \endverbatim 00076 00077 \copydetails forward_unary1_op_0 00078 */ 00079 template <class Base> 00080 inline void forward_log_op_0( 00081 size_t i_z , 00082 size_t i_x , 00083 size_t nc_taylor , 00084 Base* taylor ) 00085 { 00086 00087 // check assumptions 00088 CPPAD_ASSERT_UNKNOWN( NumArg(LogOp) == 1 ); 00089 CPPAD_ASSERT_UNKNOWN( NumRes(LogOp) == 1 ); 00090 CPPAD_ASSERT_UNKNOWN( i_x < i_z ); 00091 CPPAD_ASSERT_UNKNOWN( 0 < nc_taylor ); 00092 00093 // Taylor coefficients corresponding to argument and result 00094 Base* x = taylor + i_x * nc_taylor; 00095 Base* z = taylor + i_z * nc_taylor; 00096 00097 z[0] = log( x[0] ); 00098 } 00099 00100 /*! 00101 Compute reverse mode partial derivatives for result of op = LogOp. 00102 00103 The C++ source code corresponding to this operation is 00104 \verbatim 00105 z = log(x) 00106 \endverbatim 00107 00108 \copydetails reverse_unary1_op 00109 */ 00110 00111 template <class Base> 00112 inline void reverse_log_op( 00113 size_t d , 00114 size_t i_z , 00115 size_t i_x , 00116 size_t nc_taylor , 00117 const Base* taylor , 00118 size_t nc_partial , 00119 Base* partial ) 00120 { size_t j, k; 00121 00122 // check assumptions 00123 CPPAD_ASSERT_UNKNOWN( NumArg(LogOp) == 1 ); 00124 CPPAD_ASSERT_UNKNOWN( NumRes(LogOp) == 1 ); 00125 CPPAD_ASSERT_UNKNOWN( i_x < i_z ); 00126 CPPAD_ASSERT_UNKNOWN( d < nc_taylor ); 00127 CPPAD_ASSERT_UNKNOWN( d < nc_partial ); 00128 00129 // Taylor coefficients and partials corresponding to argument 00130 const Base* x = taylor + i_x * nc_taylor; 00131 Base* px = partial + i_x * nc_partial; 00132 00133 // Taylor coefficients and partials corresponding to result 00134 const Base* z = taylor + i_z * nc_taylor; 00135 Base* pz = partial + i_z * nc_partial; 00136 00137 j = d; 00138 while(j) 00139 { // scale partial w.r.t z[j] 00140 pz[j] /= x[0]; 00141 00142 px[0] -= pz[j] * z[j]; 00143 px[j] += pz[j]; 00144 00145 // further scale partial w.r.t. z[j] 00146 pz[j] /= Base(j); 00147 00148 for(k = 1; k < j; k++) 00149 { pz[k] -= pz[j] * Base(k) * x[j-k]; 00150 px[j-k] -= pz[j] * Base(k) * z[k]; 00151 } 00152 --j; 00153 } 00154 px[0] += pz[0] / x[0]; 00155 } 00156 00157 /*! \} */ 00158 CPPAD_END_NAMESPACE 00159 # endif