1 # ifndef CPPAD_LOCAL_SPARSE_UNARY_OP_HPP
2 # define CPPAD_LOCAL_SPARSE_UNARY_OP_HPP
3 /* --------------------------------------------------------------------------
4 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-20 Bradley M. Bell
5
6 CppAD is distributed under the terms of the
7 Eclipse Public License Version 2.0.
8
9 This Source Code may also be made available under the following
10 Secondary License when the conditions for such availability set forth
11 in the Eclipse Public License, Version 2.0 are satisfied:
12 GNU General Public License, Version 2.0 or later.
13 ---------------------------------------------------------------------------- */
14
15 // BEGIN_CPPAD_LOCAL_SPARSE_NAMESPACE
16 namespace CppAD { namespace local { namespace sparse {
17 /*!
18 \file sparse_unary_op.hpp
19 Forward and reverse mode sparsity patterns for unary operators.
20 */
21
22
23 /*!
24 Forward mode Jacobian sparsity pattern for all unary operators.
25
26 The C++ source code corresponding to a unary operation has the form
27 \verbatim
28 z = fun(x)
29 \endverbatim
30 where fun is a C++ unary function, or it has the form
31 \verbatim
32 z = x op q
33 \endverbatim
34 where op is a C++ binary unary operator and q is a parameter.
35
36 \tparam Vector_set
37 is the type used for vectors of sets. It can be either
38 sparse::pack_setvec or sparse::list_setvec.
39
40 \param i_z
41 variable index corresponding to the result for this operation;
42 i.e., z.
43
44 \param i_x
45 variable index corresponding to the argument for this operator;
46 i.e., x.
47
48
49 \param sparsity
50 \b Input: The set with index arg[0] in sparsity
51 is the sparsity bit pattern for x.
52 This identifies which of the independent variables the variable x
53 depends on.
54 \n
55 \n
56 \b Output: The set with index i_z in sparsity
57 is the sparsity bit pattern for z.
58 This identifies which of the independent variables the variable z
59 depends on.
60 \n
61
62 \par Checked Assertions:
63 \li i_x < i_z
64 */
65
66 template <class Vector_set>
for_jac_unary_op(size_t i_z,size_t i_x,Vector_set & sparsity)67 void for_jac_unary_op(
68 size_t i_z ,
69 size_t i_x ,
70 Vector_set& sparsity )
71 {
72 // check assumptions
73 CPPAD_ASSERT_UNKNOWN( i_x < i_z );
74
75 sparsity.assignment(i_z, i_x, sparsity);
76 }
77 /*!
78 Reverse mode Jacobian sparsity pattern for all unary operators.
79
80 The C++ source code corresponding to a unary operation has the form
81 \verbatim
82 z = fun(x)
83 \endverbatim
84 where fun is a C++ unary function, or it has the form
85 \verbatim
86 z = x op q
87 \endverbatim
88 where op is a C++ bianry operator and q is a parameter.
89
90 This routine is given the sparsity patterns
91 for a function G(z, y, ... )
92 and it uses them to compute the sparsity patterns for
93 \verbatim
94 H( x , w , u , ... ) = G[ z(x) , x , w , u , ... ]
95 \endverbatim
96
97 \tparam Vector_set
98 is the type used for vectors of sets. It can be either
99 sparse::pack_setvec or sparse::list_setvec.
100
101
102 \param i_z
103 variable index corresponding to the result for this operation;
104 i.e. the row index in sparsity corresponding to z.
105
106 \param i_x
107 variable index corresponding to the argument for this operator;
108 i.e. the row index in sparsity corresponding to x.
109
110 \param sparsity
111 \b Input:
112 The set with index i_z in sparsity
113 is the sparsity bit pattern for G with respect to the variable z.
114 \n
115 \b Input:
116 The set with index i_x in sparsity
117 is the sparsity bit pattern for G with respect to the variable x.
118 \n
119 \b Output:
120 The set with index i_x in sparsity
121 is the sparsity bit pattern for H with respect to the variable x.
122
123 \par Checked Assertions:
124 \li i_x < i_z
125 */
126
127 template <class Vector_set>
rev_jac_unary_op(size_t i_z,size_t i_x,Vector_set & sparsity)128 void rev_jac_unary_op(
129 size_t i_z ,
130 size_t i_x ,
131 Vector_set& sparsity )
132 {
133 // check assumptions
134 CPPAD_ASSERT_UNKNOWN( i_x < i_z );
135
136 sparsity.binary_union(i_x, i_x, i_z, sparsity);
137
138 return;
139 }
140 // ---------------------------------------------------------------------------
141 /*!
142 Reverse mode Hessian sparsity pattern for linear unary operators.
143
144 The C++ source code corresponding to this operation is
145 \verbatim
146 z = fun(x)
147 \endverbatim
148 where fun is a linear functions; e.g. abs, or
149 \verbatim
150 z = x op q
151 \endverbatim
152 where op is a C++ binary operator and q is a parameter.
153
154 \copydetails CppAD::local::reverse_sparse_hessian_unary_op
155 */
156 template <class Vector_set>
rev_hes_lin_unary_op(size_t i_z,size_t i_x,bool * rev_jacobian,const Vector_set & for_jac_sparsity,Vector_set & rev_hes_sparsity)157 void rev_hes_lin_unary_op(
158 size_t i_z ,
159 size_t i_x ,
160 bool* rev_jacobian ,
161 const Vector_set& for_jac_sparsity ,
162 Vector_set& rev_hes_sparsity )
163 {
164 // check assumptions
165 CPPAD_ASSERT_UNKNOWN( i_x < i_z );
166
167 // check for no effect
168 if( ! rev_jacobian[i_z] )
169 return;
170
171 rev_hes_sparsity.binary_union(i_x, i_x, i_z, rev_hes_sparsity);
172
173 rev_jacobian[i_x] = true;
174 return;
175 }
176
177 /*!
178 Reverse mode Hessian sparsity pattern for non-linear unary operators.
179
180 The C++ source code corresponding to this operation is
181 \verbatim
182 z = fun(x)
183 \endverbatim
184 where fun is a non-linear functions; e.g. sin. or
185 \verbatim
186 z = q / x
187 \endverbatim
188 where q is a parameter.
189
190
191 \copydetails CppAD::local::reverse_sparse_hessian_unary_op
192 */
193 template <class Vector_set>
rev_hes_nl_unary_op(size_t i_z,size_t i_x,bool * rev_jacobian,const Vector_set & for_jac_sparsity,Vector_set & rev_hes_sparsity)194 void rev_hes_nl_unary_op(
195 size_t i_z ,
196 size_t i_x ,
197 bool* rev_jacobian ,
198 const Vector_set& for_jac_sparsity ,
199 Vector_set& rev_hes_sparsity )
200 {
201 // check assumptions
202 CPPAD_ASSERT_UNKNOWN( i_x < i_z );
203
204 // check for no effect
205 if( ! rev_jacobian[i_z] )
206 return;
207
208 rev_hes_sparsity.binary_union(i_x, i_x, i_z, rev_hes_sparsity);
209 rev_hes_sparsity.binary_union(i_x, i_x, i_x, for_jac_sparsity);
210
211 rev_jacobian[i_x] = true;
212 return;
213 }
214 // ---------------------------------------------------------------------------
215 /*
216 $begin for_hes_nl_unary_op$$
217 $spell
218 hes
219 nl
220 op
221 np
222 numvar
223 Jacobian
224 $$
225
226 $section Forward Hessian Sparsity for Non-linear Unary Operators$$
227
228 $head Syntax$$
229 $codei%local::for_hes_nl_unary_op(
230 %np1%, %numvar%, %i_v%, %for_sparsity%
231 )%$$
232
233 $head Prototype$$
234 $srcthisfile%
235 0%// BEGIN_for_hes_nl_unary_op%// END_for_hes_nl_unary_op%1
236 %$$
237
238 $head C++ Source$$
239 The C++ source code corresponding to this operation is
240 $codei%
241 %w% = %fun%( %v% )
242 %$$
243 where $icode fun$$ is a non-linear function.
244
245 $head np1$$
246 This is the number of independent variables plus one;
247 i.e. size of $icode x$$ plus one.
248
249 $head numvar$$
250 This is the total number of variables in the tape.
251
252 $head i_w$$
253 is the index of the variable corresponding to the result $icode w$$.
254
255 $head i_v$$
256 is the index of the variable corresponding to the argument $icode v$$.
257
258 $head for_sparsity$$
259 We have the conditions $icode%np1% = %for_sparsity%.end()%$$
260 and $icode%for_sparsity%.n_set() = %np1% + %numvar%$$.
261
262 $subhead Input Jacobian Sparsity$$
263 For $icode%i%= 0, ..., %i_w%-1%$$,
264 the $icode%np1%+%i%$$ row of $icode for_sparsity$$ is the Jacobian sparsity
265 for the $th i$$ variable. These values do not change.
266 Note that $icode%i%=0%$$ corresponds to a parameter and
267 the corresponding Jacobian sparsity is empty.
268
269 $subhead Input Hessian Sparsity$$
270 For $icode%j%=1, ..., %n%$$,
271 the $th j$$ row of $icode for_sparsity$$ is the Hessian sparsity
272 before including the function $latex w(x)$$.
273
274 $subhead Output Jacobian Sparsity$$
275 the $icode i_w$$ row of $icode for_sparsity$$ is the Jacobian sparsity
276 for the variable $icode w$$.
277
278 $subhead Output Hessian Sparsity$$
279 For $icode%j%=1, ..., %n%$$,
280 the $th j$$ row of $icode for_sparsity$$ is the Hessian sparsity
281 after including the function $latex w(x)$$.
282
283 $end
284 */
285 // BEGIN_for_hes_nl_unary_op
286 template <class Vector_set>
for_hes_nl_unary_op(size_t np1,size_t numvar,size_t i_w,size_t i_v,Vector_set & for_sparsity)287 void for_hes_nl_unary_op(
288 size_t np1 ,
289 size_t numvar ,
290 size_t i_w ,
291 size_t i_v ,
292 Vector_set& for_sparsity )
293 // END_for_hes_nl_unary_op
294 { CPPAD_ASSERT_UNKNOWN( i_v < i_w );
295 CPPAD_ASSERT_UNKNOWN( i_w < numvar );
296 CPPAD_ASSERT_UNKNOWN( for_sparsity.end() == np1 );
297 CPPAD_ASSERT_UNKNOWN( for_sparsity.n_set() == np1 + numvar );
298 CPPAD_ASSERT_UNKNOWN( for_sparsity.number_elements(np1) == 0 );
299
300 // set Jacobian sparsity J(i_w)
301 for_sparsity.assignment(np1 + i_w, np1 + i_v, for_sparsity);
302
303 // set of independent variables that v depends on
304 typename Vector_set::const_iterator itr(for_sparsity, i_v + np1);
305
306 // loop over independent variables with non-zero partial for v
307 size_t i_x = *itr;
308 while( i_x < np1 )
309 { // N(i_x) = N(i_x) union J(i_v)
310 for_sparsity.binary_union(i_x, i_x, i_v + np1, for_sparsity);
311 i_x = *(++itr);
312 }
313 return;
314 }
315
316 } } } // END_CPPAD_LOCAL_SPARSE_NAMESPACE
317 # endif
318