1 | /* $Id: atomic_base.hpp 2900 2013-09-18 19:30:31Z bradbell $ */ |
---|
2 | # ifndef CPPAD_ATOMIC_BASE_INCLUDED |
---|
3 | # define CPPAD_ATOMIC_BASE_INCLUDED |
---|
4 | |
---|
5 | /* -------------------------------------------------------------------------- |
---|
6 | CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-13 Bradley M. Bell |
---|
7 | |
---|
8 | CppAD is distributed under multiple licenses. This distribution is under |
---|
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. |
---|
13 | Please visit http://www.coin-or.org/CppAD/ for information on other licenses. |
---|
14 | -------------------------------------------------------------------------- */ |
---|
15 | |
---|
16 | # include <set> |
---|
17 | # include <cppad/local/cppad_assert.hpp> |
---|
18 | // needed before one can use CPPAD_ASSERT_FIRST_CALL_NOT_PARALLEL |
---|
19 | # include <cppad/thread_alloc.hpp> |
---|
20 | |
---|
21 | CPPAD_BEGIN_NAMESPACE |
---|
22 | /*! |
---|
23 | \defgroup atomic_base.hpp atomic_base.hpp |
---|
24 | \{ |
---|
25 | \file atomic_base.hpp |
---|
26 | Base class for atomic user operations. |
---|
27 | */ |
---|
28 | |
---|
29 | template <class Base> |
---|
30 | class atomic_base { |
---|
31 | // =================================================================== |
---|
32 | public: |
---|
33 | enum option_enum { bool_sparsity_enum, set_sparsity_enum}; |
---|
34 | private: |
---|
35 | // ------------------------------------------------------ |
---|
36 | // constants |
---|
37 | /// name for this atomic funciton (used for error reporting) |
---|
38 | const std::string afun_name_; |
---|
39 | |
---|
40 | /// index of this object in class_object |
---|
41 | const size_t index_; |
---|
42 | |
---|
43 | // ----------------------------------------------------- |
---|
44 | // variables |
---|
45 | /// sparsity pattern this object is currently using |
---|
46 | /// (set by constructor and option member functions) |
---|
47 | option_enum sparsity_; |
---|
48 | |
---|
49 | /// temporary work space used afun, declared here to avoid memory |
---|
50 | /// allocation/deallocation for each call to afun |
---|
51 | vector<bool> afun_vx_[CPPAD_MAX_NUM_THREADS]; |
---|
52 | vector<bool> afun_vy_[CPPAD_MAX_NUM_THREADS]; |
---|
53 | vector<Base> afun_tx_[CPPAD_MAX_NUM_THREADS]; |
---|
54 | vector<Base> afun_ty_[CPPAD_MAX_NUM_THREADS]; |
---|
55 | |
---|
56 | // ----------------------------------------------------- |
---|
57 | // static member functions |
---|
58 | // 2DO: seems like this function needs to be public to be used by |
---|
59 | // *_sweep.hpp, but it seems to work this way. Why ? |
---|
60 | /// List of all the object in this class |
---|
61 | /// (null pointer used for objects that have been deleted) |
---|
62 | static std::vector<atomic_base *>& class_object(void) |
---|
63 | { CPPAD_ASSERT_FIRST_CALL_NOT_PARALLEL; |
---|
64 | static std::vector<atomic_base *> list_; |
---|
65 | return list_; |
---|
66 | } |
---|
67 | // ===================================================================== |
---|
68 | public: |
---|
69 | // ----------------------------------------------------- |
---|
70 | // member functions not in user API |
---|
71 | /// current sparsity setting |
---|
72 | option_enum sparsity(void) const |
---|
73 | { return sparsity_; } |
---|
74 | |
---|
75 | /// Name corresponding to a base_atomic object |
---|
76 | const std::string& afun_name(void) const |
---|
77 | { return afun_name_; } |
---|
78 | /* |
---|
79 | $begin atomic_ctor$$ |
---|
80 | $spell |
---|
81 | std |
---|
82 | afun |
---|
83 | arg |
---|
84 | CppAD |
---|
85 | bool |
---|
86 | ctor |
---|
87 | const |
---|
88 | $$ |
---|
89 | |
---|
90 | $section Atomic Function Constructor$$ |
---|
91 | $index constructor, atomic function$$ |
---|
92 | $index atomic, function constructor$$ |
---|
93 | $index function, atomic constructor$$ |
---|
94 | |
---|
95 | $head Syntax$$ |
---|
96 | $icode%atomic_user afun%(%ctor_arg_list%) |
---|
97 | %$$ |
---|
98 | $codei%atomic_base<%Base%>(%name%) |
---|
99 | %$$ |
---|
100 | |
---|
101 | $head atomic_user$$ |
---|
102 | |
---|
103 | $subhead ctor_arg_list$$ |
---|
104 | Is a list of arguments for the $icode atomic_user$$ constructor. |
---|
105 | |
---|
106 | $subhead afun$$ |
---|
107 | The object $icode afun$$ must stay in scope for as long |
---|
108 | as the corresponding atomic function is used. |
---|
109 | This includes use by any $cref/ADFun<Base>/ADFun/$$ that |
---|
110 | has this $icode atomic_user$$ operation in its |
---|
111 | $cref/operation sequence/glossary/Operation/Sequence/$$. |
---|
112 | |
---|
113 | $subhead Implementation$$ |
---|
114 | The user defined $icode atomic_user$$ class is a publicly derived class of |
---|
115 | $codei%atomic_base<%Base%>%$$. |
---|
116 | It should be declared as follows: |
---|
117 | $codei% |
---|
118 | class %atomic_user% : public CppAD::atomic_base<%Base%> { |
---|
119 | public: |
---|
120 | %atomic_user%(%ctor_arg_list%) : atomic_base<%Base%>(%name%) |
---|
121 | %...% |
---|
122 | }; |
---|
123 | %$$ |
---|
124 | where $icode ...$$ |
---|
125 | denotes the rest of the implementation of the derived class. |
---|
126 | This includes completing the constructor and |
---|
127 | all the virtual functions that have their |
---|
128 | $code atomic_base$$ implementations replaced by |
---|
129 | $icode atomic_user$$ implementations. |
---|
130 | |
---|
131 | $head atomic_base$$ |
---|
132 | |
---|
133 | $subhead Restrictions$$ |
---|
134 | The $code atomic_base$$ constructor cannot be called in |
---|
135 | $cref/parallel/ta_in_parallel/$$ mode. |
---|
136 | |
---|
137 | $subhead Base$$ |
---|
138 | The template parameter determines the |
---|
139 | $icode Base$$ type for this $codei%AD<%Base%>%$$ atomic operation. |
---|
140 | |
---|
141 | $subhead name$$ |
---|
142 | This $icode atomic_base$$ constructor argument has either of the |
---|
143 | following prototypes |
---|
144 | $codei% |
---|
145 | const char* %name% |
---|
146 | const std::string& %name% |
---|
147 | %$$ |
---|
148 | It is the name for this atomic function and is used for error reporting. |
---|
149 | The suggested value for $icode name$$ is $icode afun$$ or $icode atomic_user$$, |
---|
150 | i.e., the name of the corresponding atomic object or class. |
---|
151 | |
---|
152 | $end |
---|
153 | */ |
---|
154 | /*! |
---|
155 | Base class for atomic_user functions. |
---|
156 | |
---|
157 | \tparam Base |
---|
158 | This class is used for defining an AD<Base> atomic operation y = f(x). |
---|
159 | */ |
---|
160 | /// make sure user does not invoke the default constructor |
---|
161 | atomic_base(void) |
---|
162 | { CPPAD_ASSERT_KNOWN(false, |
---|
163 | "Attempt to use the atomic_base default constructor" |
---|
164 | ); |
---|
165 | } |
---|
166 | /*! |
---|
167 | Constructor |
---|
168 | |
---|
169 | \param name |
---|
170 | name used for error reporting |
---|
171 | */ |
---|
172 | atomic_base( const std::string& name) : |
---|
173 | afun_name_( name ) , |
---|
174 | index_( class_object().size() ) , |
---|
175 | sparsity_( set_sparsity_enum ) |
---|
176 | { CPPAD_ASSERT_KNOWN( |
---|
177 | ! thread_alloc::in_parallel() , |
---|
178 | "atomic_base: constructor cannot be called in parallel mode." |
---|
179 | ); |
---|
180 | class_object().push_back(this); |
---|
181 | } |
---|
182 | /// destructor informs CppAD that this atomic function with this index |
---|
183 | /// has dropped out of scope by setting its pointer to null |
---|
184 | virtual ~atomic_base(void) |
---|
185 | { CPPAD_ASSERT_UNKNOWN( class_object().size() > index_ ); |
---|
186 | class_object()[index_] = CPPAD_NULL; |
---|
187 | } |
---|
188 | /// atomic_base function object corresponding to a certain index |
---|
189 | static atomic_base* class_object(size_t index) |
---|
190 | { CPPAD_ASSERT_UNKNOWN( class_object().size() > index ); |
---|
191 | return class_object()[index]; |
---|
192 | } |
---|
193 | /* |
---|
194 | $begin atomic_option$$ |
---|
195 | $spell |
---|
196 | enum |
---|
197 | afun |
---|
198 | bool |
---|
199 | CppAD |
---|
200 | std |
---|
201 | typedef |
---|
202 | $$ |
---|
203 | |
---|
204 | $section Set Atomic Function Options$$ |
---|
205 | $index atomic, options$$ |
---|
206 | $index options, atomic$$ |
---|
207 | |
---|
208 | $head Syntax$$ |
---|
209 | $icode%afun%.option(%option_value%)%$$ |
---|
210 | |
---|
211 | $head atomic_sparsity$$ |
---|
212 | $index atomic_sparsity$$ |
---|
213 | $index sparsity, atomic$$ |
---|
214 | You can used this option to set to type used for |
---|
215 | $icode afun$$ sparsity patterns. |
---|
216 | This does not apply individual calls to $icode afun$$, |
---|
217 | but rather all its uses between when the sparsity pattern is set and when |
---|
218 | it is changed. |
---|
219 | If neither the $code set_sparsity_enum$$ or |
---|
220 | $code bool_sparsity_enum$$ option is set, |
---|
221 | the type for $icode atomic_sparsity$$ is one of the two choices below |
---|
222 | (and otherwise unspecified). |
---|
223 | |
---|
224 | $subhead bool_sparsity_enum$$ |
---|
225 | $index bool_sparsity_enum$$ |
---|
226 | If $icode option_value$$ is $code atomic_base::bool_sparsity_enum$$, |
---|
227 | then the type used by $icode afun$$ for |
---|
228 | $cref/sparsity patterns/glossary/Sparsity Pattern/$$, |
---|
229 | (after the option is set) will be |
---|
230 | $codei% |
---|
231 | typedef CppAD::vector<bool> %atomic_sparsity% |
---|
232 | %$$ |
---|
233 | If $icode r$$ is a sparsity pattern |
---|
234 | for a matrix $latex R \in B^{p \times q}$$: |
---|
235 | $icode%r%.size() == %p% * %q%$$. |
---|
236 | |
---|
237 | $subhead set_sparsity_enum$$ |
---|
238 | $index set_sparsity_enum$$ |
---|
239 | If $icode option_value$$ is $code atomic_base::set_sparsity_enum$$, |
---|
240 | then the type used by $icode afun$$ for |
---|
241 | $cref/sparsity patterns/glossary/Sparsity Pattern/$$, |
---|
242 | (after the option is set) will be |
---|
243 | $codei% |
---|
244 | typedef CppAD::vector< std::set<size_t> > %atomic_sparsity% |
---|
245 | %$$ |
---|
246 | If $icode r$$ is a sparsity pattern |
---|
247 | for a matrix $code R \in B^{p \times q}$$: |
---|
248 | $icode%r%.size() == %p%$$, and for $latex i = 0 , \ldots , p-1$$, |
---|
249 | the elements of $icode%r%[%i%]%$$ are between zero and $latex q-1$$ inclusive. |
---|
250 | |
---|
251 | $end |
---|
252 | */ |
---|
253 | void option(enum option_enum option_value) |
---|
254 | { switch( option_value ) |
---|
255 | { case bool_sparsity_enum: |
---|
256 | case set_sparsity_enum: |
---|
257 | sparsity_ = option_value; |
---|
258 | break; |
---|
259 | |
---|
260 | default: |
---|
261 | CPPAD_ASSERT_KNOWN( |
---|
262 | false, |
---|
263 | "atoic_base::option: option_value is not valid" |
---|
264 | ); |
---|
265 | } |
---|
266 | return; |
---|
267 | } |
---|
268 | /* |
---|
269 | ----------------------------------------------------------------------------- |
---|
270 | $begin atomic_afun$$ |
---|
271 | |
---|
272 | $spell |
---|
273 | afun |
---|
274 | const |
---|
275 | CppAD |
---|
276 | $$ |
---|
277 | |
---|
278 | $section Using an Atomic Function$$ |
---|
279 | $index atomic, use function$$ |
---|
280 | |
---|
281 | $head Syntax$$ |
---|
282 | $icode%afun%(%ax%, %ay%)%$$ |
---|
283 | |
---|
284 | $head Purpose$$ |
---|
285 | Given $icode ax$$, |
---|
286 | this call computes the corresponding value of $icode ay$$. |
---|
287 | If $codei%AD<%Base%>%$$ operations are being recorded, |
---|
288 | it enters the computation as an atomic operation in the recording; |
---|
289 | see $cref/start recording/Independent/Start Recording/$$. |
---|
290 | |
---|
291 | $head ADVector$$ |
---|
292 | The type $icode ADVector$$ must be a |
---|
293 | $cref/simple vector class/SimpleVector/$$ with elements of type |
---|
294 | $codei%AD<%Base%>%$$; see $cref/Base/atomic_ctor/atomic_base/Base/$$. |
---|
295 | |
---|
296 | $head afun$$ |
---|
297 | is a $cref/atomic_user/atomic_ctor/atomic_user/$$ object |
---|
298 | and this $icode afun$$ function call is implemented by the |
---|
299 | $cref/atomic_base/atomic_ctor/atomic_base/$$ class. |
---|
300 | |
---|
301 | $head ax$$ |
---|
302 | This argument has prototype |
---|
303 | $codei% |
---|
304 | const %ADVector%& %ax% |
---|
305 | %$$ |
---|
306 | and size must be equal to $icode n$$. |
---|
307 | It specifies vector $latex x \in B^n$$ |
---|
308 | at which an $codei%AD<%Base%>%$$ version of |
---|
309 | $latex y = f(x)$$ is to be evaluated; see |
---|
310 | $cref/Base/atomic_ctor/atomic_base/Base/$$. |
---|
311 | |
---|
312 | $head ay$$ |
---|
313 | This argument has prototype |
---|
314 | $codei% |
---|
315 | %ADVector%& %ay% |
---|
316 | %$$ |
---|
317 | and size must be equal to $icode m$$. |
---|
318 | The input values of its elements |
---|
319 | are not specified (must not matter). |
---|
320 | Upon return, it is an $codei%AD<%Base%>%$$ version of |
---|
321 | $latex y = f(x)$$. |
---|
322 | |
---|
323 | $end |
---|
324 | ----------------------------------------------------------------------------- |
---|
325 | */ |
---|
326 | /*! |
---|
327 | Implement the user call to <tt>afun(ax, ay)</tt> and old_atomic call to |
---|
328 | <tt>afun(ax, ay, id)</tt>. |
---|
329 | |
---|
330 | \tparam ADVector |
---|
331 | A simple vector class with elements of type <code>AD<Base></code>. |
---|
332 | |
---|
333 | \param id |
---|
334 | optional extra information vector that is just passed through by CppAD, |
---|
335 | and used by old_atomic derived class (not other derived classes). |
---|
336 | This is an extra parameter to the virtual callbacks for old_atomic; |
---|
337 | see the set_id member function. |
---|
338 | |
---|
339 | \param ax |
---|
340 | is the argument vector for this call, |
---|
341 | <tt>ax.size()</tt> determines the number of arguments. |
---|
342 | |
---|
343 | \param ay |
---|
344 | is the result vector for this call, |
---|
345 | <tt>ay.size()</tt> determines the number of results. |
---|
346 | */ |
---|
347 | template <class ADVector> |
---|
348 | void operator()( |
---|
349 | const ADVector& ax , |
---|
350 | ADVector& ay , |
---|
351 | size_t id = 0 ) |
---|
352 | { size_t i, j; |
---|
353 | size_t n = ax.size(); |
---|
354 | size_t m = ay.size(); |
---|
355 | # ifndef NDEBUG |
---|
356 | bool ok; |
---|
357 | std::string msg = "atomic_base: " + afun_name_ + ".eval: "; |
---|
358 | if( (n == 0) | (m == 0) ) |
---|
359 | { msg += "ax.size() or ay.size() is zero"; |
---|
360 | CPPAD_ASSERT_KNOWN(false, msg.c_str() ); |
---|
361 | } |
---|
362 | # endif |
---|
363 | size_t thread = thread_alloc::thread_num(); |
---|
364 | vector <Base>& tx = afun_tx_[thread]; |
---|
365 | vector <Base>& ty = afun_ty_[thread]; |
---|
366 | vector <bool>& vx = afun_vx_[thread]; |
---|
367 | vector <bool>& vy = afun_vy_[thread]; |
---|
368 | // |
---|
369 | if( vx.size() != n ) |
---|
370 | { vx.resize(n); |
---|
371 | tx.resize(n); |
---|
372 | } |
---|
373 | if( vy.size() != m ) |
---|
374 | { vy.resize(m); |
---|
375 | ty.resize(m); |
---|
376 | } |
---|
377 | // |
---|
378 | // Determine tape corresponding to variables in ax |
---|
379 | tape_id_t tape_id = 0; |
---|
380 | ADTape<Base>* tape = CPPAD_NULL; |
---|
381 | for(j = 0; j < n; j++) |
---|
382 | { tx[j] = ax[j].value_; |
---|
383 | vx[j] = Variable( ax[j] ); |
---|
384 | if( vx[j] ) |
---|
385 | { |
---|
386 | if( tape_id == 0 ) |
---|
387 | { tape = ax[j].tape_this(); |
---|
388 | tape_id = ax[j].tape_id_; |
---|
389 | CPPAD_ASSERT_UNKNOWN( tape != CPPAD_NULL ); |
---|
390 | } |
---|
391 | # ifndef NDEBUG |
---|
392 | if( tape_id != ax[j].tape_id_ ) |
---|
393 | { msg += afun_name_ + |
---|
394 | ": ax contains variables from different threads."; |
---|
395 | CPPAD_ASSERT_KNOWN(false, msg.c_str()); |
---|
396 | } |
---|
397 | # endif |
---|
398 | } |
---|
399 | } |
---|
400 | // Use zero order forward mode to compute values |
---|
401 | size_t q = 0, p = 0; |
---|
402 | set_id(id); |
---|
403 | # ifdef NDEBUG |
---|
404 | forward(q, p, vx, vy, tx, ty); |
---|
405 | # else |
---|
406 | ok = forward(q, p, vx, vy, tx, ty); |
---|
407 | if( ! ok ) |
---|
408 | { msg += afun_name_ + ": ok is false for " |
---|
409 | "zero order forward mode calculation."; |
---|
410 | CPPAD_ASSERT_KNOWN(false, msg.c_str()); |
---|
411 | } |
---|
412 | # endif |
---|
413 | bool record_operation = false; |
---|
414 | for(i = 0; i < m; i++) |
---|
415 | { |
---|
416 | // pass back values |
---|
417 | ay[i].value_ = ty[i]; |
---|
418 | |
---|
419 | // initialize entire vector parameters (not in tape) |
---|
420 | ay[i].tape_id_ = 0; |
---|
421 | ay[i].taddr_ = 0; |
---|
422 | |
---|
423 | // we need to record this operation if |
---|
424 | // any of the eleemnts of ay are variables, |
---|
425 | record_operation |= vy[i]; |
---|
426 | } |
---|
427 | # ifndef NDEBUG |
---|
428 | if( record_operation & (tape == CPPAD_NULL) ) |
---|
429 | { msg += |
---|
430 | "all elements of vx are false but vy contains a true element"; |
---|
431 | CPPAD_ASSERT_KNOWN(false, msg.c_str() ); |
---|
432 | } |
---|
433 | # endif |
---|
434 | // if tape is not null, ay is on the tape |
---|
435 | if( record_operation ) |
---|
436 | { |
---|
437 | // Operator that marks beginning of this atomic operation |
---|
438 | CPPAD_ASSERT_UNKNOWN( NumRes(UserOp) == 0 ); |
---|
439 | CPPAD_ASSERT_UNKNOWN( NumArg(UserOp) == 4 ); |
---|
440 | tape->Rec_.PutArg(index_, id, n, m); |
---|
441 | tape->Rec_.PutOp(UserOp); |
---|
442 | |
---|
443 | // Now put n operators, one for each element of arugment vector |
---|
444 | CPPAD_ASSERT_UNKNOWN( NumRes(UsravOp) == 0 ); |
---|
445 | CPPAD_ASSERT_UNKNOWN( NumRes(UsrapOp) == 0 ); |
---|
446 | CPPAD_ASSERT_UNKNOWN( NumArg(UsravOp) == 1 ); |
---|
447 | CPPAD_ASSERT_UNKNOWN( NumArg(UsrapOp) == 1 ); |
---|
448 | for(j = 0; j < n; j++) |
---|
449 | { if( vx[j] ) |
---|
450 | { // information for an argument that is a variable |
---|
451 | tape->Rec_.PutArg(ax[j].taddr_); |
---|
452 | tape->Rec_.PutOp(UsravOp); |
---|
453 | } |
---|
454 | else |
---|
455 | { // information for an arugment that is parameter |
---|
456 | addr_t par = tape->Rec_.PutPar(ax[j].value_); |
---|
457 | tape->Rec_.PutArg(par); |
---|
458 | tape->Rec_.PutOp(UsrapOp); |
---|
459 | } |
---|
460 | } |
---|
461 | |
---|
462 | // Now put m operators, one for each element of result vector |
---|
463 | CPPAD_ASSERT_UNKNOWN( NumArg(UsrrpOp) == 1 ); |
---|
464 | CPPAD_ASSERT_UNKNOWN( NumRes(UsrrpOp) == 0 ); |
---|
465 | CPPAD_ASSERT_UNKNOWN( NumArg(UsrrvOp) == 0 ); |
---|
466 | CPPAD_ASSERT_UNKNOWN( NumRes(UsrrvOp) == 1 ); |
---|
467 | for(i = 0; i < m; i++) |
---|
468 | { if( vy[i] ) |
---|
469 | { ay[i].taddr_ = tape->Rec_.PutOp(UsrrvOp); |
---|
470 | ay[i].tape_id_ = tape_id; |
---|
471 | } |
---|
472 | else |
---|
473 | { addr_t par = tape->Rec_.PutPar(ay[i].value_); |
---|
474 | tape->Rec_.PutArg(par); |
---|
475 | tape->Rec_.PutOp(UsrrpOp); |
---|
476 | } |
---|
477 | } |
---|
478 | |
---|
479 | // Put a duplicate UserOp at end of UserOp sequence |
---|
480 | tape->Rec_.PutArg(index_, id, n, m); |
---|
481 | tape->Rec_.PutOp(UserOp); |
---|
482 | } |
---|
483 | return; |
---|
484 | } |
---|
485 | /* |
---|
486 | ----------------------------------------------------------------------------- |
---|
487 | $begin atomic_forward$$ |
---|
488 | $spell |
---|
489 | hes |
---|
490 | afun |
---|
491 | vx |
---|
492 | vy |
---|
493 | ty |
---|
494 | Taylor |
---|
495 | const |
---|
496 | CppAD |
---|
497 | bool |
---|
498 | $$ |
---|
499 | |
---|
500 | $section Atomic Forward Mode$$ |
---|
501 | $index atomic, forward callback$$ |
---|
502 | $index forward, atomic callback$$ |
---|
503 | $index forward, atomic virtual$$ |
---|
504 | |
---|
505 | |
---|
506 | $head Syntax$$ |
---|
507 | $icode%ok% = %afun%.forward(%q%, %p%, %vx%, %vy%, %tx%, %ty%)%$$ |
---|
508 | |
---|
509 | $head Purpose$$ |
---|
510 | This virtual function is used by $cref atomic_afun$$ |
---|
511 | to evaluate function values. |
---|
512 | It is also used buy $cref/forward/Forward/$$ |
---|
513 | to compute function vales and derivatives. |
---|
514 | |
---|
515 | $head Implementation$$ |
---|
516 | This virtual function must be defined by the |
---|
517 | $cref/atomic_user/atomic_ctor/atomic_user/$$ class. |
---|
518 | It can just return $icode%ok% == false%$$ |
---|
519 | (and not compute anything) for values |
---|
520 | of $icode%p% > 0%$$ that are greater than those used by your |
---|
521 | $cref/forward/Forward/$$ mode calculations. |
---|
522 | |
---|
523 | $head q$$ |
---|
524 | The argument $icode q$$ has prototype |
---|
525 | $codei% |
---|
526 | size_t %q% |
---|
527 | %$$ |
---|
528 | It specifies the lowest order Taylor coefficient that we are evaluating. |
---|
529 | During calls to $cref atomic_afun$$, $icode%q% == 0%$$. |
---|
530 | |
---|
531 | $head p$$ |
---|
532 | The argument $icode p$$ has prototype |
---|
533 | $codei% |
---|
534 | size_t %p% |
---|
535 | %$$ |
---|
536 | It specifies the highest order Taylor coefficient that we are evaluating. |
---|
537 | During calls to $cref atomic_afun$$, $icode%p% == 0%$$. |
---|
538 | |
---|
539 | $head vx$$ |
---|
540 | The $code forward$$ argument $icode vx$$ has prototype |
---|
541 | $codei% |
---|
542 | const CppAD::vector<bool>& %vx% |
---|
543 | %$$ |
---|
544 | The case $icode%vx%.size() > 0%$$ only occurs while evaluating a call to |
---|
545 | $cref atomic_afun$$. |
---|
546 | In this case, |
---|
547 | $icode%q% == %p% == 0%$$, |
---|
548 | $icode%vx%.size() == %n%$$, and |
---|
549 | for $latex j = 0 , \ldots , n-1$$, |
---|
550 | $icode%vx%[%j%]%$$ is true if and only if |
---|
551 | $icode%ax%[%j%]%$$ is a $cref/variable/glossary/Variable/$$ |
---|
552 | in the corresponding call to |
---|
553 | $codei% |
---|
554 | %afun%(%ax%, %ay%, %id%) |
---|
555 | %$$ |
---|
556 | This is only necessary to use for second order cross terms because |
---|
557 | CppAD knows which components are variables; e.g., see |
---|
558 | $code rev_sparse_hes$$ in the $cref atomic_matrix_mul.hpp$$ example. |
---|
559 | $pre |
---|
560 | |
---|
561 | $$ |
---|
562 | If $icode%vx%.size() == 0%$$, |
---|
563 | then $icode%vy%.size() == 0%$$ and neither of these vectors |
---|
564 | should be used. |
---|
565 | |
---|
566 | $head vy$$ |
---|
567 | The $code forward$$ argument $icode vy$$ has prototype |
---|
568 | $codei% |
---|
569 | CppAD::vector<bool>& %vy% |
---|
570 | %$$ |
---|
571 | If $icode%vy%.size() == 0%$$, it should not be used. |
---|
572 | Otherwise, |
---|
573 | $icode%p% == 0%$$ and $icode%vy%.size() == %m%$$. |
---|
574 | The input values of the elements of $icode vy$$ |
---|
575 | are not specified (must not matter). |
---|
576 | Upon return, for $latex j = 0 , \ldots , m-1$$, |
---|
577 | $icode%vy%[%i%]%$$ is true if and only if |
---|
578 | $icode%ay%[%i%]%$$ is a variable |
---|
579 | (CppAD uses $icode vy$$ to reduce the necessary computations). |
---|
580 | |
---|
581 | $head tx$$ |
---|
582 | The argument $icode tx$$ has prototype |
---|
583 | $codei% |
---|
584 | const CppAD::vector<%Base%>& %tx% |
---|
585 | %$$ |
---|
586 | and $icode%tx%.size() == (%p%+1)*%n%$$. |
---|
587 | For $latex j = 0 , \ldots , n-1$$ and $latex k = 0 , \ldots , p$$, |
---|
588 | we use the Taylor coefficient notation |
---|
589 | $latex \[ |
---|
590 | \begin{array}{rcl} |
---|
591 | x_j^k & = & tx [ j * ( p + 1 ) + k ] |
---|
592 | \\ |
---|
593 | X_j (t) & = & x_j^0 + x_j^1 t^1 + \cdots + x_j^p t^p |
---|
594 | \end{array} |
---|
595 | \] $$ |
---|
596 | Note that superscripts represent an index for $latex x_j^k$$ |
---|
597 | and an exponent for $latex t^k$$. |
---|
598 | Also note that the Taylor coefficients for $latex X(t)$$ correspond |
---|
599 | to the derivatives of $latex X(t)$$ at $latex t = 0$$ in the following way: |
---|
600 | $latex \[ |
---|
601 | x_j^k = \frac{1}{ k ! } X_j^{(k)} (0) |
---|
602 | \] $$ |
---|
603 | |
---|
604 | $head ty$$ |
---|
605 | The argument $icode ty$$ has prototype |
---|
606 | $codei% |
---|
607 | CppAD::vector<%Base%>& %ty% |
---|
608 | %$$ |
---|
609 | and $icode%tx%.size() == (%p%+1)*%m%$$. |
---|
610 | Upon return, |
---|
611 | For $latex i = 0 , \ldots , m-1$$ and $latex k = 0 , \ldots , p$$, |
---|
612 | $latex \[ |
---|
613 | \begin{array}{rcl} |
---|
614 | Y_i (t) & = & f_i [ X(t) ] |
---|
615 | \\ |
---|
616 | Y_i (t) & = & y_i^0 + y_i^1 t^1 + \cdots + y_i^p t^p + o ( t^p ) |
---|
617 | \\ |
---|
618 | ty [ i * ( p + 1 ) + k ] & = & y_i^k |
---|
619 | \end{array} |
---|
620 | \] $$ |
---|
621 | where $latex o( t^p ) / t^p \rightarrow 0$$ as $latex t \rightarrow 0$$. |
---|
622 | Note that superscripts represent an index for $latex y_j^k$$ |
---|
623 | and an exponent for $latex t^k$$. |
---|
624 | Also note that the Taylor coefficients for $latex Y(t)$$ correspond |
---|
625 | to the derivatives of $latex Y(t)$$ at $latex t = 0$$ in the following way: |
---|
626 | $latex \[ |
---|
627 | y_j^k = \frac{1}{ k ! } Y_j^{(k)} (0) |
---|
628 | \] $$ |
---|
629 | If $latex q > 0$$, |
---|
630 | for $latex i = 0 , \ldots , m-1$$ and $latex k = 0 , \ldots , q-1$$, |
---|
631 | the input of $icode ty$$ satisfies |
---|
632 | $latex \[ |
---|
633 | ty [ i * ( p + 1 ) + k ] = y_i^k |
---|
634 | \]$$ |
---|
635 | and hence the corresponding elements need not be recalculated. |
---|
636 | |
---|
637 | $head ok$$ |
---|
638 | If the required results are calculated, $icode ok$$ should be true. |
---|
639 | Otherwise, it should be false. |
---|
640 | |
---|
641 | $head Example$$ |
---|
642 | For example, suppose that $icode%p% == 2%$$, |
---|
643 | and you know how to compute the function $latex f(x)$$, |
---|
644 | its first derivative $latex f^{(1)} (x)$$, |
---|
645 | and it component wise Hessian $latex f_i^{(2)} (x)$$. |
---|
646 | Then you can compute $icode ty$$ using the following formulas: |
---|
647 | $latex \[ |
---|
648 | \begin{array}{rcl} |
---|
649 | y_i^0 & = & Y(0) |
---|
650 | = f_i ( x^0 ) |
---|
651 | \\ |
---|
652 | y_i^1 & = & Y^{(1)} ( 0 ) |
---|
653 | = f_i^{(1)} ( x^0 ) X^{(1)} ( 0 ) |
---|
654 | = f_i^{(1)} ( x^0 ) x^1 |
---|
655 | \\ |
---|
656 | y_i^2 |
---|
657 | & = & \frac{1}{2 !} Y^{(2)} (0) |
---|
658 | \\ |
---|
659 | & = & \frac{1}{2} X^{(1)} (0)^\R{T} f_i^{(2)} ( x^0 ) X^{(1)} ( 0 ) |
---|
660 | + \frac{1}{2} f_i^{(1)} ( x^0 ) X^{(2)} ( 0 ) |
---|
661 | \\ |
---|
662 | & = & \frac{1}{2} (x^1)^\R{T} f_i^{(2)} ( x^0 ) x^1 |
---|
663 | + f_i^{(1)} ( x^0 ) x^2 |
---|
664 | \end{array} |
---|
665 | \] $$ |
---|
666 | For $latex i = 0 , \ldots , m-1$$, and $latex k = 0 , 1 , 2$$, |
---|
667 | $latex \[ |
---|
668 | ty [ i * (p + 1) + k ] = y_i^k |
---|
669 | \] $$ |
---|
670 | |
---|
671 | $end |
---|
672 | ----------------------------------------------------------------------------- |
---|
673 | */ |
---|
674 | /*! |
---|
675 | Link from atomic_base to forward mode |
---|
676 | |
---|
677 | \param q [in] |
---|
678 | lowerest order for this forward mode calculation. |
---|
679 | |
---|
680 | \param p [in] |
---|
681 | highest order for this forward mode calculation. |
---|
682 | |
---|
683 | \param vx [in] |
---|
684 | if size not zero, which components of \c x are variables |
---|
685 | |
---|
686 | \param vy [out] |
---|
687 | if size not zero, which components of \c y are variables |
---|
688 | |
---|
689 | \param tx [in] |
---|
690 | Taylor coefficients corresponding to \c x for this calculation. |
---|
691 | |
---|
692 | \param ty [out] |
---|
693 | Taylor coefficient corresponding to \c y for this calculation |
---|
694 | |
---|
695 | See the forward mode in user's documentation for base_atomic |
---|
696 | */ |
---|
697 | virtual bool forward( |
---|
698 | size_t q , |
---|
699 | size_t p , |
---|
700 | const vector<bool>& vx , |
---|
701 | vector<bool>& vy , |
---|
702 | const vector<Base>& tx , |
---|
703 | vector<Base>& ty ) |
---|
704 | { return false; } |
---|
705 | /* |
---|
706 | ----------------------------------------------------------------------------- |
---|
707 | $begin atomic_reverse$$ |
---|
708 | $spell |
---|
709 | afun |
---|
710 | ty |
---|
711 | px |
---|
712 | py |
---|
713 | Taylor |
---|
714 | const |
---|
715 | CppAD |
---|
716 | $$ |
---|
717 | |
---|
718 | $section Atomic Reverse Mode$$ |
---|
719 | $index atomic, reverse callback$$ |
---|
720 | $index reverse, atomic callback$$ |
---|
721 | $index reverse, atomic virtual$$ |
---|
722 | $spell |
---|
723 | bool |
---|
724 | $$ |
---|
725 | |
---|
726 | $head Syntax$$ |
---|
727 | $icode%ok% = %afun%.reverse(%p%, %tx%, %ty%, %px%, %py%)%$$ |
---|
728 | |
---|
729 | $head Purpose$$ |
---|
730 | This function is used by $cref/reverse/Reverse/$$ |
---|
731 | to compute derivatives. |
---|
732 | |
---|
733 | $head Implementation$$ |
---|
734 | If you are using |
---|
735 | $cref/reverse/Reverse/$$ mode, |
---|
736 | this virtual function must be defined by the |
---|
737 | $cref/atomic_user/atomic_ctor/atomic_user/$$ class. |
---|
738 | It can just return $icode%ok% == false%$$ |
---|
739 | (and not compute anything) for values |
---|
740 | of $icode p$$ that are greater than those used by your |
---|
741 | $cref/reverse/Reverse/$$ mode calculations. |
---|
742 | |
---|
743 | $head p$$ |
---|
744 | The argument $icode p$$ has prototype |
---|
745 | $codei% |
---|
746 | size_t %p% |
---|
747 | %$$ |
---|
748 | It specifies the highest order Taylor coefficient that |
---|
749 | computing the derivative of. |
---|
750 | |
---|
751 | $head tx$$ |
---|
752 | The argument $icode tx$$ has prototype |
---|
753 | $codei% |
---|
754 | const CppAD::vector<%Base%>& %tx% |
---|
755 | %$$ |
---|
756 | and $icode%tx%.size() == (%p%+1)*%n%$$. |
---|
757 | For $latex j = 0 , \ldots , n-1$$ and $latex k = 0 , \ldots , p$$, |
---|
758 | we use the Taylor coefficient notation |
---|
759 | $latex \[ |
---|
760 | \begin{array}{rcl} |
---|
761 | x_j^k & = & tx [ j * ( p + 1 ) + k ] |
---|
762 | \\ |
---|
763 | X_j (t) & = & x_j^0 + x_j^1 t^1 + \cdots + x_j^p t^p |
---|
764 | \end{array} |
---|
765 | \] $$ |
---|
766 | Note that superscripts represent an index for $latex x_j^k$$ |
---|
767 | and an exponent for $latex t^k$$. |
---|
768 | Also note that the Taylor coefficients for $latex X(t)$$ correspond |
---|
769 | to the derivatives of $latex X(t)$$ at $latex t = 0$$ in the following way: |
---|
770 | $latex \[ |
---|
771 | x_j^k = \frac{1}{ k ! } X_j^{(k)} (0) |
---|
772 | \] $$ |
---|
773 | |
---|
774 | $head ty$$ |
---|
775 | The argument $icode ty$$ has prototype |
---|
776 | $codei% |
---|
777 | const CppAD::vector<%Base%>& %ty% |
---|
778 | %$$ |
---|
779 | and $icode%tx%.size() == (%p%+1)*%m%$$. |
---|
780 | For $latex i = 0 , \ldots , m-1$$ and $latex k = 0 , \ldots , p$$, |
---|
781 | we use the Taylor coefficient notation |
---|
782 | $latex \[ |
---|
783 | \begin{array}{rcl} |
---|
784 | Y_i (t) & = & f_i [ X(t) ] |
---|
785 | \\ |
---|
786 | Y_i (t) & = & y_i^0 + y_i^1 t^1 + \cdots + y_i^p t^p + o ( t^p ) |
---|
787 | \\ |
---|
788 | y_i^k & = & ty [ i * ( p + 1 ) + k ] |
---|
789 | \end{array} |
---|
790 | \] $$ |
---|
791 | where $latex o( t^p ) / t^p \rightarrow 0$$ as $latex t \rightarrow 0$$. |
---|
792 | Note that superscripts represent an index for $latex y_j^k$$ |
---|
793 | and an exponent for $latex t^k$$. |
---|
794 | Also note that the Taylor coefficients for $latex Y(t)$$ correspond |
---|
795 | to the derivatives of $latex Y(t)$$ at $latex t = 0$$ in the following way: |
---|
796 | $latex \[ |
---|
797 | y_j^k = \frac{1}{ k ! } Y_j^{(k)} (0) |
---|
798 | \] $$ |
---|
799 | |
---|
800 | |
---|
801 | $head F, G, H$$ |
---|
802 | We use the notation $latex \{ x_j^k \} \in B^{n \times (p+1)}$$ for |
---|
803 | $latex \[ |
---|
804 | \{ x_j^k \W{:} j = 0 , \ldots , n-1, k = 0 , \ldots , p \} |
---|
805 | \]$$ |
---|
806 | We use the notation $latex \{ y_i^k \} \in B^{m \times (p+1)}$$ for |
---|
807 | $latex \[ |
---|
808 | \{ y_i^k \W{:} i = 0 , \ldots , m-1, k = 0 , \ldots , p \} |
---|
809 | \]$$ |
---|
810 | We define the function |
---|
811 | $latex F : B^{n \times (p+1)} \rightarrow B^{m \times (p+1)}$$ by |
---|
812 | $latex \[ |
---|
813 | y_i^k = F_i^k [ \{ x_j^k \} ] |
---|
814 | \] $$ |
---|
815 | We use $latex G : B^{m \times (p+1)} \rightarrow B$$ |
---|
816 | to denote an arbitrary scalar valued function of $latex \{ y_i^k \}$$. |
---|
817 | We use $latex H : B^{n \times (p+1)} \rightarrow B$$ |
---|
818 | defined by |
---|
819 | $latex \[ |
---|
820 | H ( \{ x_j^k \} ) = G[ F( \{ x_j^k \} ) ] |
---|
821 | \] $$ |
---|
822 | |
---|
823 | $head py$$ |
---|
824 | The argument $icode py$$ has prototype |
---|
825 | $codei% |
---|
826 | const CppAD::vector<%Base%>& %py% |
---|
827 | %$$ |
---|
828 | and $icode%py%.size() == m * (%p%+1)%$$. |
---|
829 | For $latex i = 0 , \ldots , m-1$$, $latex k = 0 , \ldots , p$$, |
---|
830 | $latex \[ |
---|
831 | py[ i * (p + 1 ) + k ] = \partial G / \partial y_i^k |
---|
832 | \] $$ |
---|
833 | |
---|
834 | $subhead px$$ |
---|
835 | The $icode px$$ has prototype |
---|
836 | $codei% |
---|
837 | CppAD::vector<%Base%>& %px% |
---|
838 | %$$ |
---|
839 | and $icode%px%.size() == n * (%p%+1)%$$. |
---|
840 | The input values of the elements of $icode px$$ |
---|
841 | are not specified (must not matter). |
---|
842 | Upon return, |
---|
843 | for $latex j = 0 , \ldots , n-1$$ and $latex \ell = 0 , \ldots , p$$, |
---|
844 | $latex \[ |
---|
845 | \begin{array}{rcl} |
---|
846 | px [ j * (p + 1) + \ell ] & = & \partial H / \partial x_j^\ell |
---|
847 | \\ |
---|
848 | & = & |
---|
849 | ( \partial G / \partial \{ y_i^k \} ) |
---|
850 | ( \partial \{ y_i^k \} / \partial x_j^\ell ) |
---|
851 | \\ |
---|
852 | & = & |
---|
853 | \sum_{i=0}^{m-1} \sum_{k=0}^p |
---|
854 | ( \partial G / \partial y_i^k ) ( \partial y_i^k / \partial x_j^\ell ) |
---|
855 | \\ |
---|
856 | & = & |
---|
857 | \sum_{i=0}^{m-1} \sum_{k=\ell}^p |
---|
858 | py[ i * (p + 1 ) + k ] ( \partial F_i^k / \partial x_j^\ell ) |
---|
859 | \end{array} |
---|
860 | \] $$ |
---|
861 | Note that we have used the fact that for $latex k < \ell$$, |
---|
862 | $latex \partial F_i^k / \partial x_j^\ell = 0$$. |
---|
863 | |
---|
864 | $head ok$$ |
---|
865 | The return value $icode ok$$ has prototype |
---|
866 | $codei% |
---|
867 | bool %ok% |
---|
868 | %$$ |
---|
869 | If it is $code true$$, the corresponding evaluation succeeded, |
---|
870 | otherwise it failed. |
---|
871 | |
---|
872 | $end |
---|
873 | ----------------------------------------------------------------------------- |
---|
874 | */ |
---|
875 | /*! |
---|
876 | Link from reverse mode sweep to users routine. |
---|
877 | |
---|
878 | \param p [in] |
---|
879 | highest order for this reverse mode calculation. |
---|
880 | |
---|
881 | \param tx [in] |
---|
882 | Taylor coefficients corresponding to \c x for this calculation. |
---|
883 | |
---|
884 | \param ty [in] |
---|
885 | Taylor coefficient corresponding to \c y for this calculation |
---|
886 | |
---|
887 | \param px [out] |
---|
888 | Partials w.r.t. the \c x Taylor coefficients. |
---|
889 | |
---|
890 | \param py [in] |
---|
891 | Partials w.r.t. the \c y Taylor coefficients. |
---|
892 | |
---|
893 | See atomic_reverse mode use documentation |
---|
894 | */ |
---|
895 | virtual bool reverse( |
---|
896 | size_t p , |
---|
897 | const vector<Base>& tx , |
---|
898 | const vector<Base>& ty , |
---|
899 | vector<Base>& px , |
---|
900 | const vector<Base>& py ) |
---|
901 | { return false; } |
---|
902 | /* |
---|
903 | -------------------------------------- --------------------------------------- |
---|
904 | $begin atomic_for_sparse_jac$$ |
---|
905 | $spell |
---|
906 | afun |
---|
907 | Jacobian |
---|
908 | jac |
---|
909 | const |
---|
910 | CppAD |
---|
911 | std |
---|
912 | bool |
---|
913 | std |
---|
914 | $$ |
---|
915 | |
---|
916 | $section Atomic Forward Jacobian Sparsity Patterns$$ |
---|
917 | $index atomic, for_sparse_jac callback$$ |
---|
918 | $index for_sparse_jac, atomic callback$$ |
---|
919 | $index for_sparse_jac, atomic virtual$$ |
---|
920 | |
---|
921 | $head Syntax$$ |
---|
922 | $icode%ok% = %afun%.for_sparse_jac(%q%, %r%, %s%)%$$ |
---|
923 | |
---|
924 | $head Purpose$$ |
---|
925 | This function is used by $cref ForSparseJac$$ to compute |
---|
926 | Jacobian sparsity patterns. |
---|
927 | For a fixed matrix $latex R \in B^{n \times q}$$, |
---|
928 | the Jacobian of $latex f( x + R * u)$$ with respect to $latex u \in B^q$$ is |
---|
929 | $latex \[ |
---|
930 | S(x) = f^{(1)} (x) * R |
---|
931 | \] $$ |
---|
932 | Given a $cref/sparsity pattern/glossary/Sparsity Pattern/$$ for $latex R$$, |
---|
933 | $code for_sparse_jac$$ computes a sparsity pattern for $latex S(x)$$. |
---|
934 | |
---|
935 | $head Implementation$$ |
---|
936 | If you are using $cref ForSparseJac$$, |
---|
937 | this virtual function must be defined by the |
---|
938 | $cref/atomic_user/atomic_ctor/atomic_user/$$ class. |
---|
939 | |
---|
940 | $subhead q$$ |
---|
941 | The argument $icode q$$ has prototype |
---|
942 | $codei% |
---|
943 | size_t %q% |
---|
944 | %$$ |
---|
945 | It specifies the number of columns in |
---|
946 | $latex R \in B^{n \times q}$$ and the Jacobian |
---|
947 | $latex S(x) \in B^{m \times q}$$. |
---|
948 | |
---|
949 | $subhead r$$ |
---|
950 | This argument has prototype |
---|
951 | $codei% |
---|
952 | const %atomic_sparsity%& %r% |
---|
953 | %$$ |
---|
954 | and is a $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
955 | $latex R \in B^{n \times q}$$. |
---|
956 | |
---|
957 | $subhead s$$ |
---|
958 | This argument has prototype |
---|
959 | $codei% |
---|
960 | %atomic_sparsity%& %s% |
---|
961 | %$$ |
---|
962 | The input values of its elements |
---|
963 | are not specified (must not matter). |
---|
964 | Upon return, $icode s$$ is a |
---|
965 | $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
966 | $latex S(x) \in B^{m \times q}$$. |
---|
967 | |
---|
968 | $head ok$$ |
---|
969 | The return value $icode ok$$ has prototype |
---|
970 | $codei% |
---|
971 | bool %ok% |
---|
972 | %$$ |
---|
973 | If it is $code true$$, the corresponding evaluation succeeded, |
---|
974 | otherwise it failed. |
---|
975 | |
---|
976 | $end |
---|
977 | ----------------------------------------------------------------------------- |
---|
978 | */ |
---|
979 | /*! |
---|
980 | Link from forward Jacobian sparsity sweep to atomic_base |
---|
981 | |
---|
982 | \param q |
---|
983 | is the column dimension for the Jacobian sparsity partterns. |
---|
984 | |
---|
985 | \param r |
---|
986 | is the Jacobian sparsity pattern for the argument vector x |
---|
987 | |
---|
988 | \param s |
---|
989 | is the Jacobian sparsity pattern for the result vector y |
---|
990 | */ |
---|
991 | virtual bool for_sparse_jac( |
---|
992 | size_t q , |
---|
993 | const vector< std::set<size_t> >& r , |
---|
994 | vector< std::set<size_t> >& s ) |
---|
995 | { return false; } |
---|
996 | virtual bool for_sparse_jac( |
---|
997 | size_t q , |
---|
998 | const vector<bool>& r , |
---|
999 | vector<bool>& s ) |
---|
1000 | { return false; } |
---|
1001 | /* |
---|
1002 | -------------------------------------- --------------------------------------- |
---|
1003 | $begin atomic_rev_sparse_jac$$ |
---|
1004 | $spell |
---|
1005 | rt |
---|
1006 | afun |
---|
1007 | Jacobian |
---|
1008 | jac |
---|
1009 | CppAD |
---|
1010 | std |
---|
1011 | bool |
---|
1012 | const |
---|
1013 | hes |
---|
1014 | $$ |
---|
1015 | |
---|
1016 | $section Atomic Reverse Jacobian Sparsity Patterns$$ |
---|
1017 | $index atomic, rev_sparse_jac callback$$ |
---|
1018 | $index rev_sparse_jac, atomic callback$$ |
---|
1019 | $index rev_sparse_jac, atomic virtual$$ |
---|
1020 | |
---|
1021 | $head Syntax$$ |
---|
1022 | $icode%ok% = %afun%.rev_sparse_jac(%q%, %rt%, %st%)%$$ |
---|
1023 | |
---|
1024 | $head Purpose$$ |
---|
1025 | This function is used by $cref RevSparseJac$$ to compute |
---|
1026 | Jacobian sparsity patterns. |
---|
1027 | For a fixed matrix $latex R \in B^{q \times m}$$, |
---|
1028 | the Jacobian of $latex R * f( x )$$ with respect to $latex x \in B^q$$ is |
---|
1029 | $latex \[ |
---|
1030 | S(x) = R * f^{(1)} (x) |
---|
1031 | \] $$ |
---|
1032 | Given a $cref/sparsity pattern/glossary/Sparsity Pattern/$$ for $latex R$$, |
---|
1033 | $code rev_sparse_jac$$ computes a sparsity pattern for $latex S(x)$$. |
---|
1034 | |
---|
1035 | $head Implementation$$ |
---|
1036 | If you are using $cref RevSparseJac$$, |
---|
1037 | this virtual function must be defined by the |
---|
1038 | $cref/atomic_user/atomic_ctor/atomic_user/$$ class. |
---|
1039 | |
---|
1040 | $subhead q$$ |
---|
1041 | The argument $icode q$$ has prototype |
---|
1042 | $codei% |
---|
1043 | size_t %q% |
---|
1044 | %$$ |
---|
1045 | It specifies the number of rows in |
---|
1046 | $latex R \in B^{q \times m}$$ and the Jacobian |
---|
1047 | $latex S(x) \in B^{q \times n}$$. |
---|
1048 | |
---|
1049 | $subhead rt$$ |
---|
1050 | This argument has prototype |
---|
1051 | $codei% |
---|
1052 | const %atomic_sparsity%& %rt% |
---|
1053 | %$$ |
---|
1054 | and is a |
---|
1055 | $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
1056 | $latex R^\R{T} \in B^{m \times q}$$. |
---|
1057 | |
---|
1058 | $subhead st$$ |
---|
1059 | This argument has prototype |
---|
1060 | $codei% |
---|
1061 | %atomic_sparsity%& %st% |
---|
1062 | %$$ |
---|
1063 | The input value of its elements |
---|
1064 | are not specified (must not matter). |
---|
1065 | Upon return, $icode s$$ is a |
---|
1066 | $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
1067 | $latex S(x)^\R{T} \in B^{n \times q}$$. |
---|
1068 | |
---|
1069 | $head ok$$ |
---|
1070 | The return value $icode ok$$ has prototype |
---|
1071 | $codei% |
---|
1072 | bool %ok% |
---|
1073 | %$$ |
---|
1074 | If it is $code true$$, the corresponding evaluation succeeded, |
---|
1075 | otherwise it failed. |
---|
1076 | |
---|
1077 | $end |
---|
1078 | ----------------------------------------------------------------------------- |
---|
1079 | */ |
---|
1080 | /*! |
---|
1081 | Link from reverse Jacobian sparsity sweep to atomic_base |
---|
1082 | |
---|
1083 | \param q [in] |
---|
1084 | is the row dimension for the Jacobian sparsity partterns |
---|
1085 | |
---|
1086 | \param rt [out] |
---|
1087 | is the tansposed Jacobian sparsity pattern w.r.t to range variables y |
---|
1088 | |
---|
1089 | \param st [in] |
---|
1090 | is the tansposed Jacobian sparsity pattern for the argument variables x |
---|
1091 | */ |
---|
1092 | virtual bool rev_sparse_jac( |
---|
1093 | size_t q , |
---|
1094 | const vector< std::set<size_t> >& rt , |
---|
1095 | vector< std::set<size_t> >& st ) |
---|
1096 | { return false; } |
---|
1097 | virtual bool rev_sparse_jac( |
---|
1098 | size_t q , |
---|
1099 | const vector<bool>& rt , |
---|
1100 | vector<bool>& st ) |
---|
1101 | { return false; } |
---|
1102 | /* |
---|
1103 | -------------------------------------- --------------------------------------- |
---|
1104 | $begin atomic_rev_sparse_hes$$ |
---|
1105 | $spell |
---|
1106 | vx |
---|
1107 | afun |
---|
1108 | Jacobian |
---|
1109 | jac |
---|
1110 | CppAD |
---|
1111 | std |
---|
1112 | bool |
---|
1113 | hes |
---|
1114 | const |
---|
1115 | $$ |
---|
1116 | |
---|
1117 | $section Atomic Reverse Hessian Sparsity Patterns$$ |
---|
1118 | $index atomic, rev_sparse_hes callback$$ |
---|
1119 | $index rev_sparse_hes, atomic callback$$ |
---|
1120 | $index rev_sparse_hes, atomic virtual$$ |
---|
1121 | |
---|
1122 | $head Syntax$$ |
---|
1123 | $icode%ok% = %afun%.rev_sparse_hes(%vx%, %s%, %t%, %q%, %r%, %u%, %v%)%$$ |
---|
1124 | |
---|
1125 | $head Purpose$$ |
---|
1126 | This function is used by $cref RevSparseHes$$ to compute |
---|
1127 | Hessian sparsity patterns. |
---|
1128 | There is an unspecified scalar valued function |
---|
1129 | $latex g : B^m \rightarrow B$$. |
---|
1130 | Given a $cref/sparsity pattern/glossary/Sparsity Pattern/$$ for |
---|
1131 | $latex R \in B^{n \times q}$$, |
---|
1132 | and information about the function $latex z = g(y)$$, |
---|
1133 | this routine computes the sparsity pattern for |
---|
1134 | $latex \[ |
---|
1135 | V(x) = (g \circ f)^{(2)}( x ) R |
---|
1136 | \] $$ |
---|
1137 | |
---|
1138 | $head Implementation$$ |
---|
1139 | If you are using and $cref RevSparseHes$$, |
---|
1140 | this virtual function must be defined by the |
---|
1141 | $cref/atomic_user/atomic_ctor/atomic_user/$$ class. |
---|
1142 | |
---|
1143 | $subhead vx$$ |
---|
1144 | The argument $icode vx$$ has prototype |
---|
1145 | $codei% |
---|
1146 | const CppAD:vector<bool>& %vx% |
---|
1147 | %$$ |
---|
1148 | $icode%vx%.size() == %n%$$, and |
---|
1149 | for $latex j = 0 , \ldots , n-1$$, |
---|
1150 | $icode%vx%[%j%]%$$ is true if and only if |
---|
1151 | $icode%ax%[%j%]%$$ is a $cref/variable/glossary/Variable/$$ |
---|
1152 | in the corresponding call to |
---|
1153 | $codei% |
---|
1154 | %afun%(%ax%, %ay%, %id%) |
---|
1155 | %$$ |
---|
1156 | |
---|
1157 | $subhead s$$ |
---|
1158 | The argument $icode s$$ has prototype |
---|
1159 | $codei% |
---|
1160 | const CppAD:vector<bool>& %s% |
---|
1161 | %$$ |
---|
1162 | and its size is $icode m$$. |
---|
1163 | It is a sparsity pattern for |
---|
1164 | $latex S(x) = g^{(1)} (y) \in B^{1 \times m}$$. |
---|
1165 | |
---|
1166 | $subhead t$$ |
---|
1167 | This argument has prototype |
---|
1168 | $codei% |
---|
1169 | CppAD:vector<bool>& %t% |
---|
1170 | %$$ |
---|
1171 | and its size is $icode m$$. |
---|
1172 | The input values of its elements |
---|
1173 | are not specified (must not matter). |
---|
1174 | Upon return, $icode t$$ is a |
---|
1175 | sparsity pattern for |
---|
1176 | $latex T(x) \in B^{1 \times n}$$ where |
---|
1177 | $latex \[ |
---|
1178 | T(x) = (g \circ f)^{(1)} (x) = S(x) * f^{(1)} (x) |
---|
1179 | \]$$ |
---|
1180 | |
---|
1181 | $subhead q$$ |
---|
1182 | The argument $icode q$$ has prototype |
---|
1183 | $codei% |
---|
1184 | size_t %q% |
---|
1185 | %$$ |
---|
1186 | It specifies the number of columns in |
---|
1187 | $latex R \in B^{n \times q}$$, |
---|
1188 | $latex U(x) \in B^{m \times q}$$, and |
---|
1189 | $latex V(x) \in B^{n \times q}$$. |
---|
1190 | |
---|
1191 | $subhead r$$ |
---|
1192 | This argument has prototype |
---|
1193 | $codei% |
---|
1194 | const %atomic_sparsity%& %r% |
---|
1195 | %$$ |
---|
1196 | and is a $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
1197 | $latex R \in B^{n \times q}$$. |
---|
1198 | |
---|
1199 | $head u$$ |
---|
1200 | This argument has prototype |
---|
1201 | $codei% |
---|
1202 | const %atomic_sparsity%& %u% |
---|
1203 | %$$ |
---|
1204 | and is a $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
1205 | $latex U(x) \in B^{m \times q}$$ which is defined by |
---|
1206 | $latex \[ |
---|
1207 | \begin{array}{rcl} |
---|
1208 | U(x) |
---|
1209 | & = & |
---|
1210 | \partial_u \{ \partial_y g[ y + f^{(1)} (x) R u ] \}_{u=0} |
---|
1211 | \\ |
---|
1212 | & = & |
---|
1213 | \partial_u \{ g^{(1)} [ y + f^{(1)} (x) R u ] \}_{u=0} |
---|
1214 | \\ |
---|
1215 | & = & |
---|
1216 | g^{(2)} (y) f^{(1)} (x) R |
---|
1217 | \end{array} |
---|
1218 | \] $$ |
---|
1219 | |
---|
1220 | $subhead v$$ |
---|
1221 | This argument has prototype |
---|
1222 | $codei% |
---|
1223 | %atomic_sparsity%& %v% |
---|
1224 | %$$ |
---|
1225 | The input value of its elements |
---|
1226 | are not specified (must not matter). |
---|
1227 | Upon return, $icode v$$ is a |
---|
1228 | $cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for |
---|
1229 | $latex V(x) \in B^{n \times q}$$ which is defined by |
---|
1230 | $latex \[ |
---|
1231 | \begin{array}{rcl} |
---|
1232 | V(x) |
---|
1233 | & = & |
---|
1234 | \partial_u [ \partial_x (g \circ f) ( x + R u ) ]_{u=0} |
---|
1235 | \\ |
---|
1236 | & = & |
---|
1237 | \partial_u [ (g \circ f)^{(1)}( x + R u ) ]_{u=0} |
---|
1238 | \\ |
---|
1239 | & = & |
---|
1240 | (g \circ f)^{(2)}( x ) R |
---|
1241 | \\ |
---|
1242 | & = & |
---|
1243 | f^{(1)} (x)^\R{T} g^{(2)} ( y ) f^{(1)} (x) R |
---|
1244 | + |
---|
1245 | \sum_{i=1}^m g_i^{(1)} (y) \; f_i^{(2)} (x) R |
---|
1246 | \\ |
---|
1247 | & = & |
---|
1248 | f^{(1)} (x)^\R{T} U(x) |
---|
1249 | + |
---|
1250 | \sum_{i=1}^m S_i (x) \; f_i^{(2)} (x) R |
---|
1251 | \end{array} |
---|
1252 | \] $$ |
---|
1253 | |
---|
1254 | $end |
---|
1255 | ----------------------------------------------------------------------------- |
---|
1256 | */ |
---|
1257 | /*! |
---|
1258 | Link from reverse Hessian sparsity sweep to base_atomic |
---|
1259 | |
---|
1260 | \param vx [in] |
---|
1261 | which componens of x are variables. |
---|
1262 | |
---|
1263 | \param s [in] |
---|
1264 | is the reverse Jacobian sparsity pattern w.r.t the result vector y. |
---|
1265 | |
---|
1266 | \param t [out] |
---|
1267 | is the reverse Jacobian sparsity pattern w.r.t the argument vector x. |
---|
1268 | |
---|
1269 | \param q [in] |
---|
1270 | is the column dimension for the sparsity partterns. |
---|
1271 | |
---|
1272 | \param r [in] |
---|
1273 | is the forward Jacobian sparsity pattern w.r.t the argument vector x |
---|
1274 | |
---|
1275 | \param u [in] |
---|
1276 | is the Hessian sparsity pattern w.r.t the result vector y. |
---|
1277 | |
---|
1278 | \param v [out] |
---|
1279 | is the Hessian sparsity pattern w.r.t the argument vector x. |
---|
1280 | */ |
---|
1281 | virtual bool rev_sparse_hes( |
---|
1282 | const vector<bool>& vx , |
---|
1283 | const vector<bool>& s , |
---|
1284 | vector<bool>& t , |
---|
1285 | size_t q , |
---|
1286 | const vector< std::set<size_t> >& r , |
---|
1287 | const vector< std::set<size_t> >& u , |
---|
1288 | vector< std::set<size_t> >& v ) |
---|
1289 | { return false; } |
---|
1290 | virtual bool rev_sparse_hes( |
---|
1291 | const vector<bool>& vx , |
---|
1292 | const vector<bool>& s , |
---|
1293 | vector<bool>& t , |
---|
1294 | size_t q , |
---|
1295 | const vector<bool>& r , |
---|
1296 | const vector<bool>& u , |
---|
1297 | vector<bool>& v ) |
---|
1298 | { return false; } |
---|
1299 | /* |
---|
1300 | ------------------------------------------------------------------------------ |
---|
1301 | $begin atomic_base_clear$$ |
---|
1302 | |
---|
1303 | $section Free Static Variables$$ |
---|
1304 | $index free, atomic static$$ |
---|
1305 | $index atomic, free static$$ |
---|
1306 | $index static, free atomic$$ |
---|
1307 | $index clear, atomic static$$ |
---|
1308 | |
---|
1309 | $head Syntax$$ |
---|
1310 | $codei%atomic_base<%Base%>::clear()%$$ |
---|
1311 | |
---|
1312 | $head Purpose$$ |
---|
1313 | The $code atomic_base$$ class holds onto static work space in order to |
---|
1314 | increase speed by avoiding system memory allocation calls. |
---|
1315 | This call makes to work space $cref/available/ta_available/$$ to |
---|
1316 | for other uses by the same thread. |
---|
1317 | This should be called when you are done using the |
---|
1318 | user atomic functions for a specific value of $icode Base$$. |
---|
1319 | |
---|
1320 | $head Restriction$$ |
---|
1321 | This routine cannot be called |
---|
1322 | while in $cref/parallel/ta_in_parallel/$$ execution mode. |
---|
1323 | |
---|
1324 | $end |
---|
1325 | ------------------------------------------------------------------------------ |
---|
1326 | */ |
---|
1327 | /*! |
---|
1328 | Free all thread_alloc static memory held by atomic_base (avoids reallocations). |
---|
1329 | (This does not include class_object() which is an std::vector.) |
---|
1330 | */ |
---|
1331 | /// Free vector memory used by this class (work space) |
---|
1332 | static void clear(void) |
---|
1333 | { CPPAD_ASSERT_KNOWN( |
---|
1334 | ! thread_alloc::in_parallel() , |
---|
1335 | "cannot use atomic_base clear during parallel execution" |
---|
1336 | ); |
---|
1337 | size_t i = class_object().size(); |
---|
1338 | while(i--) |
---|
1339 | { size_t thread = CPPAD_MAX_NUM_THREADS; |
---|
1340 | while(thread--) |
---|
1341 | { |
---|
1342 | atomic_base* op = class_object()[i]; |
---|
1343 | if( op != CPPAD_NULL ) |
---|
1344 | { op->afun_vx_[thread].clear(); |
---|
1345 | op->afun_vy_[thread].clear(); |
---|
1346 | op->afun_tx_[thread].clear(); |
---|
1347 | op->afun_ty_[thread].clear(); |
---|
1348 | } |
---|
1349 | } |
---|
1350 | } |
---|
1351 | return; |
---|
1352 | } |
---|
1353 | // ------------------------------------------------------------------------- |
---|
1354 | /*! |
---|
1355 | Set value of id (used by deprecated old_atomic class) |
---|
1356 | |
---|
1357 | This function is called just before calling any of the virtual funcitons |
---|
1358 | and has the corresponding id of the corresponding virtual call. |
---|
1359 | */ |
---|
1360 | virtual void set_id(size_t id) |
---|
1361 | { } |
---|
1362 | // --------------------------------------------------------------------------- |
---|
1363 | }; |
---|
1364 | /*! \} */ |
---|
1365 | CPPAD_END_NAMESPACE |
---|
1366 | # endif |
---|