source: trunk/cppad/local/optimize.hpp @ 3607

Last change on this file since 3607 was 3607, checked in by bradbell, 5 years ago

merge to branch: trunk
from repository: https://github.com/coin-or/CppAD
start hash code: f4baab0ea7c2679ba9351b0439b6efebfa77ba04
end hash code: 894d554a00ceec8a3545f05eca62708a7b5cb43d

commit 894d554a00ceec8a3545f05eca62708a7b5cb43d
Author: Brad Bell <bradbell@…>
Date: Tue Jan 20 09:06:10 2015 -0700

Fix copyright end date.


whats_new_15.omh: Add comment about date of deprecation.

commit 611e982000168db91aba22b763c14bb78d57da47
Author: Brad Bell <bradbell@…>
Date: Tue Jan 20 08:53:00 2015 -0700

Squashed commit from old/compare_op to master:


In addition, fix copyright end date for some files, and add note about
deprecated date in whats_new_15.omh.


commit 6e46df5c850ecd58d7a886db4043bc3f2d4579d1
Author: Brad Bell <bradbell@…>
Date: Tue Jan 20 08:16:57 2015 -0700


Always return f.compare_change_op_index() == 0 after f.optimize().


checkpoint.hpp: ignore comparison operators.
fun_construct.hpp: remove double initilaization of values.
compare_change.cpp: demonstrate more features of new interface.
whats_new_15.omh: prepare for merging this branch into master.
wish_list.omh: update wish list item to new compare_change interface.


commit 45315907c70e5b383d984fb9498b54a474001af0
Author: Brad Bell <bradbell@…>
Date: Tue Jan 20 05:04:37 2015 -0700


Use the new f.compare_change_count(0) option in speed tests.


commit bb6e72befd6d01f1fb62c43b9b19471ffaa7cc2c
Author: Brad Bell <bradbell@…>
Date: Tue Jan 20 04:51:16 2015 -0700


Move CompareChange? -> deprecated and plan -> compare_change.


forward0sweep.hpp: skip comparison operator when count is zero.
forward1sweep.hpp: skip comparison operator when count is zero.
compare_change.cpp: demonstrate count == 0 case.


commit 622a13c568c612d9dfe9ccd1a01f4ac5f74ba824
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 23:17:42 2015 -0700


Add example and test of new compare change user API.


ad_fun.hpp: fix f.compare_change_op_index.
compare_change.cpp: Change example to use new API.
compare_change.cpp: Move old example here (just test now).


commit ec4c1613eae8df56fbf31e7b8711ce69cc41df83
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 21:12:11 2015 -0700


Implement the compare change plan, still needs an example and test.
Also have the change from 'plan' to just plain documentation.


commit a81a46f27011bee08ba072551044dc9f4a99a906
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 17:49:05 2015 -0700


Change name of compare_change functions and partially implement this new plan.


commit 146faad48a700a56362e74f9c3a3c39144a79738
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 14:22:40 2015 -0700


Branch: compare_op
plan.omh: change name of count function.


commit 35d91d126765d1a0ab4bfe9e2b006bbf535cd648
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 13:19:07 2015 -0700


Add deprecation date for some more cases.


commit 5bb65a8c48fae4263b66fcd04520e10e66febc11
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 13:13:51 2015 -0700


Add date of deprecation for some more cases.


commit e95ee6b209601cd9a075d2e37c602e73c32fb6ab
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 12:58:44 2015 -0700


Add date of deprecation for some more cases.


commit 0ea84ccd87383edc62a6ae1711da104b12e8c444
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 12:47:01 2015 -0700


Add date of deprecation for some cases.


commit 17755e609ea8e03472b08dcc2fb5ad347eb723cb
Author: Brad Bell <bradbell@…>
Date: Mon Jan 19 08:20:45 2015 -0700


plan.omh: changs some names.


commit 29f369c06d4d0ee284c4c668d52d8461613066dc
Author: Brad Bell <bradbell@…>
Date: Fri Jan 16 06:39:45 2015 -0700


Document the plan for compare_change user API.


compare_change.omh: fix minor typo.
plan.txt: change to the omhelp file plan.omh.


commit a3a2f4dedd202a722812b6eb2714851b40726e6e
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 21:03:44 2015 -0700


new_branch.sh: remove unused variable.
push_git2svn.py: move repository from github/bradbell to coin-or/bradbell.


commit 3751a197ab2897e76616f9d9b0915148bd855356
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 20:56:17 2015 -0700


plan.txt: plan for this branches API will go here.


commit 76013ec2ad7baacdeab5e761812d542867910174
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 18:04:33 2015 -0700


Store the operator index for the first comparision change in the ADFun object.


commit 9caf25014079a60df5de17bcac76775daf8ee201
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 12:45:56 2015 -0700


Make compare_change a parameter (so will be easy to add compare_first).


commit 2246d22fe82b8909d432f82ab0d783ce3351a02f
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 09:12:40 2015 -0700


speed_branch.sh: fix directory (before cd).


commit b3910de86558a97749741bfb728e45c5a86d1c73
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 05:14:01 2015 -0700


search.sh: use git to get list of source files.
ad_fun.hpp: imporve doxygen doc for compare_change_.
ad_tape.hpp: remove RecordCompare? (no longer used).
atomic_base.hpp: minor edit to user documentation.


commit dd74f331386cadc9cc272c264296e575691aa3f8
Author: Brad Bell <bradbell@…>
Date: Thu Jan 15 04:12:34 2015 -0700


Change Leq..Op -> Le..Op and leq.._op -> le.._op.


commit ae729296323eb7f4f4a7c0e90a303a8d7f4ed42a
Author: Brad Bell <bradbell@…>
Date: Wed Jan 14 21:19:55 2015 -0700


comp_op.hpp: Add doxygen documentaiton for compare operators.
compare.hpp: avoid an extra branch.


commit b064d59a5ad01dff5c708cc8c02f628f58c863ec
Author: Brad Bell <bradbell@…>
Date: Wed Jan 14 16:11:25 2015 -0700


Remove ComOp?, replaced by specific cases Eq..Op, Ne..Op, Leq..Op, Lt..Op,
which should run faster.


forward0sweep.hpp: Remove out-dated comment about CompOp?.
forward1sweep.hpp: Remove out-dated comment about CompOp?.


commit 5bb0a70d1151e9086b88024050cea6cf11e83aa7
Author: Brad Bell <bradbell@…>
Date: Wed Jan 14 09:02:17 2015 -0700


Use Eq..Op, Ne..Op to implement == and !=.


commit 0ebeec61bc040a00c50db41ca5da31fb87194f93
Author: Brad Bell <bradbell@…>
Date: Wed Jan 14 07:19:26 2015 -0700


compare.hpp: Finish folding < <= > >= into Lt..Op and Leq..Op.


commit c949e5b72158b98cbab61c8aea98e76008e9c2f4
Author: Brad Bell <bradbell@…>
Date: Wed Jan 14 06:28:41 2015 -0700


  1. Change plan to fold all compare operations into Leq..Op and Lt..Op cases.
  2. Fix alphabetical order between Ld and Leq.


commit 6ffee88b68b682359d62bc75a8c2ba3e28d012ac
Author: Brad Bell <bradbell@…>
Date: Tue Jan 13 22:40:18 2015 -0700


Splite parameter <= variable and parameter > variable as separate compare
operators.


commit 0841014db4ead690d1c2358f5e09494030ae1e5f
Author: Brad Bell <bradbell@…>
Date: Tue Jan 13 21:12:57 2015 -0700


Attempting to recording and playback of compare operators faster:
Split variable <= variable and variable > variable out as separate cases.


compare.hpp: remove white space at end of lines.


commit cfd3afceebd4b3383b3042cbca98caf82ff77670
Author: Brad Bell <bradbell@…>
Date: Tue Jan 13 08:39:12 2015 -0700


speed_branch.sh: compare speed_cppad between two git branches.
speed_new.sh: correct list of options, remove unused variable.
wish_list.omh: correct discussion of effect of keeping compare operators.


commit 9ed03e1ee2c258ca38561ad083067e235d032e14
Author: Brad Bell <bradbell@…>
Date: Tue Jan 13 05:54:30 2015 -0700


Change test so it corresponds to optimization keeping compare operators.


commit 5c418d477d58984b094bba027ebb0794e759e557
Author: Brad Bell <bradbell@…>
Date: Tue Jan 13 05:12:50 2015 -0700


Include example and test of using CompareChange? with optimized tape.


commit c1b48edfa56ca096ce8c2db1dbadb2658cb13fe3
Author: Brad Bell <bradbell@…>
Date: Tue Jan 13 04:24:54 2015 -0700


Fix optimization so variables used in compare opertions are alwasy available.


commit 6fe607ca30db07b356fd3a9fe9779fa2dfd382d8
Author: Brad Bell <bradbell@…>
Date: Mon Jan 12 07:56:41 2015 -0700


Keep CompareChange? in optimized versions of tape and when NDEBUG is defined.

commit b4c0e51489cc0878499a331b4af4875b2781d8d8
Author: Brad Bell <bradbell@…>
Date: Sun Jan 11 17:01:59 2015 -0700

new_branch.sh: final procedure (only hand tested).

  • Property svn:keywords set to Id
File size: 83.0 KB
Line 
1/* $Id: optimize.hpp 3607 2015-01-20 16:20:41Z bradbell $ */
2# ifndef CPPAD_OPTIMIZE_INCLUDED
3# define CPPAD_OPTIMIZE_INCLUDED
4
5/* --------------------------------------------------------------------------
6CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-15 Bradley M. Bell
7
8CppAD is distributed under multiple licenses. This distribution is under
9the terms of the
10                    Eclipse Public License Version 1.0.
11
12A copy of this license is included in the COPYING file of this distribution.
13Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
14-------------------------------------------------------------------------- */
15
16/*
17$begin optimize$$
18$spell
19        jac
20        bool
21        Taylor
22        var
23        CppAD
24        cppad
25        std
26        CondExpEq
27$$
28
29$section Optimize an ADFun Object Tape$$
30
31$index optimize$$
32$index tape, optimize$$
33$index sequence, optimize operations$$
34$index operations, optimize sequence$$
35$index speed, optimize$$
36$index memory, optimize$$
37
38$head Syntax$$
39$icode%f%.optimize()%$$
40
41
42$head Purpose$$
43The operation sequence corresponding to an $cref ADFun$$ object can
44be very large and involve many operations; see the
45size functions in $cref seq_property$$.
46The $icode%f%.optimize%$$ procedure reduces the number of operations,
47and thereby the time and the memory, required to
48compute function and derivative values.
49
50$head f$$
51The object $icode f$$ has prototype
52$codei%
53        ADFun<%Base%> %f%
54%$$
55
56$head Improvements$$
57You can see the reduction in number of variables in the operation sequence
58by calling the function $cref/f.size_var()/seq_property/size_var/$$
59before and after the optimization procedure.
60Given that the optimization procedure takes time,
61it may be faster to skip this optimize procedure and just compute
62derivatives using the original operation sequence.
63
64$subhead Testing$$
65You can run the CppAD $cref/speed/speed_main/$$ tests and see
66the corresponding changes in number of variables and execution time;
67see $cref cmake_check$$.
68
69$head Efficiency$$
70The $code optimize$$ member function
71may greatly reduce the number of variables
72in the operation sequence; see $cref/size_var/seq_property/size_var/$$.
73If a $cref/zero order forward/forward_zero/$$ calculation is done during
74the construction of $icode f$$, it will require more memory
75and time than required after the optimization procedure.
76In addition, it will need to be redone.
77For this reason, it is more efficient to use
78$codei%
79        ADFun<%Base%> %f%;
80        %f%.Dependent(%x%, %y%);
81        %f%.optimize();
82%$$
83instead of
84$codei%
85        ADFun<%Base%> %f%(%x%, %y%)
86        %f%.optimize();
87%$$
88See the discussion about
89$cref/sequence constructors/FunConstruct/Sequence Constructor/$$.
90
91$head Atomic Functions$$
92There are some subtitle issue with optimized $cref atomic$$ functions
93$latex v = g(u)$$:
94
95$subhead rev_sparse_jac$$
96The $cref atomic_rev_sparse_jac$$ function is be used to determine
97which components of $icode u$$ affect the dependent variables of $icode f$$.
98The current setting of the
99$cref/atomic_sparsity/atomic_option/atomic_sparsity/$$ pattern for each
100atomic function is used to determine if the $code bool$$ or
101$code std::set<size_t>$$ version of $cref atomic_rev_sparse_jac$$ is used.
102
103$subhead nan$$
104If $icode%u%[%i%]%$$ does not affect the value of
105the dependent variables for $icode f$$,
106the value of $icode%u%[%i%]%$$ is set to $cref nan$$.
107
108
109$head Checking Optimization$$
110$index NDEBUG$$
111If $cref/NDEBUG/Faq/Speed/NDEBUG/$$ is not defined,
112and $cref/f.size_order()/size_order/$$ is greater than zero,
113a $cref forward_zero$$ calculation is done using the optimized version
114of $icode f$$ and the results are checked to see that they are
115the same as before.
116If they are not the same, the
117$cref ErrorHandler$$ is called with a known error message
118related to $icode%f%.optimize()%$$.
119
120$head Example$$
121$children%
122        example/optimize.cpp
123%$$
124The file
125$cref optimize.cpp$$
126contains an example and test of this operation.
127It returns true if it succeeds and false otherwise.
128
129$end
130-----------------------------------------------------------------------------
131*/
132# include <stack>
133
134namespace CppAD { // BEGIN_CPPAD_NAMESPACE
135namespace optimize { // BEGIN_CPPAD_OPTIMIZE_NAMESPACE
136/*!
137\file optimize.hpp
138Routines for optimizing a tape
139*/
140
141
142/*!
143State for this variable set during reverse sweep.
144*/
145enum enum_connect_type {
146        /// There is no operation that connects this variable to the
147        /// independent variables.
148        not_connected        ,
149
150        /// There is one or more operations that connects this variable to the
151        /// independent variables.
152        yes_connected        ,
153
154        /// There is only one parrent that connects this variable to the
155        /// independent variables and the parent is a summation operation; i.e.,
156        /// AddvvOp, AddpvOp, SubpvOp, SubvpOp, or SubvvOp.
157        sum_connected        ,
158
159        /// Satisfies the sum_connected assumptions above and in addition
160        /// this variable is the result of summation operator.
161        csum_connected       ,
162
163        /// This node is only connected in the case where the comparision is
164        /// true for the conditional expression with index \c connect_index.
165        cexp_connected
166
167};
168
169/*!
170Class used to hold information about one conditional expression.
171*/
172class class_cexp_pair {
173public:
174        /// If this is true (false) this connection is only for the case where
175        /// the comparision in the conditional expression is true (false)
176        bool compare_; 
177        /// This is the index of the conditional expression (in cksip_info)
178        /// for this connection
179        size_t index_;
180        /// constructor
181        class_cexp_pair(const bool& compare, const size_t& index)
182        : compare_(compare), index_(index)
183        { }
184        /// assignment operator
185        void operator=(const class_cexp_pair& right)
186        {       index_   = right.index_;
187                compare_ = right.compare_;
188        }
189        /// not equal operator
190        bool operator!=(const class_cexp_pair& right)
191        {       return (index_ != right.index_) | (compare_ != right.compare_);
192        }
193        /// Less than operator
194        /// (required for intersection of two sets of class_cexp_pair elements).
195        bool operator<(const class_cexp_pair& right) const
196        {       if( index_ == right.index_ )
197                        return size_t(compare_) < size_t(right.compare_);       
198                return index_ < right.index_;
199        }
200};
201/*!
202Compute intersection of two sets of class_cexp_pair elements.
203
204\param left
205first operand of the intersection
206
207\param right
208second operand of the intersection
209
210\result
211the intersection of left and right
212*/
213inline std::set<class_cexp_pair> intersection(
214        std::set<class_cexp_pair>& left  ,
215        std::set<class_cexp_pair>& right )
216{       std::set<class_cexp_pair> result;
217        std::set_intersection(
218                left.begin()  ,
219                left.end()    ,
220                right.begin() ,
221                right.end()   ,
222                std::inserter(result, result.begin())
223        );
224        return result; 
225}
226
227/*!
228Structure used by \c optimize to hold information about one variable.
229in the old operation seqeunce.
230*/
231struct struct_old_variable {
232        /// Operator for which this variable is the result, \c NumRes(op) > 0.
233        /// Set by the reverse sweep at beginning of optimization.
234        OpCode              op;       
235
236        /// Pointer to first argument (child) for this operator.
237        /// Set by the reverse sweep at beginning of optimization.
238        const addr_t*       arg;
239
240        /// How is this variable connected to the independent variables
241        enum_connect_type connect_type; 
242
243        /*!
244        If \c connect_type is \c cexp_connected,
245        this is the corresponding infromation for the conditional connections.
246        */
247        std::set<class_cexp_pair> cexp_set;
248
249        /// New operation sequence corresponding to this old varable.
250        /// Set during forward sweep to the index in the new tape
251        addr_t new_var;
252
253        /// New operator index for this varable.
254        /// Set during forward sweep to the index in the new tape
255        size_t new_op;
256
257        /// Did this variable match another variable in the operation sequence
258        bool match;
259};
260
261struct struct_size_pair {
262        size_t i_op;  // an operator index
263        size_t i_var; // a variable index
264};
265
266/*!
267Structures used by \c record_csum
268to hold information about one variable.
269*/
270struct struct_csum_variable {
271        /// Operator for which this variable is the result, \c NumRes(op) > 0.
272        OpCode              op;       
273
274        /// Pointer to first argument (child) for this operator.
275        /// Set by the reverse sweep at beginning of optimization.
276        const addr_t*       arg;
277
278        /// Is this variable added to the summation
279        /// (if not it is subtracted)
280        bool                add;
281};
282
283/*!
284Structure used to pass work space from \c optimize to \c record_csum
285(so that stacks do not start from zero size every time).
286*/
287struct struct_csum_stacks {
288        /// stack of operations in the cummulative summation
289        std::stack<struct struct_csum_variable>   op_stack;
290        /// stack of variables to be added
291        std::stack<size_t >                         add_stack;
292        /// stack of variables to be subtracted
293        std::stack<size_t >                         sub_stack;
294};
295
296/*!
297CExpOp information that is copied to corresponding CSkipOp
298*/
299struct struct_cskip_info {
300        /// comparision operator
301        CompareOp cop;
302        /// (flag & 1) is true if and only if left is a variable
303        /// (flag & 2) is true if and only if right is a variable
304        size_t flag;
305        /// index for left comparison operand
306        size_t left; 
307        /// index for right comparison operand
308        size_t right; 
309        /// maximum variable index between left and right
310        size_t max_left_right;
311        /// set of variables to skip on true
312        CppAD::vector<size_t> skip_var_true;
313        /// set of variables to skip on false
314        CppAD::vector<size_t> skip_var_false;
315        /// set of operations to skip on true
316        CppAD::vector<size_t> skip_op_true;
317        /// set of operations to skip on false
318        CppAD::vector<size_t> skip_op_false;
319        /// size of skip_op_true
320        size_t n_op_true;
321        /// size of skip_op_false
322        size_t n_op_false;
323        /// index in the argument recording of first argument for this CSkipOp
324        size_t i_arg;
325};
326/*!
327Connection information for a user atomic function
328*/
329struct struct_user_info {
330        /// type of connection for this atomic function
331        enum_connect_type connect_type;
332        /// If this is an conditional connection, this is the information
333        /// of the correpsonding CondExpOp operators
334        std::set<class_cexp_pair> cexp_set;
335        /// If this is a conditional connection, this is the operator
336        /// index of the beginning of the atomic call sequence; i.e.,
337        /// the first UserOp.
338        size_t op_begin;
339        /// If this is a conditional connection, this is one more than the
340        ///  operator index of the ending of the atomic call sequence; i.e.,
341        /// the second UserOp.
342        size_t op_end;
343};
344
345/*!
346Shared documentation for optimization helper functions (not called).
347
348<!-- define prototype -->
349\param tape
350is a vector that maps a variable index, in the old operation sequence,
351to an <tt>struct_old_variable</tt> information record.
352Note that the index for this vector must be greater than or equal zero and
353less than <tt>tape.size()</tt>.
354
355\li <tt>tape[i].op</tt>
356is the operator in the old operation sequence
357corresponding to the old variable index \c i.
358Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
359
360\li <tt>tape[i].arg</tt>
361for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
362is the j-th the argument, in the old operation sequence,
363corresponding to the old variable index \c i.
364Assertion: <tt>tape[i].arg[j] < i</tt>.
365
366\li <tt>tape[i].new_var</tt>
367Suppose
368<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
369and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
370It follows that <tt>tape[k].new_var</tt>
371has alread been set to the variable in the new operation sequence
372corresponding to the old variable index \c k.
373This means that the \c new_var value has been set
374for all the possible arguments that come before \a current.
375
376\param current
377is the index in the old operation sequence for
378the variable corresponding to the result for the current operator.
379Assertions:
380<tt>
381current < tape.size(),
382NumRes( tape[current].op ) > 0.
383</tt>
384
385\param npar
386is the number of parameters corresponding to this operation sequence.
387
388\param par
389is a vector of length \a npar containing the parameters
390for this operation sequence; i.e.,
391given a parameter index \c i, the corresponding parameter value is
392<tt>par[i]</tt>.
393<!-- end prototype -->
394*/
395template <class Base>
396void prototype(
397        const CppAD::vector<struct struct_old_variable>& tape           ,
398        size_t                                             current        ,
399        size_t                                             npar           ,
400        const Base*                                        par            )
401{       CPPAD_ASSERT_UNKNOWN(false); }
402
403/*!
404Check a unary operator for a complete match with a previous operator.
405
406A complete match means that the result of the previous operator
407can be used inplace of the result for current operator.
408
409<!-- replace prototype -->
410\param tape
411is a vector that maps a variable index, in the old operation sequence,
412to an <tt>struct_old_variable</tt> information record.
413Note that the index for this vector must be greater than or equal zero and
414less than <tt>tape.size()</tt>.
415
416\li <tt>tape[i].op</tt>
417is the operator in the old operation sequence
418corresponding to the old variable index \c i.
419Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
420
421\li <tt>tape[i].arg</tt>
422for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
423is the j-th the argument, in the old operation sequence,
424corresponding to the old variable index \c i.
425Assertion: <tt>tape[i].arg[j] < i</tt>.
426
427\li <tt>tape[i].new_var</tt>
428Suppose
429<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
430and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
431It follows that <tt>tape[k].new_var</tt>
432has alread been set to the variable in the new operation sequence
433corresponding to the old variable index \c k.
434This means that the \c new_var value has been set
435for all the possible arguments that come before \a current.
436
437\param current
438is the index in the old operation sequence for
439the variable corresponding to the result for the current operator.
440Assertions:
441<tt>
442current < tape.size(),
443NumRes( tape[current].op ) > 0.
444</tt>
445
446\param npar
447is the number of parameters corresponding to this operation sequence.
448
449\param par
450is a vector of length \a npar containing the parameters
451for this operation sequence; i.e.,
452given a parameter index \c i, the corresponding parameter value is
453<tt>par[i]</tt>.
454<!-- end prototype -->
455
456\param hash_table_var
457is a vector with size CPPAD_HASH_TABLE_SIZE
458that maps a hash code to the corresponding
459variable index in the old operation sequence.
460All the values in this table must be less than \a current.
461
462\param code
463The input value of code does not matter.
464The output value of code is the hash code corresponding to
465this operation in the new operation sequence.
466
467\return
468If the return value is zero,
469no match was found.
470If the return value is greater than zero,
471it is the old operation sequence index of a variable,
472that comes before current and can be used to replace the current variable.
473
474\par Restrictions:
475NumArg( tape[current].op ) == 1
476*/
477template <class Base>
478addr_t unary_match(
479        const CppAD::vector<struct struct_old_variable>& tape           ,
480        size_t                                             current        ,
481        size_t                                             npar           ,
482        const Base*                                        par            ,
483        const CppAD::vector<size_t>&                       hash_table_var ,
484        unsigned short&                                    code           )
485{       const addr_t* arg = tape[current].arg;
486        OpCode        op  = tape[current].op;
487        addr_t new_arg[1];
488       
489        CPPAD_ASSERT_UNKNOWN( NumArg(op) == 1 );
490        CPPAD_ASSERT_UNKNOWN( NumRes(op) > 0  );
491        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < current );
492        new_arg[0] = tape[arg[0]].new_var;
493        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < current );
494        code = hash_code(
495                op                  , 
496                new_arg             ,
497                npar                ,
498                par
499        );
500        size_t  i_var  = hash_table_var[code];
501        CPPAD_ASSERT_UNKNOWN( i_var < current );
502        if( op == tape[i_var].op )
503        {       size_t k = tape[i_var].arg[0];
504                CPPAD_ASSERT_UNKNOWN( k < i_var );
505                if (new_arg[0] == tape[k].new_var )
506                        return i_var;
507        }
508        return 0;
509} 
510
511/*!
512Check a binary operator for a complete match with a previous operator,
513
514<!-- replace prototype -->
515\param tape
516is a vector that maps a variable index, in the old operation sequence,
517to an <tt>struct_old_variable</tt> information record.
518Note that the index for this vector must be greater than or equal zero and
519less than <tt>tape.size()</tt>.
520
521\li <tt>tape[i].op</tt>
522is the operator in the old operation sequence
523corresponding to the old variable index \c i.
524Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
525
526\li <tt>tape[i].arg</tt>
527for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
528is the j-th the argument, in the old operation sequence,
529corresponding to the old variable index \c i.
530Assertion: <tt>tape[i].arg[j] < i</tt>.
531
532\li <tt>tape[i].new_var</tt>
533Suppose
534<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
535and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
536It follows that <tt>tape[k].new_var</tt>
537has alread been set to the variable in the new operation sequence
538corresponding to the old variable index \c k.
539This means that the \c new_var value has been set
540for all the possible arguments that come before \a current.
541
542\param current
543is the index in the old operation sequence for
544the variable corresponding to the result for the current operator.
545Assertions:
546<tt>
547current < tape.size(),
548NumRes( tape[current].op ) > 0.
549</tt>
550
551\param npar
552is the number of parameters corresponding to this operation sequence.
553
554\param par
555is a vector of length \a npar containing the parameters
556for this operation sequence; i.e.,
557given a parameter index \c i, the corresponding parameter value is
558<tt>par[i]</tt>.
559<!-- end prototype -->
560
561\param hash_table_var
562is a vector with size CPPAD_HASH_TABLE_SIZE
563that maps a hash code to the corresponding
564variable index in the old operation sequence.
565All the values in this table must be less than \a current.
566
567\param code
568The input value of code does not matter.
569The output value of code is the hash code corresponding to
570this operation in the new operation sequence.
571
572\return
573If the return value is zero,
574no match was found.
575If the return value is greater than zero,
576it is the index of a new variable that can be used to replace the
577old variable.
578
579
580\par Restrictions:
581The binary operator must be an addition, subtraction, multiplication, division
582or power operator.  NumArg( tape[current].op ) == 1.
583*/
584template <class Base>
585inline addr_t binary_match(
586        const CppAD::vector<struct struct_old_variable>&   tape           ,
587        size_t                                             current        ,
588        size_t                                             npar           ,
589        const Base*                                        par            ,
590        const CppAD::vector<size_t>&                       hash_table_var ,
591        unsigned short&                                    code           )
592{       OpCode        op         = tape[current].op;
593        const addr_t* arg        = tape[current].arg;
594        addr_t        new_arg[2];
595        bool          parameter[2];
596
597        // initialize return value
598        addr_t  match_var = 0;
599
600        CPPAD_ASSERT_UNKNOWN( NumArg(op) == 2 );
601        CPPAD_ASSERT_UNKNOWN( NumRes(op) >  0 );
602        switch(op)
603        {       // index op variable
604                case DisOp:
605                // parameter not defined for this case
606                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < current );
607                new_arg[0]   = arg[0];
608                new_arg[1]   = tape[arg[1]].new_var;
609                break;
610
611                // parameter op variable ----------------------------------
612                case AddpvOp:
613                case MulpvOp:
614                case DivpvOp:
615                case PowpvOp:
616                case SubpvOp:
617                // arg[0]
618                parameter[0] = true;
619                new_arg[0]   = arg[0];
620                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < npar );
621                // arg[1]
622                parameter[1] = false;
623                new_arg[1]   = tape[arg[1]].new_var;
624                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < current );
625                break;
626
627                // variable op parameter -----------------------------------
628                case DivvpOp:
629                case PowvpOp:
630                case SubvpOp:
631                // arg[0]
632                parameter[0] = false;
633                new_arg[0]   = tape[arg[0]].new_var;
634                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < current );
635                // arg[1]
636                parameter[1] = true;
637                new_arg[1]   = arg[1];
638                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < npar );
639                break;
640
641                // variable op variable -----------------------------------
642                case AddvvOp:
643                case MulvvOp:
644                case DivvvOp:
645                case PowvvOp:
646                case SubvvOp:
647                // arg[0]
648                parameter[0] = false;
649                new_arg[0]   = tape[arg[0]].new_var;
650                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < current );
651                // arg[1]
652                parameter[1] = false;
653                new_arg[1]   = tape[arg[1]].new_var;
654                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < current );
655                break;
656
657                // must be one of the cases above
658                default:
659                CPPAD_ASSERT_UNKNOWN(false);
660        }
661        code = hash_code(
662                op                  , 
663                new_arg             ,
664                npar                ,
665                par
666        );
667        size_t  i_var  = hash_table_var[code];
668        CPPAD_ASSERT_UNKNOWN( i_var < current );
669        if( op == tape[i_var].op )
670        {       bool match = true;
671                if( op == DisOp )
672                {       match   &= new_arg[0] == tape[i_var].arg[0];
673                        size_t k = tape[i_var].arg[1];
674                        match   &= new_arg[1] == tape[k].new_var;
675                }
676                else
677                {       for(size_t j = 0; j < 2; j++)
678                        {       size_t k = tape[i_var].arg[j];
679                                if( parameter[j] )
680                                {       CPPAD_ASSERT_UNKNOWN( k < npar );
681                                        match &= IdenticalEqualPar(
682                                                par[ arg[j] ], par[k]
683                                        );
684                                }
685                                else
686                                {       CPPAD_ASSERT_UNKNOWN( k < i_var );
687                                        match &= (new_arg[j] == tape[k].new_var);
688                                }
689                        }
690                }
691                if( match )
692                        match_var = i_var;
693        }
694        if( (match_var > 0) | ( (op != AddvvOp) & (op != MulvvOp ) ) )
695                return match_var;
696
697        // check for match with argument order switched ----------------------
698        CPPAD_ASSERT_UNKNOWN( op == AddvvOp || op == MulvvOp );
699        std::swap(new_arg[0], new_arg[1]);
700        unsigned short code_switch = hash_code(
701                op                  , 
702                new_arg             ,
703                npar                ,
704                par
705        );
706        i_var  = hash_table_var[code_switch];
707        CPPAD_ASSERT_UNKNOWN( i_var < current );
708        if( op == tape[i_var].op )
709        {       bool match = true;
710                size_t j;
711                for(j = 0; j < 2; j++)
712                {       size_t k = tape[i_var].arg[j];
713                        CPPAD_ASSERT_UNKNOWN( k < i_var );
714                        match &= (new_arg[j] == tape[k].new_var);
715                }
716                if( match )
717                        match_var = i_var;
718        }
719        return match_var;
720} 
721
722/*!
723Record an operation of the form (parameter op variable).
724
725<!-- replace prototype -->
726\param tape
727is a vector that maps a variable index, in the old operation sequence,
728to an <tt>struct_old_variable</tt> information record.
729Note that the index for this vector must be greater than or equal zero and
730less than <tt>tape.size()</tt>.
731
732\li <tt>tape[i].op</tt>
733is the operator in the old operation sequence
734corresponding to the old variable index \c i.
735Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
736
737\li <tt>tape[i].arg</tt>
738for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
739is the j-th the argument, in the old operation sequence,
740corresponding to the old variable index \c i.
741Assertion: <tt>tape[i].arg[j] < i</tt>.
742
743\li <tt>tape[i].new_var</tt>
744Suppose
745<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
746and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
747It follows that <tt>tape[k].new_var</tt>
748has alread been set to the variable in the new operation sequence
749corresponding to the old variable index \c k.
750This means that the \c new_var value has been set
751for all the possible arguments that come before \a current.
752
753\param current
754is the index in the old operation sequence for
755the variable corresponding to the result for the current operator.
756Assertions:
757<tt>
758current < tape.size(),
759NumRes( tape[current].op ) > 0.
760</tt>
761
762\param npar
763is the number of parameters corresponding to this operation sequence.
764
765\param par
766is a vector of length \a npar containing the parameters
767for this operation sequence; i.e.,
768given a parameter index \c i, the corresponding parameter value is
769<tt>par[i]</tt>.
770<!-- end prototype -->
771
772\param rec
773is the object that will record the operations.
774
775\param op
776is the operator that we are recording which must be one of the following:
777AddpvOp, DivpvOp, MulpvOp, PowvpOp, SubpvOp.
778 
779\param arg
780is the vector of arguments for this operator.
781
782\return
783the result is the operaiton and variable index corresponding to the current
784operation in the new operation sequence.
785*/
786template <class Base>
787struct_size_pair record_pv(
788        const CppAD::vector<struct struct_old_variable>& tape           ,
789        size_t                                             current        ,
790        size_t                                             npar           ,
791        const Base*                                        par            ,
792        recorder<Base>*                                    rec            ,
793        OpCode                                             op             ,
794        const addr_t*                                      arg            )
795{
796# ifndef NDEBUG
797        switch(op)
798        {       case AddpvOp:
799                case DivpvOp:
800                case MulpvOp:
801                case PowpvOp:
802                case SubpvOp:
803                break;
804
805                default:
806                CPPAD_ASSERT_UNKNOWN(false);
807        }
808# endif
809        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < npar    );
810        CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < current );
811        addr_t new_arg[2];
812        new_arg[0]   = rec->PutPar( par[arg[0]] );
813        new_arg[1]   = tape[ arg[1] ].new_var;
814        rec->PutArg( new_arg[0], new_arg[1] );
815
816        struct_size_pair ret;
817        ret.i_op  = rec->num_op_rec();
818        ret.i_var = rec->PutOp(op);
819        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < ret.i_var );
820        return ret;
821}
822
823
824/*!
825Record an operation of the form (variable op parameter).
826
827<!-- replace prototype -->
828\param tape
829is a vector that maps a variable index, in the old operation sequence,
830to an <tt>struct_old_variable</tt> information record.
831Note that the index for this vector must be greater than or equal zero and
832less than <tt>tape.size()</tt>.
833
834\li <tt>tape[i].op</tt>
835is the operator in the old operation sequence
836corresponding to the old variable index \c i.
837Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
838
839\li <tt>tape[i].arg</tt>
840for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
841is the j-th the argument, in the old operation sequence,
842corresponding to the old variable index \c i.
843Assertion: <tt>tape[i].arg[j] < i</tt>.
844
845\li <tt>tape[i].new_var</tt>
846Suppose
847<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
848and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
849It follows that <tt>tape[k].new_var</tt>
850has alread been set to the variable in the new operation sequence
851corresponding to the old variable index \c k.
852This means that the \c new_var value has been set
853for all the possible arguments that come before \a current.
854
855\param current
856is the index in the old operation sequence for
857the variable corresponding to the result for the current operator.
858Assertions:
859<tt>
860current < tape.size(),
861NumRes( tape[current].op ) > 0.
862</tt>
863
864\param npar
865is the number of parameters corresponding to this operation sequence.
866
867\param par
868is a vector of length \a npar containing the parameters
869for this operation sequence; i.e.,
870given a parameter index \c i, the corresponding parameter value is
871<tt>par[i]</tt>.
872<!-- end prototype -->
873
874\param rec
875is the object that will record the operations.
876
877\param op
878is the operator that we are recording which must be one of the following:
879DivvpOp, PowvpOp, SubvpOp.
880 
881\param arg
882is the vector of arguments for this operator.
883
884\return
885the result operation and variable index corresponding to the current
886operation in the new operation sequence.
887*/
888template <class Base>
889struct_size_pair record_vp(
890        const CppAD::vector<struct struct_old_variable>& tape           ,
891        size_t                                             current        ,
892        size_t                                             npar           ,
893        const Base*                                        par            ,
894        recorder<Base>*                                    rec            ,
895        OpCode                                             op             ,
896        const addr_t*                                      arg            )
897{
898# ifndef NDEBUG
899        switch(op)
900        {       case DivvpOp:
901                case PowvpOp:
902                case SubvpOp:
903                break;
904
905                default:
906                CPPAD_ASSERT_UNKNOWN(false);
907        }
908# endif
909        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < current );
910        CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < npar    );
911        addr_t new_arg[2];
912        new_arg[0]   = tape[ arg[0] ].new_var;
913        new_arg[1]   = rec->PutPar( par[arg[1]] );
914        rec->PutArg( new_arg[0], new_arg[1] );
915
916        struct_size_pair ret;
917        ret.i_op  = rec->num_op_rec();
918        ret.i_var = rec->PutOp(op);
919        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < ret.i_var );
920        return ret;
921}
922
923/*!
924Record an operation of the form (variable op variable).
925
926<!-- replace prototype -->
927\param tape
928is a vector that maps a variable index, in the old operation sequence,
929to an <tt>struct_old_variable</tt> information record.
930Note that the index for this vector must be greater than or equal zero and
931less than <tt>tape.size()</tt>.
932
933\li <tt>tape[i].op</tt>
934is the operator in the old operation sequence
935corresponding to the old variable index \c i.
936Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
937
938\li <tt>tape[i].arg</tt>
939for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
940is the j-th the argument, in the old operation sequence,
941corresponding to the old variable index \c i.
942Assertion: <tt>tape[i].arg[j] < i</tt>.
943
944\li <tt>tape[i].new_var</tt>
945Suppose
946<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
947and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
948It follows that <tt>tape[k].new_var</tt>
949has alread been set to the variable in the new operation sequence
950corresponding to the old variable index \c k.
951This means that the \c new_var value has been set
952for all the possible arguments that come before \a current.
953
954\param current
955is the index in the old operation sequence for
956the variable corresponding to the result for the current operator.
957Assertions:
958<tt>
959current < tape.size(),
960NumRes( tape[current].op ) > 0.
961</tt>
962
963\param npar
964is the number of parameters corresponding to this operation sequence.
965
966\param par
967is a vector of length \a npar containing the parameters
968for this operation sequence; i.e.,
969given a parameter index \c i, the corresponding parameter value is
970<tt>par[i]</tt>.
971<!-- end prototype -->
972
973\param rec
974is the object that will record the operations.
975
976\param op
977is the operator that we are recording which must be one of the following:
978AddvvOp, DivvvOp, MulvvOp, PowvpOp, SubvvOp.
979 
980\param arg
981is the vector of arguments for this operator.
982
983\return
984the result is the operation and variable index corresponding to the current
985operation in the new operation sequence.
986*/
987template <class Base>
988struct_size_pair record_vv(
989        const CppAD::vector<struct struct_old_variable>& tape           ,
990        size_t                                             current        ,
991        size_t                                             npar           ,
992        const Base*                                        par            ,
993        recorder<Base>*                                    rec            ,
994        OpCode                                             op             ,
995        const addr_t*                                      arg            )
996{
997# ifndef NDEBUG
998        switch(op)
999        {       case AddvvOp:
1000                case DivvvOp:
1001                case MulvvOp:
1002                case PowvvOp:
1003                case SubvvOp:
1004                break;
1005
1006                default:
1007                CPPAD_ASSERT_UNKNOWN(false);
1008        }
1009# endif
1010        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < current );
1011        CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < current );
1012        addr_t new_arg[2];
1013        new_arg[0]   = tape[ arg[0] ].new_var;
1014        new_arg[1]   = tape[ arg[1] ].new_var;
1015        rec->PutArg( new_arg[0], new_arg[1] );
1016
1017        struct_size_pair ret;
1018        ret.i_op  = rec->num_op_rec();
1019        ret.i_var = rec->PutOp(op);
1020        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < ret.i_var );
1021        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[1]) < ret.i_var );
1022        return ret;
1023}
1024
1025// ==========================================================================
1026
1027/*!
1028Recording a cummulative cummulative summation starting at its highest parrent.
1029
1030<!-- replace prototype -->
1031\param tape
1032is a vector that maps a variable index, in the old operation sequence,
1033to an <tt>struct_old_variable</tt> information record.
1034Note that the index for this vector must be greater than or equal zero and
1035less than <tt>tape.size()</tt>.
1036
1037\li <tt>tape[i].op</tt>
1038is the operator in the old operation sequence
1039corresponding to the old variable index \c i.
1040Assertion: <tt>NumRes(tape[i].op) > 0</tt>.
1041
1042\li <tt>tape[i].arg</tt>
1043for <tt>j < NumArg( tape[i].op ), tape[i].arg[j]</tt>
1044is the j-th the argument, in the old operation sequence,
1045corresponding to the old variable index \c i.
1046Assertion: <tt>tape[i].arg[j] < i</tt>.
1047
1048\li <tt>tape[i].new_var</tt>
1049Suppose
1050<tt>i <= current, j < NumArg( tape[i].op ), and k = tape[i].arg[j]</tt>,
1051and \c j corresponds to a variable for operator <tt>tape[i].op</tt>.
1052It follows that <tt>tape[k].new_var</tt>
1053has alread been set to the variable in the new operation sequence
1054corresponding to the old variable index \c k.
1055This means that the \c new_var value has been set
1056for all the possible arguments that come before \a current.
1057
1058\param current
1059is the index in the old operation sequence for
1060the variable corresponding to the result for the current operator.
1061Assertions:
1062<tt>
1063current < tape.size(),
1064NumRes( tape[current].op ) > 0.
1065</tt>
1066
1067\param npar
1068is the number of parameters corresponding to this operation sequence.
1069
1070\param par
1071is a vector of length \a npar containing the parameters
1072for this operation sequence; i.e.,
1073given a parameter index \c i, the corresponding parameter value is
1074<tt>par[i]</tt>.
1075<!-- end prototype -->
1076
1077\param rec
1078is the object that will record the operations.
1079
1080\param work
1081Is used for computation. On input and output,
1082<tt>work.op_stack.empty()</tt>,
1083<tt>work.add_stack.empty()</tt>, and
1084<tt>work.sub_stack.empty()</tt>,
1085are all true true.
1086These stacks are passed in so that elements can be allocated once
1087and then the elements can be reused with calls to \c record_csum.
1088
1089\par Exception
1090<tt>tape[i].new_var</tt>
1091is not yet defined for any node \c i that is \c csum_connected
1092to the \a current node
1093(or that is \c sum_connected to a node that is \c csum_connected).
1094For example; suppose that index \c j corresponds to a variable
1095in the current operator,
1096<tt>i = tape[current].arg[j]</tt>,
1097and
1098<tt>tape[arg[j]].connect_type == csum_connected</tt>.
1099It then follows that
1100<tt>tape[i].new_var == tape.size()</tt>.
1101
1102\par Restriction:
1103\li <tt>tape[current].op</tt>
1104must be one of <tt>AddpvOp, AddvvOp, SubpvOp, SubvpOp, SubvvOp</tt>.
1105
1106\li <tt>tape[current].connect_type</tt> must be \c yes_connected.
1107
1108\li <tt>tape[j].connect_type == csum_connected</tt> for some index
1109j that is a variable operand for the current operation.
1110*/
1111
1112
1113template <class Base>
1114struct_size_pair record_csum(
1115        const CppAD::vector<struct struct_old_variable>& tape           ,
1116        size_t                                             current        ,
1117        size_t                                             npar           ,
1118        const Base*                                        par            ,
1119        recorder<Base>*                                    rec            ,
1120        struct_csum_stacks&                              work           )
1121{
1122       
1123        CPPAD_ASSERT_UNKNOWN( work.op_stack.empty() );
1124        CPPAD_ASSERT_UNKNOWN( work.add_stack.empty() );
1125        CPPAD_ASSERT_UNKNOWN( work.sub_stack.empty() );
1126        CPPAD_ASSERT_UNKNOWN( tape[current].connect_type == yes_connected );
1127
1128        size_t                        i;
1129        OpCode                        op;
1130        const addr_t*                 arg;
1131        bool                          add;
1132        struct struct_csum_variable var;
1133
1134        var.op  = tape[current].op;
1135        var.arg = tape[current].arg;
1136        var.add = true; 
1137        work.op_stack.push( var );
1138        Base sum_par(0);
1139
1140# ifndef NDEBUG
1141        bool ok = false;
1142        if( var.op == SubvpOp )
1143                ok = tape[ tape[current].arg[0] ].connect_type == csum_connected;
1144        if( var.op == AddpvOp || var.op == SubpvOp )
1145                ok = tape[ tape[current].arg[1] ].connect_type == csum_connected;
1146        if( var.op == AddvvOp || var.op == SubvvOp )
1147        {       ok  = tape[ tape[current].arg[0] ].connect_type == csum_connected;
1148                ok |= tape[ tape[current].arg[1] ].connect_type == csum_connected;
1149        }
1150        CPPAD_ASSERT_UNKNOWN( ok );
1151# endif
1152        while( ! work.op_stack.empty() )
1153        {       var     = work.op_stack.top();
1154                work.op_stack.pop();
1155                op      = var.op;
1156                arg     = var.arg;
1157                add     = var.add;
1158                // process first argument to this operator
1159                switch(op)
1160                {       case AddpvOp:
1161                        case SubpvOp:
1162                        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < npar );
1163                        if( add )
1164                                sum_par += par[arg[0]];
1165                        else    sum_par -= par[arg[0]];
1166                        break;
1167
1168                        case AddvvOp:
1169                        case SubvpOp:
1170                        case SubvvOp:
1171                        if( tape[arg[0]].connect_type == csum_connected )
1172                        {       CPPAD_ASSERT_UNKNOWN(
1173                                        size_t(tape[arg[0]].new_var) == tape.size()
1174                                );
1175                                var.op  = tape[arg[0]].op;
1176                                var.arg = tape[arg[0]].arg;
1177                                var.add = add; 
1178                                work.op_stack.push( var );
1179                        }
1180                        else if( add )
1181                                work.add_stack.push(arg[0]);
1182                        else    work.sub_stack.push(arg[0]);
1183                        break;
1184
1185                        default:
1186                        CPPAD_ASSERT_UNKNOWN(false);
1187                }
1188                // process second argument to this operator
1189                switch(op)
1190                {
1191                        case SubvpOp:
1192                        CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < npar );
1193                        if( add )
1194                                sum_par -= par[arg[1]];
1195                        else    sum_par += par[arg[1]];
1196                        break;
1197
1198                        case SubvvOp:
1199                        case SubpvOp:
1200                        add = ! add;
1201
1202                        case AddvvOp:
1203                        case AddpvOp:
1204                        if( tape[arg[1]].connect_type == csum_connected )
1205                        {       CPPAD_ASSERT_UNKNOWN(
1206                                        size_t(tape[arg[1]].new_var) == tape.size()
1207                                );
1208                                var.op   = tape[arg[1]].op;
1209                                var.arg  = tape[arg[1]].arg;
1210                                var.add  = add;
1211                                work.op_stack.push( var );
1212                        }
1213                        else if( add )
1214                                work.add_stack.push(arg[1]);
1215                        else    work.sub_stack.push(arg[1]);
1216                        break;
1217
1218                        default:
1219                        CPPAD_ASSERT_UNKNOWN(false);
1220                }
1221        }
1222        // number of variables in this cummulative sum operator
1223        size_t n_add = work.add_stack.size();
1224        size_t n_sub = work.sub_stack.size();
1225        size_t old_arg, new_arg;
1226        rec->PutArg(n_add);                // arg[0]
1227        rec->PutArg(n_sub);                // arg[1]
1228        new_arg = rec->PutPar( sum_par );
1229        rec->PutArg(new_arg);              // arg[2]
1230        for(i = 0; i < n_add; i++)
1231        {       CPPAD_ASSERT_UNKNOWN( ! work.add_stack.empty() );
1232                old_arg = work.add_stack.top();
1233                new_arg = tape[old_arg].new_var;
1234                CPPAD_ASSERT_UNKNOWN( new_arg < tape.size() );
1235                rec->PutArg(new_arg);      // arg[3+i]
1236                work.add_stack.pop();
1237        }
1238        for(i = 0; i < n_sub; i++)
1239        {       CPPAD_ASSERT_UNKNOWN( ! work.sub_stack.empty() );
1240                old_arg = work.sub_stack.top();
1241                new_arg = tape[old_arg].new_var;
1242                CPPAD_ASSERT_UNKNOWN( new_arg < tape.size() );
1243                rec->PutArg(new_arg);      // arg[3 + arg[0] + i]
1244                work.sub_stack.pop();
1245        }
1246        rec->PutArg(n_add + n_sub);        // arg[3 + arg[0] + arg[1]]
1247
1248
1249        struct_size_pair ret;
1250        ret.i_op  = rec->num_op_rec();
1251        ret.i_var = rec->PutOp(CSumOp);
1252        CPPAD_ASSERT_UNKNOWN( new_arg < ret.i_var );
1253        return ret;
1254}
1255// ==========================================================================
1256/*!
1257Convert a player object to an optimized recorder object
1258
1259\tparam Base
1260base type for the operator; i.e., this operation was recorded
1261using AD< \a Base > and computations by this routine are done using type
1262\a Base.
1263
1264\param n
1265is the number of independent variables on the tape.
1266
1267\param dep_taddr
1268On input this vector contains the indices for each of the dependent
1269variable values in the operation sequence corresponding to \a play.
1270Upon return it contains the indices for the same variables but in
1271the operation sequence corresponding to \a rec.
1272
1273\param play
1274This is the operation sequence that we are optimizing.
1275It is essentially const, except for play back state which
1276changes while it plays back the operation seqeunce.
1277
1278\param rec
1279The input contents of this recording does not matter.
1280Upon return, it contains an optimized verison of the
1281operation sequence corresponding to \a play.
1282*/
1283
1284template <class Base>
1285void optimize_run(
1286        size_t                       n         ,
1287        CppAD::vector<size_t>&       dep_taddr ,
1288        player<Base>*                play      ,
1289        recorder<Base>*              rec       ) 
1290{
1291        // temporary indices
1292        size_t i, j, k;
1293
1294        // temporary variables
1295        OpCode        op;   // current operator
1296        const addr_t* arg;  // operator arguments
1297        size_t        i_var;  // index of first result for current operator
1298
1299        // range and domain dimensions for F
1300        size_t m = dep_taddr.size();
1301
1302        // number of variables in the player
1303        const size_t num_var = play->num_var_rec(); 
1304
1305# ifndef NDEBUG
1306        // number of parameters in the player
1307        const size_t num_par = play->num_par_rec();
1308# endif
1309
1310        // number of  VecAD indices
1311        size_t num_vecad_ind   = play->num_vec_ind_rec();
1312
1313        // number of VecAD vectors
1314        size_t num_vecad_vec   = play->num_vecad_vec_rec();
1315
1316        // -------------------------------------------------------------
1317        // data structure that maps variable index in original operation
1318        // sequence to corresponding operator information
1319        CppAD::vector<struct struct_old_variable> tape(num_var);
1320        // -------------------------------------------------------------
1321        // Determine how each variable is connected to the dependent variables
1322
1323        // initialize all variables has having no connections
1324        for(i = 0; i < num_var; i++)
1325                tape[i].connect_type = not_connected;
1326
1327        for(j = 0; j < m; j++)
1328        {       // mark dependent variables as having one or more connections
1329                tape[ dep_taddr[j] ].connect_type = yes_connected;
1330        }
1331
1332        // vecad_connect contains a value for each VecAD object.
1333        // vecad maps a VecAD index (which corresponds to the beginning of the
1334        // VecAD object) to the vecad_connect falg for the VecAD object.
1335        CppAD::vector<enum_connect_type>   vecad_connect(num_vecad_vec);
1336        CppAD::vector<size_t> vecad(num_vecad_ind);
1337        j = 0;
1338        for(i = 0; i < num_vecad_vec; i++)
1339        {       vecad_connect[i] = not_connected;
1340                // length of this VecAD
1341                size_t length = play->GetVecInd(j);
1342                // set to proper index for this VecAD
1343                vecad[j] = i; 
1344                for(k = 1; k <= length; k++)
1345                        vecad[j+k] = num_vecad_vec; // invalid index
1346                // start of next VecAD
1347                j       += length + 1;
1348        }
1349        CPPAD_ASSERT_UNKNOWN( j == num_vecad_ind );
1350
1351        // work space used by UserOp.
1352        typedef std::set<size_t> size_set;
1353        vector<size_set> user_r_set;   // set sparsity pattern for result
1354        vector<size_set> user_s_set;   // set sparisty pattern for argument
1355        vector<bool>     user_r_bool;  // bool sparsity pattern for result
1356        vector<bool>     user_s_bool;  // bool sparisty pattern for argument
1357        //
1358        size_t user_q     = 0;       // column dimension for sparsity patterns
1359        size_t user_index = 0;       // indentifier for this user_atomic operation
1360        size_t user_id    = 0;       // user identifier for this call to operator
1361        size_t user_i     = 0;       // index in result vector
1362        size_t user_j     = 0;       // index in argument vector
1363        size_t user_m     = 0;       // size of result vector
1364        size_t user_n     = 0;       // size of arugment vector
1365        //
1366        atomic_base<Base>* user_atom = CPPAD_NULL; // current user atomic function
1367        bool               user_set  = true;       // use set sparsity (or bool)
1368
1369        // next expected operator in a UserOp sequence
1370        enum { user_start, user_arg, user_ret, user_end } user_state;
1371
1372        // During reverse mode, compute type of connection for each call to
1373        // a user atomic function.
1374        CppAD::vector<struct_user_info>    user_info;
1375        size_t                             user_curr = 0;
1376
1377        /// During reverse mode, information for each CSkip operation
1378        CppAD::vector<struct_cskip_info>   cskip_info;
1379
1380        // Initialize a reverse mode sweep through the operation sequence
1381        size_t i_op;
1382        play->reverse_start(op, arg, i_op, i_var);
1383        CPPAD_ASSERT_UNKNOWN( op == EndOp );
1384        size_t mask;
1385        user_state = user_end;
1386        while(op != BeginOp)
1387        {       // next op
1388                play->reverse_next(op, arg, i_op, i_var);
1389
1390                // Store the operator corresponding to each variable
1391                if( NumRes(op) > 0 )
1392                {       tape[i_var].op = op;
1393                        tape[i_var].arg = arg;
1394                }
1395# ifndef NDEBUG
1396                if( i_op <= n )
1397                {       CPPAD_ASSERT_UNKNOWN((op == InvOp) | (op == BeginOp));
1398                }
1399                else    CPPAD_ASSERT_UNKNOWN((op != InvOp) & (op != BeginOp));
1400# endif
1401                enum_connect_type connect_type      = tape[i_var].connect_type;
1402                std::set<class_cexp_pair>& cexp_set = tape[i_var].cexp_set;
1403                switch( op )
1404                {
1405                        // One variable corresponding to arg[0]
1406                        case AbsOp:
1407                        case AcosOp:
1408                        case AsinOp:
1409                        case AtanOp:
1410                        case CosOp:
1411                        case CoshOp:
1412                        case DivvpOp:
1413                        case ExpOp:
1414                        case LogOp:
1415                        case PowvpOp:
1416                        case SignOp:
1417                        case SinOp:
1418                        case SinhOp:
1419                        case SqrtOp:
1420                        case TanOp:
1421                        case TanhOp:
1422                        switch( connect_type )
1423                        {       case not_connected:
1424                                break;
1425       
1426                                case yes_connected:
1427                                case sum_connected:
1428                                case csum_connected: 
1429                                tape[arg[0]].connect_type = yes_connected;
1430                                break;
1431
1432                                case cexp_connected:
1433                                if( tape[arg[0]].connect_type == not_connected )
1434                                {       tape[arg[0]].connect_type = cexp_connected;
1435                                        tape[arg[0]].cexp_set     = cexp_set;
1436                                }
1437                                else if( tape[arg[0]].connect_type == cexp_connected )
1438                                {       tape[arg[0]].cexp_set = intersection(
1439                                                tape[arg[0]].cexp_set, cexp_set
1440                                        );
1441                                        if( tape[arg[0]].cexp_set.empty() )
1442                                                tape[arg[0]].connect_type = yes_connected;
1443                                }
1444                                else    tape[arg[0]].connect_type = yes_connected;
1445                                break;
1446
1447                                default:
1448                                CPPAD_ASSERT_UNKNOWN(false);
1449                        }
1450                        break; // --------------------------------------------
1451
1452                        // One variable corresponding to arg[1]
1453                        case DisOp:
1454                        case DivpvOp:
1455                        case MulpvOp:
1456                        case PowpvOp:
1457                        switch( connect_type )
1458                        {       case not_connected:
1459                                break;
1460       
1461                                case yes_connected:
1462                                case sum_connected:
1463                                case csum_connected: 
1464                                tape[arg[1]].connect_type = yes_connected;
1465                                break;
1466
1467                                case cexp_connected:
1468                                if( tape[arg[1]].connect_type == not_connected )
1469                                {       tape[arg[1]].connect_type = cexp_connected;
1470                                        tape[arg[1]].cexp_set     = cexp_set;
1471                                }
1472                                else if( tape[arg[1]].connect_type == cexp_connected )
1473                                {       tape[arg[1]].cexp_set = intersection(
1474                                                tape[arg[1]].cexp_set, cexp_set
1475                                        );
1476                                        if( tape[arg[1]].cexp_set.empty() )
1477                                                tape[arg[1]].connect_type = yes_connected;
1478                                }
1479                                else    tape[arg[1]].connect_type = yes_connected;
1480                                break;
1481
1482                                default:
1483                                CPPAD_ASSERT_UNKNOWN(false);
1484                        }
1485                        break; // --------------------------------------------
1486               
1487                        // Special case for SubvpOp
1488                        case SubvpOp:
1489                        switch( connect_type )
1490                        {       case not_connected:
1491                                break;
1492       
1493                                case yes_connected:
1494                                case sum_connected:
1495                                case csum_connected: 
1496                                if( tape[arg[0]].connect_type == not_connected )
1497                                        tape[arg[0]].connect_type = sum_connected;
1498                                else    tape[arg[0]].connect_type = yes_connected;
1499                                break;
1500
1501                                case cexp_connected:
1502                                if( tape[arg[0]].connect_type == not_connected )
1503                                {       tape[arg[0]].connect_type = cexp_connected;
1504                                        tape[arg[0]].cexp_set     = cexp_set;
1505                                }
1506                                else if( tape[arg[0]].connect_type == cexp_connected )
1507                                {       tape[arg[0]].cexp_set = intersection(
1508                                                tape[arg[0]].cexp_set, cexp_set
1509                                        );
1510                                        if( tape[arg[0]].cexp_set.empty() )
1511                                                tape[arg[0]].connect_type = yes_connected;
1512                                }
1513                                else    tape[arg[0]].connect_type = yes_connected;
1514                                break;
1515
1516                                default:
1517                                CPPAD_ASSERT_UNKNOWN(false);
1518                        }
1519                        if( connect_type == sum_connected )
1520                        {       // convert sum to csum connection for this variable
1521                                tape[i_var].connect_type = connect_type = csum_connected;
1522                        }
1523                        break; // --------------------------------------------
1524               
1525                        // Special case for AddpvOp and SubpvOp
1526                        case AddpvOp:
1527                        case SubpvOp:
1528                        switch( connect_type )
1529                        {       case not_connected:
1530                                break;
1531       
1532                                case yes_connected:
1533                                case sum_connected:
1534                                case csum_connected:
1535                                if( tape[arg[1]].connect_type == not_connected )
1536                                        tape[arg[1]].connect_type = sum_connected;
1537                                else    tape[arg[1]].connect_type = yes_connected;
1538                                break;
1539
1540                                case cexp_connected:
1541                                if( tape[arg[1]].connect_type == not_connected )
1542                                {       tape[arg[1]].connect_type = cexp_connected;
1543                                        tape[arg[1]].cexp_set     = cexp_set;
1544                                }
1545                                else if( tape[arg[1]].connect_type == cexp_connected )
1546                                {       tape[arg[1]].cexp_set = intersection(
1547                                                tape[arg[1]].cexp_set, cexp_set
1548                                        );
1549                                        if( tape[arg[1]].cexp_set.empty() )
1550                                                tape[arg[1]].connect_type = yes_connected;
1551                                }
1552                                else    tape[arg[1]].connect_type = yes_connected;
1553                                break;
1554
1555                                default:
1556                                CPPAD_ASSERT_UNKNOWN(false);
1557                        }
1558                        if( connect_type == sum_connected )
1559                        {       // convert sum to csum connection for this variable
1560                                tape[i_var].connect_type = connect_type = csum_connected;
1561                        }
1562                        break; // --------------------------------------------
1563
1564               
1565                        // Special case for AddvvOp and SubvvOp
1566                        case AddvvOp:
1567                        case SubvvOp:
1568                        for(i = 0; i < 2; i++) switch( connect_type )
1569                        {       case not_connected:
1570                                break;
1571
1572                                case yes_connected:
1573                                case sum_connected:
1574                                case csum_connected:
1575                                if( tape[arg[i]].connect_type == not_connected )
1576                                        tape[arg[i]].connect_type = sum_connected;
1577                                else    tape[arg[i]].connect_type = yes_connected;
1578                                break;
1579
1580                                case cexp_connected:
1581                                if( tape[arg[i]].connect_type == not_connected )
1582                                {       tape[arg[i]].connect_type = cexp_connected;
1583                                        tape[arg[i]].cexp_set     = cexp_set;
1584                                }
1585                                else if( tape[arg[i]].connect_type == cexp_connected )
1586                                {       tape[arg[i]].cexp_set = intersection(
1587                                                tape[arg[i]].cexp_set, cexp_set
1588                                        );
1589                                        if( tape[arg[i]].cexp_set.empty() )
1590                                                tape[arg[i]].connect_type = yes_connected;
1591                                }
1592                                else    tape[arg[i]].connect_type = yes_connected;
1593                                break;
1594
1595                                default:
1596                                CPPAD_ASSERT_UNKNOWN(false);
1597                        }
1598                        if( connect_type == sum_connected )
1599                        {       // convert sum to csum connection for this variable
1600                                tape[i_var].connect_type = connect_type = csum_connected;
1601                        }
1602                        break; // --------------------------------------------
1603
1604                        // Other binary operators
1605                        // where operands are arg[0], arg[1]
1606                        case DivvvOp:
1607                        case MulvvOp:
1608                        case PowvvOp:
1609                        for(i = 0; i < 2; i++) switch( connect_type )
1610                        {       case not_connected:
1611                                break;
1612
1613                                case yes_connected:
1614                                case sum_connected:
1615                                case csum_connected:
1616                                tape[arg[i]].connect_type = yes_connected;
1617                                break;
1618
1619                                case cexp_connected:
1620                                if( tape[arg[i]].connect_type == not_connected )
1621                                {       tape[arg[i]].connect_type = cexp_connected;
1622                                        tape[arg[i]].cexp_set     = cexp_set;
1623                                }
1624                                else if( tape[arg[i]].connect_type == cexp_connected )
1625                                {       tape[arg[i]].cexp_set = intersection(
1626                                                tape[arg[i]].cexp_set, cexp_set
1627                                        );
1628                                        if( tape[arg[i]].cexp_set.empty() )
1629                                                tape[arg[i]].connect_type = yes_connected;
1630                                }
1631                                else    tape[arg[i]].connect_type = yes_connected;
1632                                break;
1633
1634                                default:
1635                                CPPAD_ASSERT_UNKNOWN(false);
1636                        }
1637                        break; // --------------------------------------------
1638
1639                        // Conditional expression operators
1640                        case CExpOp:
1641                        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
1642                        if( connect_type != not_connected )
1643                        {       struct_cskip_info info;
1644                                info.cop        = CompareOp( arg[0] );
1645                                info.flag       = arg[1];
1646                                info.left       = arg[2];
1647                                info.right      = arg[3];
1648                                info.n_op_true  = 0;
1649                                info.n_op_false = 0;
1650                                info.i_arg      = 0; // case where no CSkipOp for this CExpOp
1651                                //
1652                                size_t index    = 0;
1653                                if( arg[1] & 1 )
1654                                {       index = std::max(index, info.left);
1655                                        tape[info.left].connect_type = yes_connected;
1656                                }
1657                                if( arg[1] & 2 )
1658                                {       index = std::max(index, info.right);
1659                                        tape[info.right].connect_type = yes_connected;
1660                                }
1661                                CPPAD_ASSERT_UNKNOWN( index > 0 );
1662                                info.max_left_right = index;
1663                                //
1664                                index = cskip_info.size();
1665                                cskip_info.push_back(info);
1666                                //
1667                                if( arg[1] & 4 )
1668                                {       if( tape[arg[4]].connect_type == not_connected )       
1669                                        {       tape[arg[4]].connect_type = cexp_connected;
1670                                                tape[arg[4]].cexp_set     = cexp_set;
1671                                                tape[arg[4]].cexp_set.insert(
1672                                                        class_cexp_pair(true, index)
1673                                                );
1674                                        }
1675                                        else
1676                                        {       // if arg[4] is cexp_connected, it could be
1677                                                // connected for both the true and false case
1678                                                // 2DO: if previously cexp_connected
1679                                                // and the true/false sense is the same, should
1680                                                // keep this conditional connnection.
1681                                                tape[arg[4]].cexp_set.clear();
1682                                                tape[arg[4]].connect_type = yes_connected;
1683                                        }
1684                                }
1685                                if( arg[1] & 8 )
1686                                {       if( tape[arg[5]].connect_type == not_connected )       
1687                                        {       tape[arg[5]].connect_type = cexp_connected;
1688                                                tape[arg[5]].cexp_set     = cexp_set;
1689                                                tape[arg[5]].cexp_set.insert(
1690                                                        class_cexp_pair(false, index)
1691                                                );
1692                                        }
1693                                        else
1694                                        {       tape[arg[5]].cexp_set.clear();
1695                                                tape[arg[5]].connect_type = yes_connected;
1696                                        }
1697                                }
1698                        }
1699                        break;  // --------------------------------------------
1700
1701                        // Operations where there is nothing to do
1702                        case EndOp:
1703                        case ParOp:
1704                        case PriOp:
1705                        break;  // --------------------------------------------
1706
1707                        // Operators that never get removed
1708                        case BeginOp:
1709                        case InvOp:
1710                        tape[i_var].connect_type = yes_connected;
1711                        break;
1712
1713                        // Compare operators never get removed -----------------
1714                        case LepvOp:
1715                        case LtpvOp:
1716                        case EqpvOp:
1717                        case NepvOp:
1718                        tape[arg[1]].connect_type = yes_connected;
1719                        break;
1720
1721                        case LevpOp:
1722                        case LtvpOp:
1723                        tape[arg[0]].connect_type = yes_connected;
1724                        break;
1725
1726                        case LevvOp:
1727                        case LtvvOp:
1728                        case EqvvOp:
1729                        case NevvOp:
1730                        tape[arg[0]].connect_type = yes_connected;
1731                        tape[arg[1]].connect_type = yes_connected;
1732                        break;
1733
1734                        // Load using a parameter index ----------------------
1735                        case LdpOp:
1736                        if( tape[i_var].connect_type != not_connected )
1737                        {
1738                                i                = vecad[ arg[0] - 1 ];
1739                                vecad_connect[i] = yes_connected;
1740                        }
1741                        break; // --------------------------------------------
1742
1743                        // Load using a variable index
1744                        case LdvOp:
1745                        if( tape[i_var].connect_type != not_connected )
1746                        {
1747                                i                    = vecad[ arg[0] - 1 ];
1748                                vecad_connect[i]     = yes_connected;
1749                                tape[arg[1]].connect_type = yes_connected;
1750                        }
1751                        break; // --------------------------------------------
1752
1753                        // Store a variable using a parameter index
1754                        case StpvOp:
1755                        i = vecad[ arg[0] - 1 ];
1756                        if( vecad_connect[i] != not_connected )
1757                                tape[arg[2]].connect_type = yes_connected;
1758                        break; // --------------------------------------------
1759
1760                        // Store a variable using a variable index
1761                        case StvvOp:
1762                        i = vecad[ arg[0] - 1 ];
1763                        if( vecad_connect[i] )
1764                        {       tape[arg[1]].connect_type = yes_connected;
1765                                tape[arg[2]].connect_type = yes_connected;
1766                        }
1767                        break; 
1768                        // ============================================================
1769                        case UserOp:
1770                        // start or end atomic operation sequence
1771                        CPPAD_ASSERT_UNKNOWN( NumRes( UserOp ) == 0 );
1772                        CPPAD_ASSERT_UNKNOWN( NumArg( UserOp ) == 4 );
1773                        if( user_state == user_end )
1774                        {       user_index = arg[0];
1775                                user_id    = arg[1];
1776                                user_n     = arg[2];
1777                                user_m     = arg[3];
1778                                user_q     = 1;
1779                                user_atom  = atomic_base<Base>::class_object(user_index);
1780                                user_set   = user_atom->sparsity() ==
1781                                        atomic_base<Base>::set_sparsity_enum;
1782                                //
1783                                if(user_s_set.size() != user_n )
1784                                        user_s_set.resize(user_n);
1785                                if(user_r_set.size() != user_m )
1786                                        user_r_set.resize(user_m);
1787                                //
1788                                // Note user_q is 1, but use it for clarity of code
1789                                if(user_s_bool.size() != user_n * user_q )
1790                                        user_s_bool.resize(user_n * user_q);
1791                                if(user_r_bool.size() != user_m * user_q )
1792                                        user_r_bool.resize(user_m * user_q);
1793                                //
1794                                user_j     = user_n;
1795                                user_i     = user_m;
1796                                user_state = user_ret;
1797                                //
1798                                struct_user_info info;
1799                                info.connect_type = not_connected;
1800                                info.op_end       = i_op + 1;
1801                                user_info.push_back(info);
1802                               
1803                        }
1804                        else
1805                        {       CPPAD_ASSERT_UNKNOWN( user_state == user_start );
1806                                CPPAD_ASSERT_UNKNOWN( user_index == size_t(arg[0]) );
1807                                CPPAD_ASSERT_UNKNOWN( user_id    == size_t(arg[1]) );
1808                                CPPAD_ASSERT_UNKNOWN( user_n     == size_t(arg[2]) );
1809                                CPPAD_ASSERT_UNKNOWN( user_m     == size_t(arg[3]) );
1810                                user_state = user_end;
1811                                //
1812                                CPPAD_ASSERT_UNKNOWN( user_curr + 1 == user_info.size() );
1813                                user_info[user_curr].op_begin = i_op;
1814                                user_curr                     = user_info.size();
1815               }
1816                        break;
1817
1818                        case UsrapOp:
1819                        // parameter argument in an atomic operation sequence
1820                        CPPAD_ASSERT_UNKNOWN( user_state == user_arg );
1821                        CPPAD_ASSERT_UNKNOWN( 0 < user_j && user_j <= user_n );
1822                        CPPAD_ASSERT_UNKNOWN( NumArg(op) == 1 );
1823                        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
1824                        --user_j;
1825                        if( user_j == 0 )
1826                                user_state = user_start;
1827                        break;
1828
1829                        case UsravOp:
1830                        // variable argument in an atomic operation sequence
1831                        CPPAD_ASSERT_UNKNOWN( user_state == user_arg );
1832                        CPPAD_ASSERT_UNKNOWN( 0 < user_j && user_j <= user_n );
1833                        CPPAD_ASSERT_UNKNOWN( NumArg(op) == 1 );
1834                        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) <= i_var );
1835                        CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
1836                        --user_j;
1837                        if( user_set )
1838                        {       if( ! user_s_set[user_j].empty() )
1839                                        tape[arg[0]].connect_type = 
1840                                                user_info[user_curr].connect_type;
1841                        }
1842                        else
1843                        {       if( user_s_bool[user_j] )
1844                                        tape[arg[0]].connect_type = 
1845                                                user_info[user_curr].connect_type;
1846                        }
1847                        if( user_j == 0 )
1848                                user_state = user_start;
1849                        break;
1850
1851                        case UsrrvOp:
1852                        // variable result in an atomic operation sequence
1853                        CPPAD_ASSERT_UNKNOWN( user_state == user_ret );
1854                        CPPAD_ASSERT_UNKNOWN( 0 < user_i && user_i <= user_m );
1855                        --user_i;
1856                        user_r_set[user_i].clear();
1857                        user_r_bool[user_i] = false;
1858                        switch( connect_type )
1859                        {       case not_connected:
1860                                break;
1861
1862                                case yes_connected:
1863                                case sum_connected:
1864                                case csum_connected:
1865                                user_info[user_curr].connect_type = yes_connected;
1866                                user_r_set[user_i].insert(0);
1867                                user_r_bool[user_i] = true;
1868                                break;
1869
1870                                case cexp_connected:
1871                                if( user_info[user_curr].connect_type == not_connected )
1872                                {       user_info[user_curr].connect_type  = connect_type;
1873                                        user_info[user_curr].cexp_set      = cexp_set;
1874                                }
1875                                else if(user_info[user_curr].connect_type==cexp_connected)
1876                                {       user_info[user_curr].cexp_set = intersection(
1877                                                user_info[user_curr].cexp_set, cexp_set
1878                                        );
1879                                        if( user_info[user_curr].cexp_set.empty() )
1880                                                user_info[user_curr].connect_type = yes_connected;
1881                                }
1882                                else    user_info[user_curr].connect_type = yes_connected;
1883                                user_r_set[user_i].insert(0);
1884                                user_r_bool[user_i] = true;
1885                                break;
1886
1887                                default:
1888                                CPPAD_ASSERT_UNKNOWN(false);
1889                        }
1890                        // drop into op = UsrrpOp code to handle case where user_i == 0
1891                        // for both UsrrvOp and UsrrpOp together.
1892
1893                        case UsrrpOp:
1894                        if( op == UsrrpOp )
1895                        {       // parameter result in an atomic operation sequence
1896                                CPPAD_ASSERT_UNKNOWN( user_state == user_ret );
1897                                CPPAD_ASSERT_UNKNOWN( 0 < user_i && user_i <= user_m );
1898                                CPPAD_ASSERT_UNKNOWN( NumArg(op) == 1 );
1899                                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
1900                                --user_i;
1901                                user_r_set[user_i].clear();
1902                                user_r_bool[user_i] = false;
1903                        }
1904                        if( user_i == 0 )
1905                        {       // call users function for this operation
1906                                user_atom->set_id(user_id);
1907# ifdef NDEBUG
1908                                if( user_set )
1909                                {       user_atom->
1910                                                rev_sparse_jac(user_q, user_r_set, user_s_set);
1911                                }
1912                                else
1913                                {       user_atom->
1914                                                rev_sparse_jac(user_q, user_r_bool, user_s_bool);
1915                                }
1916# else
1917                                bool flag;
1918                                if( user_set )
1919                                {       flag = user_atom->
1920                                                rev_sparse_jac(user_q, user_r_set, user_s_set);
1921                                }
1922                                else
1923                                {       flag = user_atom->
1924                                                rev_sparse_jac(user_q, user_r_bool, user_s_bool);
1925                                }
1926                                if( ! flag )
1927                                {       std::string s =
1928                                                "Optimizing an ADFun object"
1929                                                " that contains the atomic function\n\t";
1930                                        s += user_atom->afun_name();
1931                                        s += "\nCurrent atomic_sparsity is set to";
1932                                        //
1933                                        if( user_set )
1934                                                s += " std::set\nand std::set";
1935                                        else    s += " bool\nand bool";
1936                                        //
1937                                        s += " version of rev_sparse_jac returned false";
1938                                        CPPAD_ASSERT_KNOWN(false, s.c_str() );
1939                                }
1940# endif
1941                                user_state = user_arg;
1942                        }
1943                        break;
1944                        // ============================================================
1945
1946                        // all cases should be handled above
1947                        default:
1948                        CPPAD_ASSERT_UNKNOWN(0);
1949                }
1950        }
1951        // values corresponding to BeginOp
1952        CPPAD_ASSERT_UNKNOWN( i_op == 0 && i_var == 0 && op == BeginOp );
1953        tape[i_var].op = op;
1954        // -------------------------------------------------------------
1955
1956        // Determine which variables can be conditionally skipped
1957        for(i = 0; i < num_var; i++)
1958        {       if( tape[i].connect_type == cexp_connected )
1959                {       std::set<class_cexp_pair>::const_iterator itr = 
1960                                tape[i].cexp_set.begin();
1961                        while( itr != tape[i].cexp_set.end() )
1962                        {       j = itr->index_;
1963                                if( itr->compare_ == true )
1964                                        cskip_info[j].skip_var_false.push_back(i);
1965                                else cskip_info[j].skip_var_true.push_back(i);
1966                                itr++;
1967                        }
1968                }
1969        }
1970        // Determine size of skip information in user_info
1971        for(i = 0; i < user_info.size(); i++)
1972        {       if( user_info[i].connect_type == cexp_connected )
1973                {       std::set<class_cexp_pair>::const_iterator itr = 
1974                                user_info[i].cexp_set.begin();
1975                        while( itr != user_info[i].cexp_set.end() )
1976                        {       j = itr->index_;
1977                                if( itr->compare_ == true )
1978                                        cskip_info[j].n_op_false = 
1979                                                user_info[i].op_end - user_info[i].op_begin;
1980                                else
1981                                        cskip_info[j].n_op_true = 
1982                                                user_info[i].op_end - user_info[i].op_begin;
1983                                itr++;
1984                        }
1985                }
1986        }
1987
1988        // Sort the conditional skip information by the maximum of the
1989        // index for the left and right comparision operands
1990        CppAD::vector<size_t> cskip_info_order( cskip_info.size() );
1991        {       CppAD::vector<size_t> keys( cskip_info.size() );
1992                for(i = 0; i < cskip_info.size(); i++)
1993                        keys[i] = std::max( cskip_info[i].left, cskip_info[i].right );
1994                CppAD::index_sort(keys, cskip_info_order);
1995        }
1996        // index in sorted order
1997        size_t cskip_order_next = 0;
1998        // index in order during reverse sweep
1999        size_t cskip_info_index = cskip_info.size();
2000
2001
2002        // Initilaize table mapping hash code to variable index in tape
2003        // as pointing to the BeginOp at the beginning of the tape
2004        CppAD::vector<size_t>  hash_table_var(CPPAD_HASH_TABLE_SIZE);
2005        for(i = 0; i < CPPAD_HASH_TABLE_SIZE; i++)
2006                hash_table_var[i] = 0;
2007        CPPAD_ASSERT_UNKNOWN( tape[0].op == BeginOp );
2008
2009        // initialize mapping from old variable index to new
2010        // operator and variable index
2011        for(i = 0; i < num_var; i++)
2012        {       tape[i].new_op  = 0;       // invalid index (except for BeginOp)
2013                tape[i].new_var = num_var; // invalid index
2014        }
2015
2016        // Erase all information in the old recording
2017        rec->free();
2018
2019        // initialize mapping from old VecAD index to new VecAD index
2020        CppAD::vector<size_t> new_vecad_ind(num_vecad_ind);
2021        for(i = 0; i < num_vecad_ind; i++)
2022                new_vecad_ind[i] = num_vecad_ind; // invalid index
2023
2024        j = 0;     // index into the old set of indices
2025        for(i = 0; i < num_vecad_vec; i++)
2026        {       // length of this VecAD
2027                size_t length = play->GetVecInd(j);
2028                if( vecad_connect[i] != not_connected )
2029                {       // Put this VecAD vector in new recording
2030                        CPPAD_ASSERT_UNKNOWN(length < num_vecad_ind);
2031                        new_vecad_ind[j] = rec->PutVecInd(length);
2032                        for(k = 1; k <= length; k++) new_vecad_ind[j+k] =
2033                                rec->PutVecInd(
2034                                        rec->PutPar(
2035                                                play->GetPar( 
2036                                                        play->GetVecInd(j+k)
2037                        ) ) );
2038                }
2039                // start of next VecAD
2040                j       += length + 1;
2041        }
2042        CPPAD_ASSERT_UNKNOWN( j == num_vecad_ind );
2043
2044        // start playing the operations in the forward direction
2045        play->forward_start(op, arg, i_op, i_var);
2046        CPPAD_ASSERT_UNKNOWN( user_curr == user_info.size() );
2047
2048        // playing forward skips BeginOp at the beginning, but not EndOp at
2049        // the end.  Put BeginOp at beginning of recording
2050        CPPAD_ASSERT_UNKNOWN( op == BeginOp );
2051        CPPAD_ASSERT_NARG_NRES(BeginOp, 1, 1);
2052        tape[i_var].new_op  = rec->num_op_rec();
2053        tape[i_var].new_var = rec->PutOp(BeginOp);
2054        rec->PutArg(0);
2055
2056
2057        // temporary buffer for new argument values
2058        addr_t new_arg[6];
2059
2060        // temporary work space used by record_csum
2061        // (decalared here to avoid realloaction of memory)
2062        struct_csum_stacks csum_work;
2063
2064        // tempory used to hold a size_pair
2065        struct_size_pair size_pair;
2066
2067        user_state = user_start;
2068        while(op != EndOp)
2069        {       // next op
2070                play->forward_next(op, arg, i_op, i_var);
2071                CPPAD_ASSERT_UNKNOWN( (i_op > n)  | (op == InvOp) );
2072                CPPAD_ASSERT_UNKNOWN( (i_op <= n) | (op != InvOp) );
2073
2074                // determine if we should insert a conditional skip here
2075                bool skip = cskip_order_next < cskip_info.size();
2076                skip     &= op != BeginOp;
2077                skip     &= op != InvOp;
2078                skip     &= user_state == user_start;
2079                if( skip )
2080                {       j     = cskip_info_order[cskip_order_next];
2081                        if( NumRes(op) > 0 )
2082                                skip &= cskip_info[j].max_left_right < i_var;
2083                        else
2084                                skip &= cskip_info[j].max_left_right <= i_var;
2085                }
2086                if( skip )
2087                {       cskip_order_next++;
2088                        struct_cskip_info info = cskip_info[j];
2089                        size_t n_true  = info.skip_var_true.size() + info.n_op_true;
2090                        size_t n_false = info.skip_var_false.size() + info.n_op_false;
2091                        skip &= n_true > 0 || n_false > 0;
2092                        if( skip )
2093                        {       CPPAD_ASSERT_UNKNOWN( NumRes(CSkipOp) == 0 );
2094                                size_t n_arg   = 7 + n_true + n_false; 
2095                                // reserve space for the arguments to this operator but
2096                                // delay setting them until we have all the new addresses
2097                                cskip_info[j].i_arg = rec->ReserveArg(n_arg);
2098                                CPPAD_ASSERT_UNKNOWN( cskip_info[j].i_arg > 0 );
2099                                rec->PutOp(CSkipOp);
2100                        }
2101                }
2102
2103                // determine if we should keep this operation in the new
2104                // operation sequence
2105                bool keep;
2106                switch( op )
2107                {       // see wish_list/Optimize/CompareChange entry.
2108                        case EqpvOp:
2109                        case EqvvOp:
2110                        case LepvOp:
2111                        case LevpOp:
2112                        case LevvOp:
2113                        case LtpvOp:
2114                        case LtvpOp:
2115                        case LtvvOp:
2116                        case NepvOp:
2117                        case NevvOp:
2118                        keep = true;
2119                        break;
2120
2121                        case PriOp:
2122                        keep = false;
2123                        break;
2124
2125                        case InvOp:
2126                        case EndOp:
2127                        keep = true;
2128                        break;
2129
2130                        case StppOp:
2131                        case StvpOp:
2132                        case StpvOp:
2133                        case StvvOp:
2134                        CPPAD_ASSERT_UNKNOWN( NumRes(op) == 0 );
2135                        i = vecad[ arg[0] - 1 ];
2136                        keep = vecad_connect[i] != not_connected;
2137                        break;
2138
2139                        case AddpvOp:
2140                        case AddvvOp:
2141                        case SubpvOp:
2142                        case SubvpOp:
2143                        case SubvvOp:
2144                        keep  = tape[i_var].connect_type != not_connected;
2145                        keep &= tape[i_var].connect_type != csum_connected;
2146                        break; 
2147
2148                        case UserOp:
2149                        case UsrapOp:
2150                        case UsravOp:
2151                        case UsrrpOp:
2152                        case UsrrvOp:
2153                        keep = true;
2154                        break;
2155
2156                        default:
2157                        keep = tape[i_var].connect_type != not_connected;
2158                        break;
2159                }
2160
2161                unsigned short code         = 0;
2162                bool           replace_hash = false;
2163                addr_t         match_var;
2164                tape[i_var].match = false;
2165                if( keep ) switch( op )
2166                {
2167                        // Unary operator where operand is arg[0]
2168                        case AbsOp:
2169                        case AcosOp:
2170                        case AsinOp:
2171                        case AtanOp:
2172                        case CosOp:
2173                        case CoshOp:
2174                        case ExpOp:
2175                        case LogOp:
2176                        case SignOp:
2177                        case SinOp:
2178                        case SinhOp:
2179                        case SqrtOp:
2180                        case TanOp:
2181                        case TanhOp:
2182                        match_var = unary_match(
2183                                tape                ,  // inputs
2184                                i_var               ,
2185                                play->num_par_rec() ,
2186                                play->GetPar()      ,
2187                                hash_table_var      ,
2188                                code                  // outputs
2189                        );
2190                        if( match_var > 0 )
2191                        {       tape[i_var].match     = true;
2192                                tape[match_var].match = true;
2193                                tape[i_var].new_var   = tape[match_var].new_var;
2194                        }
2195                        else
2196                        {
2197                                replace_hash = true;
2198                                new_arg[0]   = tape[ arg[0] ].new_var;
2199                                rec->PutArg( new_arg[0] );
2200                                tape[i_var].new_op  = rec->num_op_rec();
2201                                tape[i_var].new_var = i = rec->PutOp(op);
2202                                CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < i );
2203                        }
2204                        break;
2205                        // ---------------------------------------------------
2206                        // Binary operators where
2207                        // left is a variable and right is a parameter
2208                        case SubvpOp:
2209                        if( tape[arg[0]].connect_type == csum_connected )
2210                        {
2211                                // convert to a sequence of summation operators
2212                                size_pair = record_csum(
2213                                        tape                , // inputs
2214                                        i_var               ,
2215                                        play->num_par_rec() ,
2216                                        play->GetPar()      ,
2217                                        rec                 ,
2218                                        csum_work
2219                                );
2220                                tape[i_var].new_op  = size_pair.i_op;
2221                                tape[i_var].new_var = size_pair.i_var;
2222                                // abort rest of this case
2223                                break;
2224                        }
2225                        case DivvpOp:
2226                        case PowvpOp:
2227                        match_var = binary_match(
2228                                tape                ,  // inputs
2229                                i_var               ,
2230                                play->num_par_rec() ,
2231                                play->GetPar()      ,
2232                                hash_table_var      ,
2233                                code                  // outputs
2234                        );
2235                        if( match_var > 0 )
2236                        {       tape[i_var].match     = true;
2237                                tape[match_var].match = true;
2238                                tape[i_var].new_var   = tape[match_var].new_var;
2239                        }
2240                        else
2241                        {       size_pair = record_vp(
2242                                        tape                , // inputs
2243                                        i_var               ,
2244                                        play->num_par_rec() ,
2245                                        play->GetPar()      ,
2246                                        rec                 ,
2247                                        op                  ,
2248                                        arg
2249                                );
2250                                tape[i_var].new_op  = size_pair.i_op;
2251                                tape[i_var].new_var = size_pair.i_var;
2252                                replace_hash = true;
2253                        }
2254                        break;
2255                        // ---------------------------------------------------
2256                        // Binary operators where
2257                        // left is an index and right is a variable
2258                        case DisOp:
2259                        match_var = binary_match(
2260                                tape                ,  // inputs
2261                                i_var               ,
2262                                play->num_par_rec() ,
2263                                play->GetPar()      ,
2264                                hash_table_var      ,
2265                                code                  // outputs
2266                        );
2267                        if( match_var > 0 )
2268                        {       tape[i_var].match     = true;
2269                                tape[match_var].match = true;
2270                                tape[i_var].new_var   = tape[match_var].new_var;
2271                        }
2272                        else
2273                        {       new_arg[0] = arg[0];
2274                                new_arg[1] = tape[ arg[1] ].new_var;
2275                                rec->PutArg( new_arg[0], new_arg[1] ); 
2276                                tape[i_var].new_op  = rec->num_op_rec();
2277                                tape[i_var].new_var = rec->PutOp(op);
2278                                CPPAD_ASSERT_UNKNOWN( 
2279                                        size_t(new_arg[1]) < tape[i_var].new_var
2280                                );
2281                                replace_hash = true;
2282                        }
2283                        break;
2284
2285                        // ---------------------------------------------------
2286                        // Binary operators where
2287                        // left is a parameter and right is a variable
2288                        case SubpvOp:
2289                        case AddpvOp:
2290                        if( tape[arg[1]].connect_type == csum_connected )
2291                        {
2292                                // convert to a sequence of summation operators
2293                                size_pair = record_csum(
2294                                        tape                , // inputs
2295                                        i_var               ,
2296                                        play->num_par_rec() ,
2297                                        play->GetPar()      ,
2298                                        rec                 ,
2299                                        csum_work
2300                                );
2301                                tape[i_var].new_op  = size_pair.i_op;
2302                                tape[i_var].new_var = size_pair.i_var;
2303                                // abort rest of this case
2304                                break;
2305                        }
2306                        case DivpvOp:
2307                        case MulpvOp:
2308                        case PowpvOp:
2309                        match_var = binary_match(
2310                                tape                ,  // inputs
2311                                i_var               ,
2312                                play->num_par_rec() ,
2313                                play->GetPar()      ,
2314                                hash_table_var      ,
2315                                code                  // outputs
2316                        );
2317                        if( match_var > 0 )
2318                        {       tape[i_var].match     = true;
2319                                tape[match_var].match = true;
2320                                tape[i_var].new_var   = tape[match_var].new_var;
2321                        }
2322                        else
2323                        {       size_pair = record_pv(
2324                                        tape                , // inputs
2325                                        i_var               ,
2326                                        play->num_par_rec() ,
2327                                        play->GetPar()      ,
2328                                        rec                 ,
2329                                        op                  ,
2330                                        arg
2331                                );
2332                                tape[i_var].new_op  = size_pair.i_op;
2333                                tape[i_var].new_var = size_pair.i_var;
2334                                replace_hash = true;
2335                        }
2336                        break;
2337                        // ---------------------------------------------------
2338                        // Binary operator where
2339                        // both operators are variables
2340                        case AddvvOp:
2341                        case SubvvOp:
2342                        if( (tape[arg[0]].connect_type == csum_connected) |
2343                            (tape[arg[1]].connect_type == csum_connected)
2344                        )
2345                        {
2346                                // convert to a sequence of summation operators
2347                                size_pair = record_csum(
2348                                        tape                , // inputs
2349                                        i_var               ,
2350                                        play->num_par_rec() ,
2351                                        play->GetPar()      ,
2352                                        rec                 ,
2353                                        csum_work
2354                                );
2355                                tape[i_var].new_op  = size_pair.i_op;
2356                                tape[i_var].new_var = size_pair.i_var;
2357                                // abort rest of this case
2358                                break;
2359                        }
2360                        case DivvvOp:
2361                        case MulvvOp:
2362                        case PowvvOp:
2363                        match_var = binary_match(
2364                                tape                ,  // inputs
2365                                i_var               ,
2366                                play->num_par_rec() ,
2367                                play->GetPar()      ,
2368                                hash_table_var      ,
2369                                code                  // outputs
2370                        );
2371                        if( match_var > 0 )
2372                        {       tape[i_var].match     = true;
2373                                tape[match_var].match = true;
2374                                tape[i_var].new_var   = tape[match_var].new_var;
2375                        }
2376                        else
2377                        {       size_pair = record_vv(
2378                                        tape                , // inputs
2379                                        i_var               ,
2380                                        play->num_par_rec() ,
2381                                        play->GetPar()      ,
2382                                        rec                 ,
2383                                        op                  ,
2384                                        arg
2385                                );
2386                                tape[i_var].new_op  = size_pair.i_op;
2387                                tape[i_var].new_var = size_pair.i_var;
2388                                replace_hash = true;
2389                        }
2390                        break;
2391                        // ---------------------------------------------------
2392                        // Conditional expression operators
2393                        case CExpOp:
2394                        CPPAD_ASSERT_NARG_NRES(op, 6, 1);
2395                        new_arg[0] = arg[0];
2396                        new_arg[1] = arg[1];
2397                        mask = 1;
2398                        for(i = 2; i < 6; i++)
2399                        {       if( arg[1] & mask )
2400                                {       new_arg[i] = tape[arg[i]].new_var;
2401                                        CPPAD_ASSERT_UNKNOWN( 
2402                                                size_t(new_arg[i]) < num_var
2403                                        );
2404                                }
2405                                else    new_arg[i] = rec->PutPar( 
2406                                                play->GetPar( arg[i] )
2407                                );
2408                                mask = mask << 1;
2409                        }
2410                        rec->PutArg(
2411                                new_arg[0] ,
2412                                new_arg[1] ,
2413                                new_arg[2] ,
2414                                new_arg[3] ,
2415                                new_arg[4] ,
2416                                new_arg[5] 
2417                        );
2418                        tape[i_var].new_op  = rec->num_op_rec();
2419                        tape[i_var].new_var = rec->PutOp(op);
2420                        //
2421                        // The new addresses for left and right are used during
2422                        // fill in the arguments for the CSkip operations. This does not
2423                        // affect max_left_right which is used during this sweep.
2424                        CPPAD_ASSERT_UNKNOWN( cskip_info_index > 0 );
2425                        cskip_info_index--;
2426                        cskip_info[ cskip_info_index ].left  = new_arg[2];
2427                        cskip_info[ cskip_info_index ].right = new_arg[3];
2428                        break;
2429                        // ---------------------------------------------------
2430                        // Operations with no arguments and no results
2431                        case EndOp:
2432                        CPPAD_ASSERT_NARG_NRES(op, 0, 0);
2433                        rec->PutOp(op);
2434                        break;
2435                        // ---------------------------------------------------
2436                        // Operations with two arguments and no results
2437                        case LepvOp:
2438                        case LtpvOp:
2439                        case EqpvOp:
2440                        case NepvOp:
2441                        CPPAD_ASSERT_NARG_NRES(op, 2, 0);
2442                        new_arg[0] = rec->PutPar( play->GetPar(arg[0]) );
2443                        new_arg[1] = tape[arg[1]].new_var;
2444                        rec->PutArg(new_arg[0], new_arg[1]);
2445                        rec->PutOp(op);
2446                        break;
2447                        //
2448                        case LevpOp:
2449                        case LtvpOp:
2450                        CPPAD_ASSERT_NARG_NRES(op, 2, 0);
2451                        new_arg[0] = tape[arg[0]].new_var;
2452                        new_arg[1] = rec->PutPar( play->GetPar(arg[1]) );
2453                        rec->PutArg(new_arg[0], new_arg[1]);
2454                        rec->PutOp(op);
2455                        break;
2456                        //
2457                        case LevvOp:
2458                        case LtvvOp:
2459                        case EqvvOp:
2460                        case NevvOp:
2461                        CPPAD_ASSERT_NARG_NRES(op, 2, 0);
2462                        new_arg[0] = tape[arg[0]].new_var;
2463                        new_arg[1] = tape[arg[1]].new_var;
2464                        rec->PutArg(new_arg[0], new_arg[1]);
2465                        rec->PutOp(op);
2466                        break;
2467
2468                        // ---------------------------------------------------
2469                        // Operations with no arguments and one result
2470                        case InvOp:
2471                        CPPAD_ASSERT_NARG_NRES(op, 0, 1);
2472                        tape[i_var].new_op  = rec->num_op_rec();
2473                        tape[i_var].new_var = rec->PutOp(op);
2474                        break;
2475                        // ---------------------------------------------------
2476                        // Operations with one argument that is a parameter
2477                        case ParOp:
2478                        CPPAD_ASSERT_NARG_NRES(op, 1, 1);
2479                        new_arg[0] = rec->PutPar( play->GetPar(arg[0] ) );
2480
2481                        rec->PutArg( new_arg[0] );
2482                        tape[i_var].new_op  = rec->num_op_rec();
2483                        tape[i_var].new_var = rec->PutOp(op);
2484                        break;
2485                        // ---------------------------------------------------
2486                        // Load using a parameter index
2487                        case LdpOp:
2488                        CPPAD_ASSERT_NARG_NRES(op, 3, 1);
2489                        new_arg[0] = new_vecad_ind[ arg[0] ];
2490                        new_arg[1] = arg[1];
2491                        new_arg[2] = rec->num_load_op_rec();
2492                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < num_vecad_ind );
2493                        rec->PutArg( 
2494                                new_arg[0], 
2495                                new_arg[1], 
2496                                new_arg[2]
2497                        );
2498                        tape[i_var].new_op  = rec->num_op_rec();
2499                        tape[i_var].new_var = rec->PutLoadOp(op);
2500                        break;
2501                        // ---------------------------------------------------
2502                        // Load using a variable index
2503                        case LdvOp:
2504                        CPPAD_ASSERT_NARG_NRES(op, 3, 1);
2505                        new_arg[0] = new_vecad_ind[ arg[0] ];
2506                        new_arg[1] = tape[arg[1]].new_var;
2507                        new_arg[2] = rec->num_load_op_rec();
2508                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < num_vecad_ind );
2509                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[1]) < num_var );
2510                        rec->PutArg( 
2511                                new_arg[0], 
2512                                new_arg[1], 
2513                                new_arg[2]
2514                        );
2515                        tape[i_var].new_var = rec->num_op_rec();
2516                        tape[i_var].new_var = rec->PutLoadOp(op);
2517                        break;
2518                        // ---------------------------------------------------
2519                        // Store a parameter using a parameter index
2520                        case StppOp:
2521                        CPPAD_ASSERT_NARG_NRES(op, 3, 0);
2522                        new_arg[0] = new_vecad_ind[ arg[0] ];
2523                        new_arg[1] = rec->PutPar( play->GetPar(arg[1]) );
2524                        new_arg[2] = rec->PutPar( play->GetPar(arg[2]) );
2525                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < num_vecad_ind );
2526                        rec->PutArg(
2527                                new_arg[0], 
2528                                new_arg[1], 
2529                                new_arg[2]
2530                        );
2531                        rec->PutOp(op);
2532                        break;
2533                        // ---------------------------------------------------
2534                        // Store a parameter using a variable index
2535                        case StvpOp:
2536                        CPPAD_ASSERT_NARG_NRES(op, 3, 0);
2537                        new_arg[0] = new_vecad_ind[ arg[0] ];
2538                        new_arg[1] = tape[arg[1]].new_var;
2539                        new_arg[2] = rec->PutPar( play->GetPar(arg[2]) );
2540                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < num_vecad_ind );
2541                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[1]) < num_var );
2542                        rec->PutArg(
2543                                new_arg[0], 
2544                                new_arg[1], 
2545                                new_arg[2]
2546                        );
2547                        rec->PutOp(op);
2548                        break;
2549                        // ---------------------------------------------------
2550                        // Store a variable using a parameter index
2551                        case StpvOp:
2552                        CPPAD_ASSERT_NARG_NRES(op, 3, 0);
2553                        new_arg[0] = new_vecad_ind[ arg[0] ];
2554                        new_arg[1] = rec->PutPar( play->GetPar(arg[1]) );
2555                        new_arg[2] = tape[arg[2]].new_var;
2556                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < num_vecad_ind );
2557                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[1]) < num_var );
2558                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[2]) < num_var );
2559                        rec->PutArg(
2560                                new_arg[0], 
2561                                new_arg[1], 
2562                                new_arg[2]
2563                        );
2564                        rec->PutOp(op);
2565                        break;
2566                        // ---------------------------------------------------
2567                        // Store a variable using a variable index
2568                        case StvvOp:
2569                        CPPAD_ASSERT_NARG_NRES(op, 3, 0);
2570                        new_arg[0] = new_vecad_ind[ arg[0] ];
2571                        new_arg[1] = tape[arg[1]].new_var;
2572                        new_arg[2] = tape[arg[2]].new_var;
2573                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[0]) < num_vecad_ind );
2574                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[1]) < num_var );
2575                        CPPAD_ASSERT_UNKNOWN( size_t(new_arg[2]) < num_var );
2576                        rec->PutArg(
2577                                new_arg[0], 
2578                                new_arg[1], 
2579                                new_arg[2]
2580                        );
2581                        rec->PutOp(op);
2582                        break;
2583
2584                        // -----------------------------------------------------------
2585                        case UserOp:
2586                        CPPAD_ASSERT_NARG_NRES(op, 4, 0);
2587                        if( user_state == user_start )
2588                        {       user_state = user_arg;
2589                                CPPAD_ASSERT_UNKNOWN( user_curr > 0 );
2590                                user_curr--;
2591                                user_info[user_curr].op_begin = rec->num_op_rec();
2592                        }
2593                        else
2594                        {       user_state = user_start;
2595                                user_info[user_curr].op_end = rec->num_op_rec() + 1;
2596                        }
2597                        // user_index, user_id, user_n, user_m
2598                        if( user_info[user_curr].connect_type != not_connected )
2599                        {       rec->PutArg(arg[0], arg[1], arg[2], arg[3]);
2600                                rec->PutOp(UserOp);
2601                        }
2602                        break;
2603
2604                        case UsrapOp:
2605                        CPPAD_ASSERT_NARG_NRES(op, 1, 0);
2606                        if( user_info[user_curr].connect_type != not_connected )
2607                        {       new_arg[0] = rec->PutPar( play->GetPar(arg[0]) );
2608                                rec->PutArg(new_arg[0]);
2609                                rec->PutOp(UsrapOp);
2610                        }
2611                        break;
2612
2613                        case UsravOp:
2614                        CPPAD_ASSERT_NARG_NRES(op, 1, 0);
2615                        if( user_info[user_curr].connect_type != not_connected )
2616                        {       new_arg[0] = tape[arg[0]].new_var;
2617                                if( size_t(new_arg[0]) < num_var )
2618                                {       rec->PutArg(new_arg[0]);
2619                                        rec->PutOp(UsravOp);
2620                                }
2621                                else
2622                                {       // This argument does not affect the result and
2623                                        // has been optimized out so use nan in its place.
2624                                        new_arg[0] = rec->PutPar( nan(Base(0)) );
2625                                        rec->PutArg(new_arg[0]);
2626                                        rec->PutOp(UsrapOp);
2627                                }
2628                        }
2629                        break;
2630
2631                        case UsrrpOp:
2632                        CPPAD_ASSERT_NARG_NRES(op, 1, 0);
2633                        if( user_info[user_curr].connect_type != not_connected )
2634                        {       new_arg[0] = rec->PutPar( play->GetPar(arg[0]) );
2635                                rec->PutArg(new_arg[0]);
2636                                rec->PutOp(UsrrpOp);
2637                        }
2638                        break;
2639                       
2640                        case UsrrvOp:
2641                        CPPAD_ASSERT_NARG_NRES(op, 0, 1);
2642                        if( user_info[user_curr].connect_type != not_connected )
2643                        {       tape[i_var].new_op  = rec->num_op_rec();
2644                                tape[i_var].new_var = rec->PutOp(UsrrvOp);
2645                        }
2646                        break;
2647                        // ---------------------------------------------------
2648
2649                        // all cases should be handled above
2650                        default:
2651                        CPPAD_ASSERT_UNKNOWN(false);
2652
2653                }
2654                if( replace_hash )
2655                {       // The old variable index i_var corresponds to the
2656                        // new variable index tape[i_var].new_var. In addition
2657                        // this is the most recent variable that has this code.
2658                        hash_table_var[code] = i_var;
2659                }
2660
2661        }
2662        // modify the dependent variable vector to new indices
2663        for(i = 0; i < dep_taddr.size(); i++ )
2664        {       CPPAD_ASSERT_UNKNOWN( size_t(tape[dep_taddr[i]].new_var) < num_var );
2665                dep_taddr[i] = tape[ dep_taddr[i] ].new_var;
2666        }
2667
2668# ifndef NDEBUG
2669        size_t num_new_op = rec->num_op_rec();
2670        for(i_var = 0; i_var < tape.size(); i_var++)
2671                CPPAD_ASSERT_UNKNOWN( tape[i_var].new_op < num_new_op );
2672# endif
2673
2674        // Move skip information from user_info to cskip_info
2675        for(i = 0; i < user_info.size(); i++)
2676        {       if( user_info[i].connect_type == cexp_connected )
2677                {       std::set<class_cexp_pair>::const_iterator itr =
2678                                user_info[i].cexp_set.begin();
2679                        while( itr != user_info[i].cexp_set.end() )
2680                        {       j = itr->index_;
2681                                k = user_info[i].op_begin;
2682                                while(k < user_info[i].op_end)
2683                                {       if( itr->compare_ == true )
2684                                                cskip_info[j].skip_op_false.push_back(k++);
2685                                        else    cskip_info[j].skip_op_true.push_back(k++);
2686                                }
2687                                itr++;
2688                        }
2689                }
2690        }
2691
2692        // fill in the arguments for the CSkip operations
2693        CPPAD_ASSERT_UNKNOWN( cskip_order_next == cskip_info.size() );
2694        for(i = 0; i < cskip_info.size(); i++)
2695        {       struct_cskip_info info = cskip_info[i];
2696                if( info.i_arg > 0 )
2697                {       CPPAD_ASSERT_UNKNOWN( info.n_op_true==info.skip_op_true.size() );
2698                        CPPAD_ASSERT_UNKNOWN(info.n_op_false==info.skip_op_false.size());
2699                        size_t n_true  = 
2700                                info.skip_var_true.size() + info.skip_op_true.size();
2701                        size_t n_false = 
2702                                info.skip_var_false.size() + info.skip_op_false.size();
2703                        size_t i_arg   = info.i_arg;
2704                        rec->ReplaceArg(i_arg++, info.cop   );
2705                        rec->ReplaceArg(i_arg++, info.flag  );
2706                        rec->ReplaceArg(i_arg++, info.left  );
2707                        rec->ReplaceArg(i_arg++, info.right );
2708                        rec->ReplaceArg(i_arg++, n_true     );
2709                        rec->ReplaceArg(i_arg++, n_false    );
2710                        for(j = 0; j < info.skip_var_true.size(); j++)
2711                        {       i_var = info.skip_var_true[j];
2712                                if( tape[i_var].match )
2713                                {       // the operation for this argument has been removed
2714                                        rec->ReplaceArg(i_arg++, 0);
2715                                }
2716                                else
2717                                {       CPPAD_ASSERT_UNKNOWN( tape[i_var].new_op > 0 );
2718                                        rec->ReplaceArg(i_arg++, tape[i_var].new_op );
2719                                }
2720                        } 
2721                        for(j = 0; j < info.skip_op_true.size(); j++)
2722                        {       i_op = info.skip_op_true[j];
2723                                rec->ReplaceArg(i_arg++, i_op);
2724                        } 
2725                        for(j = 0; j < info.skip_var_false.size(); j++)
2726                        {       i_var = info.skip_var_false[j];
2727                                if( tape[i_var].match )
2728                                {       // the operation for this argument has been removed
2729                                        rec->ReplaceArg(i_arg++, 0);
2730                                }
2731                                else
2732                                {       CPPAD_ASSERT_UNKNOWN( tape[i_var].new_op > 0 );
2733                                        rec->ReplaceArg(i_arg++, tape[i_var].new_op );
2734                                }
2735                        } 
2736                        for(j = 0; j < info.skip_op_false.size(); j++)
2737                        {       i_op = info.skip_op_false[j];
2738                                rec->ReplaceArg(i_arg++, i_op);
2739                        } 
2740                        rec->ReplaceArg(i_arg++, n_true + n_false);
2741# ifndef NDEBUG
2742                        size_t n_arg   = 7 + n_true + n_false; 
2743                        CPPAD_ASSERT_UNKNOWN( info.i_arg + n_arg == i_arg );
2744# endif
2745                }
2746        }
2747}
2748
2749} // END_CPPAD_OPTIMIZE_NAMESPACE
2750
2751/*!
2752Optimize a player object operation sequence
2753
2754The operation sequence for this object is replaced by one with fewer operations
2755but the same funcition and derivative values.
2756
2757\tparam Base
2758base type for the operator; i.e., this operation was recorded
2759using AD< \a Base > and computations by this routine are done using type
2760\a Base.
2761
2762*/
2763template <class Base>
2764void ADFun<Base>::optimize(void)
2765{       // place to store the optimized version of the recording
2766        recorder<Base> rec;
2767
2768        // number of independent variables
2769        size_t n = ind_taddr_.size();
2770
2771# ifndef NDEBUG
2772        size_t i, j, m = dep_taddr_.size();
2773        CppAD::vector<Base> x(n), y(m), check(m);
2774        Base max_taylor(0);
2775        bool check_zero_order = num_order_taylor_ > 0;
2776        if( check_zero_order )
2777        {       // zero order coefficients for independent vars
2778                for(j = 0; j < n; j++)
2779                {       CPPAD_ASSERT_UNKNOWN( play_.GetOp(j+1) == InvOp );
2780                        CPPAD_ASSERT_UNKNOWN( ind_taddr_[j]    == j+1   );
2781                        x[j] = taylor_[ ind_taddr_[j] * cap_order_taylor_ + 0];
2782                }
2783                // zero order coefficients for dependent vars
2784                for(i = 0; i < m; i++)
2785                {       CPPAD_ASSERT_UNKNOWN( dep_taddr_[i] < num_var_tape_  );
2786                        y[i] = taylor_[ dep_taddr_[i] * cap_order_taylor_ + 0];
2787                }
2788                // maximum zero order coefficient not counting BeginOp at beginning
2789                // (which is correpsonds to uninitialized memory).
2790                for(i = 1; i < num_var_tape_; i++)
2791                {       if(  abs_geq(taylor_[i*cap_order_taylor_+0] , max_taylor) )
2792                                max_taylor = taylor_[i*cap_order_taylor_+0];
2793                }
2794        }
2795# endif
2796
2797        // create the optimized recording
2798        CppAD::optimize::optimize_run<Base>(n, dep_taddr_, &play_, &rec);
2799
2800        // number of variables in the recording
2801        num_var_tape_  = rec.num_var_rec();
2802
2803        // now replace the recording
2804        play_.get(rec);
2805
2806        // set flag so this function knows it has been optimized
2807        has_been_optimized_ = true;
2808
2809        // free memory allocated for sparse Jacobian calculation
2810        // (the results are no longer valid)
2811        for_jac_sparse_pack_.resize(0, 0);
2812        for_jac_sparse_set_.resize(0,0);
2813
2814        // free old Taylor coefficient memory
2815        taylor_.free();
2816        num_order_taylor_     = 0;
2817        cap_order_taylor_     = 0;
2818
2819        // resize and initilaize conditional skip vector
2820        // (must use player size because it now has the recoreder information)
2821        cskip_op_.erase();
2822        cskip_op_.extend( play_.num_op_rec() );
2823
2824# ifndef NDEBUG
2825        if( check_zero_order )
2826        {
2827                // zero order forward calculation using new operation sequence
2828                check = Forward(0, x);
2829
2830                // check results
2831                Base eps = 10. * epsilon<Base>();
2832                for(i = 0; i < m; i++) CPPAD_ASSERT_KNOWN( 
2833                        abs_geq( eps * max_taylor , check[i] - y[i] ) ,
2834                        "Error during check of f.optimize()."
2835                );
2836
2837                // Erase memory that this calculation was done so NDEBUG gives
2838                // same final state for this object (from users perspective)
2839                num_order_taylor_     = 0;
2840        }
2841# endif
2842}
2843
2844} // END_CPPAD_NAMESPACE
2845# endif
Note: See TracBrowser for help on using the repository browser.