CppAD: A C++ Algorithmic Differentiation Package  20171217
expm1_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
16
18 /*!
19 \file expm1_op.hpp
20 Forward and reverse mode calculations for z = expm1(x).
21 */
22
23
24 /*!
25 Forward mode Taylor coefficient for result of op = Expm1Op.
26
27 The C++ source code corresponding to this operation is
28 \verbatim
29  z = expm1(x)
30 \endverbatim
31
33 */
34 template <class Base>
35 inline void forward_expm1_op(
36  size_t p ,
37  size_t q ,
38  size_t i_z ,
39  size_t i_x ,
40  size_t cap_order ,
41  Base* taylor )
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  size_t k;
54  if( p == 0 )
55  { z[0] = expm1( x[0] );
56  p++;
57  }
58  for(size_t j = p; j <= q; j++)
59  {
60  z[j] = x[1] * z[j-1];
61  for(k = 2; k <= j; k++)
62  z[j] += Base(double(k)) * x[k] * z[j-k];
63  z[j] /= Base(double(j));
64  z[j] += x[j];
65  }
66 }
67
68
69 /*!
70 Multiple direction forward mode Taylor coefficient for op = Expm1Op.
71
72 The C++ source code corresponding to this operation is
73 \verbatim
74  z = expm1(x)
75 \endverbatim
76
78 */
79 template <class Base>
81  size_t q ,
82  size_t r ,
83  size_t i_z ,
84  size_t i_x ,
85  size_t cap_order ,
86  Base* taylor )
87 {
88  // check assumptions
91  CPPAD_ASSERT_UNKNOWN( q < cap_order );
92  CPPAD_ASSERT_UNKNOWN( 0 < q );
93
94  // Taylor coefficients corresponding to argument and result
95  size_t num_taylor_per_var = (cap_order-1) * r + 1;
96  Base* x = taylor + i_x * num_taylor_per_var;
97  Base* z = taylor + i_z * num_taylor_per_var;
98
99  size_t m = (q-1)*r + 1;
100  for(size_t ell = 0; ell < r; ell++)
101  { z[m+ell] = Base(double(q)) * x[m+ell] * z[0];
102  for(size_t k = 1; k < q; k++)
103  z[m+ell] += Base(double(k)) * x[(k-1)*r+ell+1] * z[(q-k-1)*r+ell+1];
104  z[m+ell] /= Base(double(q));
105  z[m+ell] += x[m+ell];
106  }
107 }
108
109 /*!
110 Zero order forward mode Taylor coefficient for result of op = Expm1Op.
111
112 The C++ source code corresponding to this operation is
113 \verbatim
114  z = expm1(x)
115 \endverbatim
116
118 */
119 template <class Base>
120 inline void forward_expm1_op_0(
121  size_t i_z ,
122  size_t i_x ,
123  size_t cap_order ,
124  Base* taylor )
125 {
126  // check assumptions
129  CPPAD_ASSERT_UNKNOWN( 0 < cap_order );
130
131  // Taylor coefficients corresponding to argument and result
132  Base* x = taylor + i_x * cap_order;
133  Base* z = taylor + i_z * cap_order;
134
135  z[0] = expm1( x[0] );
136 }
137 /*!
138 Reverse mode partial derivatives for result of op = Expm1Op.
139
140 The C++ source code corresponding to this operation is
141 \verbatim
142  z = expm1(x)
143 \endverbatim
144
146 */
147
148 template <class Base>
149 inline void reverse_expm1_op(
150  size_t d ,
151  size_t i_z ,
152  size_t i_x ,
153  size_t cap_order ,
154  const Base* taylor ,
155  size_t nc_partial ,
156  Base* partial )
157 {
158  // check assumptions
161  CPPAD_ASSERT_UNKNOWN( d < cap_order );
162  CPPAD_ASSERT_UNKNOWN( d < nc_partial );
163
164  // Taylor coefficients and partials corresponding to argument
165  const Base* x = taylor + i_x * cap_order;
166  Base* px = partial + i_x * nc_partial;
167
168  // Taylor coefficients and partials corresponding to result
169  const Base* z = taylor + i_z * cap_order;
170  Base* pz = partial + i_z * nc_partial;
171
172  // If pz is zero, make sure this operation has no effect
173  // (zero times infinity or nan would be non-zero).
174  bool skip(true);
175  for(size_t i_d = 0; i_d <= d; i_d++)
176  skip &= IdenticalZero(pz[i_d]);
177  if( skip )
178  return;
179
180  // loop through orders in reverse
181  size_t j, k;
182  j = d;
183  while(j)
184  { px[j] += pz[j];
185
186  // scale partial w.r.t z[j]
187  pz[j] /= Base(double(j));
188
189  for(k = 1; k <= j; k++)
190  { px[k] += Base(double(k)) * azmul(pz[j], z[j-k]);
191  pz[j-k] += Base(double(k)) * azmul(pz[j], x[k]);
192  }
193  --j;
194  }
195  px[0] += pz[0] + azmul(pz[0], z[0]);
196 }
197
199 # endif
200 # endif
Definition: azmul.hpp:94
void forward_expm1_op_dir(size_t q, size_t r, size_t i_z, size_t i_x, size_t cap_order, Base *taylor)
Multiple direction forward mode Taylor coefficient for op = Expm1Op.
Definition: expm1_op.hpp:80
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
bool IdenticalZero(const std::complex< double > &x)
void reverse_expm1_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)
Reverse mode partial derivatives for result of op = Expm1Op.
Definition: expm1_op.hpp:149
std::complex< double > expm1(const std::complex< double > &x)