1 | /* $Id: fun_construct.hpp 2970 2013-10-19 02:16:02Z bradbell $ */ |
---|
2 | # ifndef CPPAD_FUN_CONSTRUCT_INCLUDED |
---|
3 | # define CPPAD_FUN_CONSTRUCT_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 | $begin FunConstruct$$ |
---|
17 | $spell |
---|
18 | alloc |
---|
19 | num |
---|
20 | Jac |
---|
21 | bool |
---|
22 | taylor_ |
---|
23 | var |
---|
24 | ADvector |
---|
25 | const |
---|
26 | Jacobian |
---|
27 | $$ |
---|
28 | |
---|
29 | $spell |
---|
30 | $$ |
---|
31 | |
---|
32 | $section Construct an ADFun Object and Stop Recording$$ |
---|
33 | |
---|
34 | $index ADFun, construct$$ |
---|
35 | $index construct, ADFun$$ |
---|
36 | $index tape, stop recording$$ |
---|
37 | $index recording, stop tape$$ |
---|
38 | |
---|
39 | $head Syntax$$ |
---|
40 | $codei%ADFun<%Base%> %f%, %g% |
---|
41 | %$$ |
---|
42 | $codei%ADFun<%Base%> %f%(%x%, %y%) |
---|
43 | %$$ |
---|
44 | $icode%g% = %f% |
---|
45 | %$$ |
---|
46 | |
---|
47 | |
---|
48 | $head Purpose$$ |
---|
49 | The $codei%AD<%Base%>%$$ object $icode f$$ can |
---|
50 | store an AD of $icode Base$$ |
---|
51 | $cref/operation sequence/glossary/Operation/Sequence/$$. |
---|
52 | It can then be used to calculate derivatives of the corresponding |
---|
53 | $cref/AD function/glossary/AD Function/$$ |
---|
54 | $latex \[ |
---|
55 | F : B^n \rightarrow B^m |
---|
56 | \] $$ |
---|
57 | where $latex B$$ is the space corresponding to objects of type $icode Base$$. |
---|
58 | |
---|
59 | $head x$$ |
---|
60 | If the argument $icode x$$ is present, it has prototype |
---|
61 | $codei% |
---|
62 | const %VectorAD% &%x% |
---|
63 | %$$ |
---|
64 | It must be the vector argument in the previous call to |
---|
65 | $cref Independent$$. |
---|
66 | Neither its size, or any of its values, are allowed to change |
---|
67 | between calling |
---|
68 | $codei% |
---|
69 | Independent(%x%) |
---|
70 | %$$ |
---|
71 | and |
---|
72 | $codei% |
---|
73 | ADFun<%Base%> %f%(%x%, %y%) |
---|
74 | %$$. |
---|
75 | |
---|
76 | $head y$$ |
---|
77 | If the argument $icode y$$ is present, it has prototype |
---|
78 | $codei% |
---|
79 | const %VectorAD% &%y% |
---|
80 | %$$ |
---|
81 | The sequence of operations that map $icode x$$ |
---|
82 | to $icode y$$ are stored in the ADFun object $icode f$$. |
---|
83 | |
---|
84 | $head VectorAD$$ |
---|
85 | The type $icode VectorAD$$ must be a $cref SimpleVector$$ class with |
---|
86 | $cref/elements of type/SimpleVector/Elements of Specified Type/$$ |
---|
87 | $codei%AD<%Base%>%$$. |
---|
88 | The routine $cref CheckSimpleVector$$ will generate an error message |
---|
89 | if this is not the case. |
---|
90 | |
---|
91 | $head Default Constructor$$ |
---|
92 | $index default, ADFun constructor$$ |
---|
93 | $index ADFun, default constructor$$ |
---|
94 | $index constructor, ADFun constructor$$ |
---|
95 | The default constructor |
---|
96 | $codei% |
---|
97 | ADFun<%Base%> %f% |
---|
98 | %$$ |
---|
99 | creates an |
---|
100 | $codei%AD<%Base%>%$$ object with no corresponding operation sequence; i.e., |
---|
101 | $codei% |
---|
102 | %f%.size_var() |
---|
103 | %$$ |
---|
104 | returns the value zero (see $cref/size_var/seq_property/size_var/$$). |
---|
105 | |
---|
106 | $head Sequence Constructor$$ |
---|
107 | $index sequence, ADFun constructor$$ |
---|
108 | $index ADFun, sequence constructor$$ |
---|
109 | $index constructor, ADFun sequence$$ |
---|
110 | The default constructor |
---|
111 | The sequence constructor |
---|
112 | $codei% |
---|
113 | ADFun<%Base%> %f%(%x%, %y%) |
---|
114 | %$$ |
---|
115 | creates the $codei%AD<%Base%>%$$ object $icode f$$, |
---|
116 | stops the recording of AD of $icode Base$$ operations |
---|
117 | corresponding to the call |
---|
118 | $codei% |
---|
119 | Independent(%x%) |
---|
120 | %$$ |
---|
121 | and stores the corresponding operation sequence in the object $icode f$$. |
---|
122 | It then stores the first order taylor_ coefficients |
---|
123 | (corresponding to the value of $icode x$$) in $icode f$$. |
---|
124 | This is equivalent to the following steps using the default constructor: |
---|
125 | |
---|
126 | $list number$$ |
---|
127 | Create $icode f$$ with the default constructor |
---|
128 | $codei% |
---|
129 | ADFun<%Base%> %f%; |
---|
130 | %$$ |
---|
131 | $lnext |
---|
132 | Stop the tape and storing the operation sequence using |
---|
133 | $codei% |
---|
134 | %f%.Dependent(%x%, %y%); |
---|
135 | %$$ |
---|
136 | (see $cref Dependent$$). |
---|
137 | $lnext |
---|
138 | Calculating the first order taylor_ coefficients for all |
---|
139 | the variables in the operation sequence using |
---|
140 | $codei% |
---|
141 | %f%.Forward(%p%, %x_p%) |
---|
142 | %$$ |
---|
143 | with $icode p$$ equal to zero and the elements of $icode x_p$$ |
---|
144 | equal to the corresponding elements of $icode x$$ |
---|
145 | (see $cref Forward$$). |
---|
146 | $lend |
---|
147 | |
---|
148 | $head Copy Constructor$$ |
---|
149 | $index copy, ADFun constructor$$ |
---|
150 | $index ADFun, copy constructor$$ |
---|
151 | $index constructor, ADFun copy$$ |
---|
152 | It is an error to attempt to use the $codei%ADFun<%Base%>%$$ copy constructor; |
---|
153 | i.e., the following syntax is not allowed: |
---|
154 | $codei% |
---|
155 | ADFun<%Base%> %g%(%f%) |
---|
156 | %$$ |
---|
157 | where $icode f$$ is an $codei%ADFun<%Base%>%$$ object. |
---|
158 | Use its $cref/default constructor/FunConstruct/Default Constructor/$$ instead |
---|
159 | and its assignment operator. |
---|
160 | |
---|
161 | $head Assignment Operator$$ |
---|
162 | $index ADFun, assignment operator$$ |
---|
163 | $index assignment, ADFun operator$$ |
---|
164 | $index operator, ADFun assignment$$ |
---|
165 | The $codei%ADFun<%Base%>%$$ assignment operation |
---|
166 | $codei% |
---|
167 | %g% = %f% |
---|
168 | %$$. |
---|
169 | makes a copy of the operation sequence currently stored in $icode f$$ |
---|
170 | in the object $icode g$$. |
---|
171 | The object $icode f$$ is not affected by this operation and |
---|
172 | can be $code const$$. |
---|
173 | Any operation sequence or other information in $icode g$$ is lost. |
---|
174 | |
---|
175 | $subhead Taylor Coefficients$$ |
---|
176 | The Taylor coefficient information currently stored in $icode f$$ |
---|
177 | (computed by $cref/f.Forward/Forward/$$) is |
---|
178 | copied to $icode g$$. |
---|
179 | Hence, directly after this operation |
---|
180 | $codei% |
---|
181 | %g%.size_taylor() == %f%.size_taylor() |
---|
182 | %$$ |
---|
183 | |
---|
184 | $subhead Sparsity Patterns$$ |
---|
185 | The forward Jacobian sparsity pattern currently stored in $icode f$$ |
---|
186 | (computed by $cref/f.ForSparseJac/ForSparseJac/$$) is |
---|
187 | copied to $icode g$$. |
---|
188 | Hence, directly after this operation |
---|
189 | $codei% |
---|
190 | %g%.size_forward_bool() == %f%.size_forward_bool() |
---|
191 | %g%.size_forward_set() == %f%.size_forward_set() |
---|
192 | %$$ |
---|
193 | |
---|
194 | $head Parallel Mode$$ |
---|
195 | $index parallel, ADFun$$ |
---|
196 | $index ADFun, parallel$$ |
---|
197 | The call to $code Independent$$, |
---|
198 | and the corresponding call to |
---|
199 | $codei% |
---|
200 | ADFun<%Base%> %f%( %x%, %y%) |
---|
201 | %$$ |
---|
202 | or |
---|
203 | $codei% |
---|
204 | %f%.Dependent( %x%, %y%) |
---|
205 | %$$ |
---|
206 | or $cref abort_recording$$, |
---|
207 | must be preformed by the same thread; i.e., |
---|
208 | $cref/thread_alloc::thread_num/ta_thread_num/$$ must be the same. |
---|
209 | |
---|
210 | $head Example$$ |
---|
211 | |
---|
212 | $subhead Sequence Constructor$$ |
---|
213 | The file |
---|
214 | $cref independent.cpp$$ |
---|
215 | contains an example and test of the sequence constructor. |
---|
216 | It returns true if it succeeds and false otherwise. |
---|
217 | |
---|
218 | $subhead Default Constructor$$ |
---|
219 | The files |
---|
220 | $cref fun_check.cpp$$ |
---|
221 | and |
---|
222 | $cref hes_lagrangian.cpp$$ |
---|
223 | contain an examples and tests using the default constructor. |
---|
224 | They return true if they succeed and false otherwise. |
---|
225 | |
---|
226 | $children% |
---|
227 | example/fun_assign.cpp |
---|
228 | %$$ |
---|
229 | $subhead Assignment Operator$$ |
---|
230 | The file |
---|
231 | $cref fun_assign.cpp$$ |
---|
232 | contains an example and test of the $codei%ADFun<%Base%>%$$ |
---|
233 | assignment operator. |
---|
234 | It returns true if it succeeds and false otherwise. |
---|
235 | |
---|
236 | $end |
---|
237 | ---------------------------------------------------------------------------- |
---|
238 | */ |
---|
239 | |
---|
240 | namespace CppAD { // BEGIN_CPPAD_NAMESPACE |
---|
241 | /*! |
---|
242 | \defgroup fun_construct_hpp fun_construct.hpp |
---|
243 | \{ |
---|
244 | \file fun_construct.hpp |
---|
245 | ADFun function constructors and assignment operator. |
---|
246 | */ |
---|
247 | |
---|
248 | /*! |
---|
249 | ADFun default constructor |
---|
250 | |
---|
251 | The C++ syntax for this operation is |
---|
252 | \verbatim |
---|
253 | ADFun<Base> f |
---|
254 | \endverbatim |
---|
255 | An empty ADFun object is created. |
---|
256 | The Dependent member function, |
---|
257 | or the ADFun<Base> assingment operator, |
---|
258 | can then be used to put an operation sequence in this ADFun object. |
---|
259 | |
---|
260 | \tparam Base |
---|
261 | is the base for the recording that can be stored in this ADFun object; |
---|
262 | i.e., operation sequences that were recorded using the type \c AD<Base>. |
---|
263 | */ |
---|
264 | template <typename Base> |
---|
265 | ADFun<Base>::ADFun(void) : |
---|
266 | check_for_nan_(true) , |
---|
267 | total_num_var_(0) |
---|
268 | { } |
---|
269 | |
---|
270 | /*! |
---|
271 | ADFun assignment operator |
---|
272 | |
---|
273 | The C++ syntax for this operation is |
---|
274 | \verbatim |
---|
275 | g = f |
---|
276 | \endverbatim |
---|
277 | where \c g and \c f are ADFun<Base> ADFun objects. |
---|
278 | A copy of the the operation sequence currently stored in \c f |
---|
279 | is placed in this ADFun object (called \c g above). |
---|
280 | Any information currently stored in this ADFun object is lost. |
---|
281 | |
---|
282 | \tparam Base |
---|
283 | is the base for the recording that can be stored in this ADFun object; |
---|
284 | i.e., operation sequences that were recorded using the type \c AD<Base>. |
---|
285 | |
---|
286 | \param f |
---|
287 | ADFun object containing the operation sequence to be copied. |
---|
288 | */ |
---|
289 | template <typename Base> |
---|
290 | void ADFun<Base>::operator=(const ADFun<Base>& f) |
---|
291 | { size_t m = f.Range(); |
---|
292 | size_t n = f.Domain(); |
---|
293 | |
---|
294 | // go through member variables in order |
---|
295 | // (see ad_fun.hpp for meaning of each variable) |
---|
296 | check_for_nan_ = true; |
---|
297 | compare_change_ = 0; |
---|
298 | |
---|
299 | taylor_.erase(); |
---|
300 | taylor_per_var_ = 0; |
---|
301 | taylor_col_dim_ = 0; |
---|
302 | |
---|
303 | total_num_var_ = f.total_num_var_; |
---|
304 | ind_taddr_.resize(n); |
---|
305 | ind_taddr_ = f.ind_taddr_; |
---|
306 | dep_taddr_.resize(m); |
---|
307 | dep_taddr_ = f.dep_taddr_; |
---|
308 | dep_parameter_.resize(m); |
---|
309 | dep_parameter_ = f.dep_parameter_; |
---|
310 | play_ = f.play_; |
---|
311 | for_jac_sparse_pack_.resize(0, 0); |
---|
312 | for_jac_sparse_set_.resize(0, 0); |
---|
313 | |
---|
314 | // allocate and copy the Taylor coefficients |
---|
315 | taylor_per_var_ = f.taylor_per_var_; |
---|
316 | taylor_col_dim_ = f.taylor_col_dim_; |
---|
317 | size_t length = total_num_var_ * taylor_col_dim_; |
---|
318 | if( length > 0 ) |
---|
319 | taylor_.extend(length); |
---|
320 | size_t i, j; |
---|
321 | for(i = 0; i < total_num_var_; i++) |
---|
322 | { for(j = 0; j < taylor_per_var_; j++) |
---|
323 | { taylor_[ i * taylor_col_dim_ + j ] = |
---|
324 | f.taylor_[ i * taylor_col_dim_ + j ]; |
---|
325 | } |
---|
326 | } |
---|
327 | |
---|
328 | // allocate and copy the conditional skip information |
---|
329 | cskip_var_.clear(); |
---|
330 | cskip_var_ = f.cskip_var_; |
---|
331 | |
---|
332 | // allocate and copy the forward sparsity information |
---|
333 | size_t n_set = f.for_jac_sparse_pack_.n_set(); |
---|
334 | size_t end = f.for_jac_sparse_pack_.end(); |
---|
335 | if( n_set > 0 ) |
---|
336 | { CPPAD_ASSERT_UNKNOWN( n_set == total_num_var_ ); |
---|
337 | CPPAD_ASSERT_UNKNOWN( f.for_jac_sparse_set_.n_set() == 0 ); |
---|
338 | for_jac_sparse_pack_.resize(n_set, end); |
---|
339 | for(i = 0; i < total_num_var_ ; i++) |
---|
340 | { for_jac_sparse_pack_.assignment( |
---|
341 | i , |
---|
342 | i , |
---|
343 | f.for_jac_sparse_pack_ |
---|
344 | ); |
---|
345 | } |
---|
346 | } |
---|
347 | n_set = f.for_jac_sparse_set_.n_set(); |
---|
348 | end = f.for_jac_sparse_set_.end(); |
---|
349 | if( n_set > 0 ) |
---|
350 | { CPPAD_ASSERT_UNKNOWN( n_set == total_num_var_ ); |
---|
351 | CPPAD_ASSERT_UNKNOWN( f.for_jac_sparse_pack_.n_set() == 0 ); |
---|
352 | for_jac_sparse_set_.resize(n_set, end); |
---|
353 | for(i = 0; i < total_num_var_; i++) |
---|
354 | { for_jac_sparse_set_.assignment( |
---|
355 | i , |
---|
356 | i , |
---|
357 | f.for_jac_sparse_set_ |
---|
358 | ); |
---|
359 | } |
---|
360 | } |
---|
361 | |
---|
362 | } |
---|
363 | |
---|
364 | /*! |
---|
365 | ADFun constructor from an operation sequence. |
---|
366 | |
---|
367 | The C++ syntax for this operation is |
---|
368 | \verbatim |
---|
369 | ADFun<Base> f(x, y) |
---|
370 | \endverbatim |
---|
371 | The operation sequence that started with the previous call |
---|
372 | \c Independent(x), and that ends with this operation, is stored |
---|
373 | in this \c ADFun<Base> object \c f. |
---|
374 | |
---|
375 | \tparam Base |
---|
376 | is the base for the recording that will be stored in the object \c f; |
---|
377 | i.e., the operations were recorded using the type \c AD<Base>. |
---|
378 | |
---|
379 | \tparam VectorAD |
---|
380 | is a simple vector class with elements of typea \c AD<Base>. |
---|
381 | |
---|
382 | \param x |
---|
383 | is the independent variable vector for this ADFun object. |
---|
384 | The domain dimension of this object will be the size of \a x. |
---|
385 | |
---|
386 | \param y |
---|
387 | is the dependent variable vector for this ADFun object. |
---|
388 | The range dimension of this object will be the size of \a y. |
---|
389 | |
---|
390 | \par Taylor Coefficients |
---|
391 | A zero order forward mode sweep is done, |
---|
392 | and if NDEBUG is not defined the resulting values for the |
---|
393 | depenedent variables are checked against the values in \a y. |
---|
394 | Thus, the zero order Taylor coefficients |
---|
395 | corresponding to the value of the \a x vector |
---|
396 | are stored in this ADFun object. |
---|
397 | */ |
---|
398 | template <typename Base> |
---|
399 | template <typename VectorAD> |
---|
400 | ADFun<Base>::ADFun(const VectorAD &x, const VectorAD &y) : |
---|
401 | check_for_nan_(true) , |
---|
402 | total_num_var_(0) |
---|
403 | { |
---|
404 | CPPAD_ASSERT_KNOWN( |
---|
405 | x.size() > 0, |
---|
406 | "ADFun<Base>: independent variable vector has size zero." |
---|
407 | ); |
---|
408 | CPPAD_ASSERT_KNOWN( |
---|
409 | Variable(x[0]), |
---|
410 | "ADFun<Base>: independent variable vector has been changed." |
---|
411 | ); |
---|
412 | ADTape<Base>* tape = AD<Base>::tape_ptr(x[0].tape_id_); |
---|
413 | CPPAD_ASSERT_KNOWN( |
---|
414 | tape->size_independent_ == size_t ( x.size() ), |
---|
415 | "ADFun<Base>: independent variable vector has been changed." |
---|
416 | ); |
---|
417 | size_t j, n = x.size(); |
---|
418 | # ifndef NDEBUG |
---|
419 | size_t i, m = y.size(); |
---|
420 | for(j = 0; j < n; j++) |
---|
421 | { CPPAD_ASSERT_KNOWN( |
---|
422 | size_t(x[j].taddr_) == (j+1), |
---|
423 | "ADFun<Base>: independent variable vector has been changed." |
---|
424 | ); |
---|
425 | CPPAD_ASSERT_KNOWN( |
---|
426 | x[j].tape_id_ == x[0].tape_id_, |
---|
427 | "ADFun<Base>: independent variable vector has been changed." |
---|
428 | ); |
---|
429 | } |
---|
430 | for(i = 0; i < m; i++) |
---|
431 | { CPPAD_ASSERT_KNOWN( |
---|
432 | CppAD::Parameter( y[i] ) | (y[i].tape_id_ == x[0].tape_id_) , |
---|
433 | "ADFun<Base>: dependent vector contains variables for" |
---|
434 | "\na different tape than the independent variables." |
---|
435 | ); |
---|
436 | } |
---|
437 | # endif |
---|
438 | |
---|
439 | // stop the tape and store the operation sequence |
---|
440 | Dependent(tape, y); |
---|
441 | |
---|
442 | // allocate memory for one zero order taylor_ coefficient |
---|
443 | taylor_per_var_ = 1; |
---|
444 | taylor_col_dim_ = 1; |
---|
445 | taylor_.extend(total_num_var_); |
---|
446 | |
---|
447 | // set zero order coefficients corresponding to indpendent variables |
---|
448 | CPPAD_ASSERT_UNKNOWN( n == ind_taddr_.size() ); |
---|
449 | for(j = 0; j < n; j++) |
---|
450 | { CPPAD_ASSERT_UNKNOWN( ind_taddr_[j] == (j+1) ); |
---|
451 | CPPAD_ASSERT_UNKNOWN( size_t(x[j].taddr_) == (j+1) ); |
---|
452 | taylor_[ ind_taddr_[j] ] = x[j].value_; |
---|
453 | } |
---|
454 | |
---|
455 | // use independent variable values to fill in values for others |
---|
456 | # if CPPAD_USE_FORWARD0SWEEP |
---|
457 | compare_change_ = forward0sweep(std::cout, false, |
---|
458 | n, total_num_var_, &play_, taylor_col_dim_, taylor_.data(), |
---|
459 | cskip_var_ |
---|
460 | ); |
---|
461 | # else |
---|
462 | size_t p = 0; |
---|
463 | compare_change_ = forward_sweep(std::cout, false, |
---|
464 | p, p, n, total_num_var_, &play_, taylor_col_dim_, taylor_.data(), |
---|
465 | cskip_var_ |
---|
466 | ); |
---|
467 | # endif |
---|
468 | CPPAD_ASSERT_UNKNOWN( compare_change_ == 0 ); |
---|
469 | |
---|
470 | # ifndef NDEBUG |
---|
471 | for(i = 0; i < m; i++) |
---|
472 | if( taylor_[dep_taddr_[i]] != y[i].value_ || isnan( y[i].value_ ) ) |
---|
473 | { using std::endl; |
---|
474 | std::ostringstream buf; |
---|
475 | buf << "A dependent variable value is not equal to " |
---|
476 | << "its tape evaluation value," << endl |
---|
477 | << "perhaps it is nan." << endl |
---|
478 | << "Dependent variable value = " |
---|
479 | << y[i].value_ << endl |
---|
480 | << "Tape evaluation value = " |
---|
481 | << taylor_[dep_taddr_[i]] << endl |
---|
482 | << "Difference = " |
---|
483 | << y[i].value_ - taylor_[dep_taddr_[i]] << endl |
---|
484 | ; |
---|
485 | CPPAD_ASSERT_KNOWN( |
---|
486 | 0, |
---|
487 | buf.str().c_str() |
---|
488 | ); |
---|
489 | } |
---|
490 | # endif |
---|
491 | } |
---|
492 | |
---|
493 | /*! \} */ |
---|
494 | } // END_CPPAD_NAMESPACE |
---|
495 | # endif |
---|