CppAD: A C++ Algorithmic Differentiation Package  20171217
sqrt_op.hpp
Go to the documentation of this file.
3
4 /* --------------------------------------------------------------------------
6
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.
13 -------------------------------------------------------------------------- */
14
15
17 /*!
18 \file sqrt_op.hpp
19 Forward and reverse mode calculations for z = sqrt(x).
20 */
21
22
23 /*!
24 Compute forward mode Taylor coefficient for result of op = SqrtOp.
25
26 The C++ source code corresponding to this operation is
27 \verbatim
28  z = sqrt(x)
29 \endverbatim
30
32 */
33 template <class Base>
34 inline void forward_sqrt_op(
35  size_t p ,
36  size_t q ,
37  size_t i_z ,
38  size_t i_x ,
39  size_t cap_order ,
40  Base* taylor )
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  size_t k;
53  if( p == 0 )
54  { z[0] = sqrt( x[0] );
55  p++;
56  }
57  for(size_t j = p; j <= q; j++)
58  {
59  z[j] = Base(0.0);
60  for(k = 1; k < j; k++)
61  z[j] -= Base(double(k)) * z[k] * z[j-k];
62  z[j] /= Base(double(j));
63  z[j] += x[j] / Base(2.0);
64  z[j] /= z[0];
65  }
66 }
67
68 /*!
69 Multiple direction forward mode Taylor coefficient for op = SqrtOp.
70
71 The C++ source code corresponding to this operation is
72 \verbatim
73  z = sqrt(x)
74 \endverbatim
75
77 */
78 template <class Base>
79 inline void forward_sqrt_op_dir(
80  size_t q ,
81  size_t r ,
82  size_t i_z ,
83  size_t i_x ,
84  size_t cap_order ,
85  Base* taylor )
86 {
87  // check assumptions
90  CPPAD_ASSERT_UNKNOWN( 0 < q );
91  CPPAD_ASSERT_UNKNOWN( q < cap_order );
92
93  // Taylor coefficients corresponding to argument and result
94  size_t num_taylor_per_var = (cap_order-1) * r + 1;
95  Base* z = taylor + i_z * num_taylor_per_var;
96  Base* x = taylor + i_x * num_taylor_per_var;
97
98  size_t m = (q-1) * r + 1;
99  for(size_t ell = 0; ell < r; ell++)
100  { z[m+ell] = Base(0.0);
101  for(size_t k = 1; k < q; k++)
102  z[m+ell] -= Base(double(k)) * z[(k-1)*r+1+ell] * z[(q-k-1)*r+1+ell];
103  z[m+ell] /= Base(double(q));
104  z[m+ell] += x[m+ell] / Base(2.0);
105  z[m+ell] /= z[0];
106  }
107 }
108
109 /*!
110 Compute zero order forward mode Taylor coefficient for result of op = SqrtOp.
111
112 The C++ source code corresponding to this operation is
113 \verbatim
114  z = sqrt(x)
115 \endverbatim
116
118 */
119 template <class Base>
120 inline void forward_sqrt_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] = sqrt( x[0] );
136 }
137 /*!
138 Compute reverse mode partial derivatives for result of op = SqrtOp.
139
140 The C++ source code corresponding to this operation is
141 \verbatim
142  z = sqrt(x)
143 \endverbatim
144
146 */
147
148 template <class Base>
149 inline void reverse_sqrt_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  Base* px = partial + i_x * nc_partial;
166
167  // Taylor coefficients and partials corresponding to result
168  const Base* z = taylor + i_z * cap_order;
169  Base* pz = partial + i_z * nc_partial;
170
171
172  Base inv_z0 = Base(1.0) / z[0];
173
174  // number of indices to access
175  size_t j = d;
176  size_t k;
177  while(j)
178  {
179
180  // scale partial w.r.t. z[j]
181  pz[j] = azmul(pz[j], inv_z0);
182
183  pz[0] -= azmul(pz[j], z[j]);
184  px[j] += pz[j] / Base(2.0);
185  for(k = 1; k < j; k++)
186  pz[k] -= azmul(pz[j], z[j-k]);
187  --j;
188  }
189  px[0] += azmul(pz[0], inv_z0) / Base(2.0);
190 }
191
193 # endif
void forward_sqrt_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 = SqrtOp.
Definition: sqrt_op.hpp:34
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