CppAD: A C++ Algorithmic Differentiation Package  20171217
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
log_op.hpp
Go to the documentation of this file.
1 # ifndef CPPAD_LOCAL_LOG_OP_HPP
2 # define CPPAD_LOCAL_LOG_OP_HPP
3 
4 /* --------------------------------------------------------------------------
5 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-17 Bradley M. Bell
6 
7 CppAD is distributed under multiple licenses. This distribution is under
8 the terms of the
9  Eclipse Public License Version 1.0.
10 
11 A copy of this license is included in the COPYING file of this distribution.
12 Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
13 -------------------------------------------------------------------------- */
14 
15 namespace CppAD { namespace local { // BEGIN_CPPAD_LOCAL_NAMESPACE
16 /*!
17 \file log_op.hpp
18 Forward and reverse mode calculations for z = log(x).
19 */
20 
21 /*!
22 Compute forward mode Taylor coefficient for result of op = LogOp.
23 
24 The C++ source code corresponding to this operation is
25 \verbatim
26  z = log(x)
27 \endverbatim
28 
29 \copydetails CppAD::local::forward_unary1_op
30 */
31 template <class Base>
32 inline void forward_log_op(
33  size_t p ,
34  size_t q ,
35  size_t i_z ,
36  size_t i_x ,
37  size_t cap_order ,
38  Base* taylor )
39 {
40  size_t k;
41 
42  // check assumptions
45  CPPAD_ASSERT_UNKNOWN( q < cap_order );
46  CPPAD_ASSERT_UNKNOWN( p <= q );
47 
48  // Taylor coefficients corresponding to argument and result
49  Base* x = taylor + i_x * cap_order;
50  Base* z = taylor + i_z * cap_order;
51 
52  if( p == 0 )
53  { z[0] = log( x[0] );
54  p++;
55  if( q == 0 )
56  return;
57  }
58  if ( p == 1 )
59  { z[1] = x[1] / x[0];
60  p++;
61  }
62  for(size_t j = p; j <= q; j++)
63  {
64  z[j] = -z[1] * x[j-1];
65  for(k = 2; k < j; k++)
66  z[j] -= Base(double(k)) * z[k] * x[j-k];
67  z[j] /= Base(double(j));
68  z[j] += x[j];
69  z[j] /= x[0];
70  }
71 }
72 
73 /*!
74 Muiltiple directions Taylor coefficient for op = LogOp.
75 
76 The C++ source code corresponding to this operation is
77 \verbatim
78  z = log(x)
79 \endverbatim
80 
81 \copydetails CppAD::local::forward_unary1_op_dir
82 */
83 template <class Base>
84 inline void forward_log_op_dir(
85  size_t q ,
86  size_t r ,
87  size_t i_z ,
88  size_t i_x ,
89  size_t cap_order ,
90  Base* taylor )
91 {
92 
93  // check assumptions
96  CPPAD_ASSERT_UNKNOWN( 0 < q );
97  CPPAD_ASSERT_UNKNOWN( q < cap_order );
98 
99  // Taylor coefficients corresponding to argument and result
100  size_t num_taylor_per_var = (cap_order-1) * r + 1;
101  Base* x = taylor + i_x * num_taylor_per_var;
102  Base* z = taylor + i_z * num_taylor_per_var;
103 
104  size_t m = (q-1) * r + 1;
105  for(size_t ell = 0; ell < r; ell++)
106  { z[m+ell] = Base(double(q)) * x[m+ell];
107  for(size_t k = 1; k < q; k++)
108  z[m+ell] -= Base(double(k)) * z[(k-1)*r+1+ell] * x[(q-k-1)*r+1+ell];
109  z[m+ell] /= (Base(double(q)) * x[0]);
110  }
111 }
112 
113 /*!
114 Compute zero order forward mode Taylor coefficient for result of op = LogOp.
115 
116 The C++ source code corresponding to this operation is
117 \verbatim
118  z = log(x)
119 \endverbatim
120 
121 \copydetails CppAD::local::forward_unary1_op_0
122 */
123 template <class Base>
124 inline void forward_log_op_0(
125  size_t i_z ,
126  size_t i_x ,
127  size_t cap_order ,
128  Base* taylor )
129 {
130 
131  // check assumptions
134  CPPAD_ASSERT_UNKNOWN( 0 < cap_order );
135 
136  // Taylor coefficients corresponding to argument and result
137  Base* x = taylor + i_x * cap_order;
138  Base* z = taylor + i_z * cap_order;
139 
140  z[0] = log( x[0] );
141 }
142 
143 /*!
144 Compute reverse mode partial derivatives for result of op = LogOp.
145 
146 The C++ source code corresponding to this operation is
147 \verbatim
148  z = log(x)
149 \endverbatim
150 
151 \copydetails CppAD::local::reverse_unary1_op
152 */
153 
154 template <class Base>
155 inline void reverse_log_op(
156  size_t d ,
157  size_t i_z ,
158  size_t i_x ,
159  size_t cap_order ,
160  const Base* taylor ,
161  size_t nc_partial ,
162  Base* partial )
163 { size_t j, k;
164 
165  // check assumptions
168  CPPAD_ASSERT_UNKNOWN( d < cap_order );
169  CPPAD_ASSERT_UNKNOWN( d < nc_partial );
170 
171  // Taylor coefficients and partials corresponding to argument
172  const Base* x = taylor + i_x * cap_order;
173  Base* px = partial + i_x * nc_partial;
174 
175  // Taylor coefficients and partials corresponding to result
176  const Base* z = taylor + i_z * cap_order;
177  Base* pz = partial + i_z * nc_partial;
178 
179  Base inv_x0 = Base(1.0) / x[0];
180 
181  j = d;
182  while(j)
183  { // scale partial w.r.t z[j]
184  pz[j] = azmul(pz[j] , inv_x0);
185 
186  px[0] -= azmul(pz[j], z[j]);
187  px[j] += pz[j];
188 
189  // further scale partial w.r.t. z[j]
190  pz[j] /= Base(double(j));
191 
192  for(k = 1; k < j; k++)
193  { pz[k] -= Base(double(k)) * azmul(pz[j], x[j-k]);
194  px[j-k] -= Base(double(k)) * azmul(pz[j], z[k]);
195  }
196  --j;
197  }
198  px[0] += azmul(pz[0], inv_x0);
199 }
200 
201 } } // END_CPPAD_LOCAL_NAMESPACE
202 # endif
AD< Base > azmul(const AD< Base > &x, const AD< Base > &y)
Definition: azmul.hpp:94
AD< Base > log(const AD< Base > &x)
size_t NumArg(OpCode op)
Number of arguments for a specified operator.
Definition: op_code.hpp:175
size_t NumRes(OpCode op)
Number of variables resulting from the specified operation.
Definition: op_code.hpp:281
void forward_log_op(size_t p, size_t q, size_t i_z, size_t i_x, size_t cap_order, Base *taylor)
Compute forward mode Taylor coefficient for result of op = LogOp.
Definition: log_op.hpp:32
void forward_log_op_0(size_t i_z, size_t i_x, size_t cap_order, Base *taylor)
Compute zero order forward mode Taylor coefficient for result of op = LogOp.
Definition: log_op.hpp:124
#define CPPAD_ASSERT_UNKNOWN(exp)
Check that exp is true, if not terminate execution.
void forward_log_op_dir(size_t q, size_t r, size_t i_z, size_t i_x, size_t cap_order, Base *taylor)
Muiltiple directions Taylor coefficient for op = LogOp.
Definition: log_op.hpp:84
void reverse_log_op(size_t d, size_t i_z, size_t i_x, size_t cap_order, const Base *taylor, size_t nc_partial, Base *partial)
Compute reverse mode partial derivatives for result of op = LogOp.
Definition: log_op.hpp:155