CppAD: A C++ Algorithmic Differentiation Package  20171217
log1p_op.hpp
Go to the documentation of this file.
4
5 /* --------------------------------------------------------------------------
7
9 the terms of the
10  Eclipse Public License Version 1.0.
11
12 A copy of this license is included in the COPYING file of this distribution.
14 -------------------------------------------------------------------------- */
15
17 /*!
18 \file log1p_op.hpp
19 Forward and reverse mode calculations for z = log1p(x).
20 */
21
22 /*!
23 Compute forward mode Taylor coefficient for result of op = Log1pOp.
24
25 The C++ source code corresponding to this operation is
26 \verbatim
27  z = log1p(x)
28 \endverbatim
29
31 */
32 template <class Base>
33 inline void forward_log1p_op(
34  size_t p ,
35  size_t q ,
36  size_t i_z ,
37  size_t i_x ,
38  size_t cap_order ,
39  Base* taylor )
40 {
41  size_t k;
42
43  // check assumptions
46  CPPAD_ASSERT_UNKNOWN( q < cap_order );
47  CPPAD_ASSERT_UNKNOWN( p <= q );
48
49  // Taylor coefficients corresponding to argument and result
50  Base* x = taylor + i_x * cap_order;
51  Base* z = taylor + i_z * cap_order;
52
53  if( p == 0 )
54  { z[0] = log1p( x[0] );
55  p++;
56  if( q == 0 )
57  return;
58  }
59  if ( p == 1 )
60  { z[1] = x[1] / (Base(1.0) + x[0]);
61  p++;
62  }
63  for(size_t j = p; j <= q; j++)
64  {
65  z[j] = -z[1] * x[j-1];
66  for(k = 2; k < j; k++)
67  z[j] -= Base(double(k)) * z[k] * x[j-k];
68  z[j] /= Base(double(j));
69  z[j] += x[j];
70  z[j] /= (Base(1.0) + x[0]);
71  }
72 }
73
74 /*!
75 Muiltiple directions Taylor coefficient for op = Log1pOp.
76
77 The C++ source code corresponding to this operation is
78 \verbatim
79  z = log1p(x)
80 \endverbatim
81
83 */
84 template <class Base>
86  size_t q ,
87  size_t r ,
88  size_t i_z ,
89  size_t i_x ,
90  size_t cap_order ,
91  Base* taylor )
92 {
93
94  // check assumptions
97  CPPAD_ASSERT_UNKNOWN( 0 < q );
98  CPPAD_ASSERT_UNKNOWN( q < cap_order );
99
100  // Taylor coefficients corresponding to argument and result
101  size_t num_taylor_per_var = (cap_order-1) * r + 1;
102  Base* x = taylor + i_x * num_taylor_per_var;
103  Base* z = taylor + i_z * num_taylor_per_var;
104
105  size_t m = (q-1) * r + 1;
106  for(size_t ell = 0; ell < r; ell++)
107  { z[m+ell] = Base(double(q)) * x[m+ell];
108  for(size_t k = 1; k < q; k++)
109  z[m+ell] -= Base(double(k)) * z[(k-1)*r+1+ell] * x[(q-k-1)*r+1+ell];
110  z[m+ell] /= (Base(double(q)) + Base(q) * x[0]);
111  }
112 }
113
114 /*!
115 Compute zero order forward mode Taylor coefficient for result of op = Log1pOp.
116
117 The C++ source code corresponding to this operation is
118 \verbatim
119  z = log1p(x)
120 \endverbatim
121
123 */
124 template <class Base>
125 inline void forward_log1p_op_0(
126  size_t i_z ,
127  size_t i_x ,
128  size_t cap_order ,
129  Base* taylor )
130 {
131
132  // check assumptions
135  CPPAD_ASSERT_UNKNOWN( 0 < cap_order );
136
137  // Taylor coefficients corresponding to argument and result
138  Base* x = taylor + i_x * cap_order;
139  Base* z = taylor + i_z * cap_order;
140
141  z[0] = log1p( x[0] );
142 }
143
144 /*!
145 Compute reverse mode partial derivatives for result of op = Log1pOp.
146
147 The C++ source code corresponding to this operation is
148 \verbatim
149  z = log1p(x)
150 \endverbatim
151
153 */
154
155 template <class Base>
156 inline void reverse_log1p_op(
157  size_t d ,
158  size_t i_z ,
159  size_t i_x ,
160  size_t cap_order ,
161  const Base* taylor ,
162  size_t nc_partial ,
163  Base* partial )
164 { size_t j, k;
165
166  // check assumptions
169  CPPAD_ASSERT_UNKNOWN( d < cap_order );
170  CPPAD_ASSERT_UNKNOWN( d < nc_partial );
171
172  // Taylor coefficients and partials corresponding to argument
173  const Base* x = taylor + i_x * cap_order;
174  Base* px = partial + i_x * nc_partial;
175
176  // Taylor coefficients and partials corresponding to result
177  const Base* z = taylor + i_z * cap_order;
178  Base* pz = partial + i_z * nc_partial;
179
180  Base inv_1px0 = Base(1.0) / (Base(1) + x[0]);
181
182  j = d;
183  while(j)
184  { // scale partial w.r.t z[j]
185  pz[j] = azmul(pz[j] , inv_1px0);
186
187  px[0] -= azmul(pz[j], z[j]);
188  px[j] += pz[j];
189
190  // further scale partial w.r.t. z[j]
191  pz[j] /= Base(double(j));
192
193  for(k = 1; k < j; k++)
194  { pz[k] -= Base(double(k)) * azmul(pz[j], x[j-k]);
195  px[j-k] -= Base(double(k)) * azmul(pz[j], z[k]);
196  }
197  --j;
198  }
199  px[0] += azmul(pz[0], inv_1px0);
200 }
201
203 # endif
204 # endif
Definition: azmul.hpp:94
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_log1p_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 = Log1pOp.
Definition: log1p_op.hpp:85
void forward_log1p_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 = Log1pOp.
Definition: log1p_op.hpp:125
void reverse_log1p_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 = Log1pOp.
Definition: log1p_op.hpp:156