source: branches/cache/cppad/local/cond_op.hpp @ 3324

Last change on this file since 3324 was 3324, checked in by bradbell, 6 years ago

merge trunk changes into cache

  • Property svn:keywords set to Id
File size: 31.5 KB
Line 
1/* $Id: cond_op.hpp 3324 2014-09-12 12:14:53Z bradbell $ */
2# ifndef CPPAD_COND_OP_INCLUDED
3# define CPPAD_COND_OP_INCLUDED
4/* --------------------------------------------------------------------------
5CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 Bradley M. Bell
6
7CppAD is distributed under multiple licenses. This distribution is under
8the terms of the
9                    Eclipse Public License Version 1.0.
10
11A copy of this license is included in the COPYING file of this distribution.
12Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
13-------------------------------------------------------------------------- */
14
15namespace CppAD { // BEGIN_CPPAD_NAMESPACE
16/*!
17\file cond_op.hpp
18Forward, reverse, and sparse operations for conditional expressions.
19*/
20
21/*!
22Shared documentation for conditional expressions (not called).
23
24<!-- define conditional_exp_op -->
25The C++ source code coresponding to this operation is
26\verbatim
27        z = CondExpRel(y_0, y_1, y_2, y_3)
28\endverbatim
29where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
30
31\tparam Base
32base type for the operator; i.e., this operation was recorded
33using AD< \a Base > and computations by this routine are done using type
34\a Base.
35
36\param i_z
37is the AD variable index corresponding to the variable z.
38
39\param arg
40\n
41\a arg[0]
42is static cast to size_t from the enum type
43\verbatim
44        enum CompareOp {
45                CompareLt,
46                CompareLe,
47                CompareEq,
48                CompareGe,
49                CompareGt,
50                CompareNe
51        }
52\endverbatim
53for this operation.
54Note that arg[0] cannot be equal to CompareNe.
55\n
56\n
57\a arg[1] & 1
58\n
59If this is zero, y_0 is a parameter. Otherwise it is a variable.
60\n
61\n
62\a arg[1] & 2
63\n
64If this is zero, y_1 is a parameter. Otherwise it is a variable.
65\n
66\n
67\a arg[1] & 4
68\n
69If this is zero, y_2 is a parameter. Otherwise it is a variable.
70\n
71\n
72\a arg[1] & 8
73\n
74If this is zero, y_3 is a parameter. Otherwise it is a variable.
75\n
76\n
77\a arg[2 + j ] for j = 0, 1, 2, 3
78\n
79is the index corresponding to y_j.
80
81\param num_par
82is the total number of values in the vector \a parameter.
83
84\param parameter
85For j = 0, 1, 2, 3,
86if y_j is a parameter, \a parameter [ arg[2 + j] ] is its value.
87
88\param cap_order
89number of columns in the matrix containing the Taylor coefficients.
90
91\par Checked Assertions
92\li NumArg(CExpOp) == 6
93\li NumRes(CExpOp) == 1
94\li arg[0] < static_cast<size_t> ( CompareNe )
95\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
96\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
97<!-- end conditional_exp_op -->
98*/
99template <class Base>
100inline void conditional_exp_op(
101        size_t         i_z         ,
102        const addr_t*  arg         , 
103        size_t         num_par     ,
104        const Base*    parameter   ,
105        size_t         cap_order   )
106{       // This routine is only for documentation, it should never be used
107        CPPAD_ASSERT_UNKNOWN( false );
108}
109
110/*!
111Shared documentation for conditional expression sparse operations (not called).
112
113<!-- define sparse_conditional_exp_op -->
114The C++ source code coresponding to this operation is
115\verbatim
116        z = CondExpRel(y_0, y_1, y_2, y_3)
117\endverbatim
118where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
119
120\tparam Vector_set
121is the type used for vectors of sets. It can be either
122\c sparse_pack, \c sparse_set, or \c sparse_list.
123
124\param i_z
125is the AD variable index corresponding to the variable z.
126
127\param arg
128\n
129\a arg[0]
130is static cast to size_t from the enum type
131\verbatim
132        enum CompareOp {
133                CompareLt,
134                CompareLe,
135                CompareEq,
136                CompareGe,
137                CompareGt,
138                CompareNe
139        }
140\endverbatim
141for this operation.
142Note that arg[0] cannot be equal to CompareNe.
143\n
144\n
145\a arg[1] & 1
146\n
147If this is zero, y_0 is a parameter. Otherwise it is a variable.
148\n
149\n
150\a arg[1] & 2
151\n
152If this is zero, y_1 is a parameter. Otherwise it is a variable.
153\n
154\n
155\a arg[1] & 4
156\n
157If this is zero, y_2 is a parameter. Otherwise it is a variable.
158\n
159\n
160\a arg[1] & 8
161\n
162If this is zero, y_3 is a parameter. Otherwise it is a variable.
163\n
164\n
165\a arg[2 + j ] for j = 0, 1, 2, 3
166\n
167is the index corresponding to y_j.
168
169\param num_par
170is the total number of values in the vector \a parameter.
171
172\par Checked Assertions
173\li NumArg(CExpOp) == 6
174\li NumRes(CExpOp) == 1
175\li arg[0] < static_cast<size_t> ( CompareNe )
176\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
177\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
178<!-- end sparse_conditional_exp_op -->
179*/
180template <class Vector_set>
181inline void sparse_conditional_exp_op(
182        size_t         i_z           ,
183        const addr_t*  arg           , 
184        size_t         num_par       )
185{       // This routine is only for documentation, it should never be used
186        CPPAD_ASSERT_UNKNOWN( false );
187}
188
189/*!
190Compute forward mode Taylor coefficients for op = CExpOp.
191
192<!-- replace conditional_exp_op -->
193The C++ source code coresponding to this operation is
194\verbatim
195        z = CondExpRel(y_0, y_1, y_2, y_3)
196\endverbatim
197where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
198
199\tparam Base
200base type for the operator; i.e., this operation was recorded
201using AD< \a Base > and computations by this routine are done using type
202\a Base.
203
204\param i_z
205is the AD variable index corresponding to the variable z.
206
207\param arg
208\n
209\a arg[0]
210is static cast to size_t from the enum type
211\verbatim
212        enum CompareOp {
213                CompareLt,
214                CompareLe,
215                CompareEq,
216                CompareGe,
217                CompareGt,
218                CompareNe
219        }
220\endverbatim
221for this operation.
222Note that arg[0] cannot be equal to CompareNe.
223\n
224\n
225\a arg[1] & 1
226\n
227If this is zero, y_0 is a parameter. Otherwise it is a variable.
228\n
229\n
230\a arg[1] & 2
231\n
232If this is zero, y_1 is a parameter. Otherwise it is a variable.
233\n
234\n
235\a arg[1] & 4
236\n
237If this is zero, y_2 is a parameter. Otherwise it is a variable.
238\n
239\n
240\a arg[1] & 8
241\n
242If this is zero, y_3 is a parameter. Otherwise it is a variable.
243\n
244\n
245\a arg[2 + j ] for j = 0, 1, 2, 3
246\n
247is the index corresponding to y_j.
248
249\param num_par
250is the total number of values in the vector \a parameter.
251
252\param parameter
253For j = 0, 1, 2, 3,
254if y_j is a parameter, \a parameter [ arg[2 + j] ] is its value.
255
256\param cap_order
257number of columns in the matrix containing the Taylor coefficients.
258
259\par Checked Assertions
260\li NumArg(CExpOp) == 6
261\li NumRes(CExpOp) == 1
262\li arg[0] < static_cast<size_t> ( CompareNe )
263\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
264\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
265<!-- end conditional_exp_op -->
266
267\param p
268is the lowest order of the Taylor coefficient of z that we are computing.
269
270\param q
271is the highest order of the Taylor coefficient of z that we are computing.
272
273\param taylor
274\b Input:
275For j = 0, 1, 2, 3 and k = 0 , ... , q,
276if y_j is a variable then
277<code>taylor [ arg[2+j] * cap_order + k ]</code>
278is the k-th order Taylor coefficient corresponding to y_j.
279\n
280\b Input: <code>taylor [ i_z * cap_order + k ]</code>
281for k = 0 , ... , p-1,
282is the k-th order Taylor coefficient corresponding to z.
283\n
284\b Output: <code>taylor [ i_z * cap_order + k ]</code>
285for k = p , ... , q,
286is the k-th order Taylor coefficient corresponding to z.
287
288*/
289template <class Base>
290inline void forward_cond_op(
291        size_t         p           ,
292        size_t         q           ,
293        size_t         i_z         ,
294        const addr_t*  arg         , 
295        size_t         num_par     ,
296        const Base*    parameter   ,
297        size_t         cap_order   ,
298        Base*          taylor      )
299{       Base y_0, y_1, y_2, y_3;
300        Base zero(0);
301        Base* z = taylor + i_z * cap_order;
302
303        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
304        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
305        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
306        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
307
308        if( arg[1] & 1 )
309        {
310                y_0 = taylor[ arg[2] * cap_order + 0 ];
311        }
312        else
313        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
314                y_0 = parameter[ arg[2] ];
315        }
316        if( arg[1] & 2 )
317        {
318                y_1 = taylor[ arg[3] * cap_order + 0 ];
319        }
320        else
321        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
322                y_1 = parameter[ arg[3] ];
323        }
324        if( p == 0 )
325        {       if( arg[1] & 4 )
326                {
327                        y_2 = taylor[ arg[4] * cap_order + 0 ];
328                }
329                else
330                {       CPPAD_ASSERT_UNKNOWN( size_t(arg[4]) < num_par );
331                        y_2 = parameter[ arg[4] ];
332                }
333                if( arg[1] & 8 )
334                {
335                        y_3 = taylor[ arg[5] * cap_order + 0 ];
336                }
337                else
338                {       CPPAD_ASSERT_UNKNOWN( size_t(arg[5]) < num_par );
339                        y_3 = parameter[ arg[5] ];
340                }
341                z[0] = CondExpOp(
342                        CompareOp( arg[0] ),
343                        y_0,
344                        y_1,
345                        y_2,
346                        y_3
347                );
348                p++;
349        }
350        for(size_t d = p; d <= q; d++)
351        {       if( arg[1] & 4 )
352                {
353                        y_2 = taylor[ arg[4] * cap_order + d];
354                }
355                else    y_2 = zero;
356                if( arg[1] & 8 )
357                {
358                        y_3 = taylor[ arg[5] * cap_order + d];
359                }
360                else    y_3 = zero;
361                z[d] = CondExpOp(
362                        CompareOp( arg[0] ),
363                        y_0,
364                        y_1,
365                        y_2,
366                        y_3
367                );
368        }
369        return;
370}
371
372/*!
373Multiple directions forward mode Taylor coefficients for op = CExpOp.
374
375<!-- replace conditional_exp_op -->
376The C++ source code coresponding to this operation is
377\verbatim
378        z = CondExpRel(y_0, y_1, y_2, y_3)
379\endverbatim
380where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
381
382\tparam Base
383base type for the operator; i.e., this operation was recorded
384using AD< \a Base > and computations by this routine are done using type
385\a Base.
386
387\param i_z
388is the AD variable index corresponding to the variable z.
389
390\param arg
391\n
392\a arg[0]
393is static cast to size_t from the enum type
394\verbatim
395        enum CompareOp {
396                CompareLt,
397                CompareLe,
398                CompareEq,
399                CompareGe,
400                CompareGt,
401                CompareNe
402        }
403\endverbatim
404for this operation.
405Note that arg[0] cannot be equal to CompareNe.
406\n
407\n
408\a arg[1] & 1
409\n
410If this is zero, y_0 is a parameter. Otherwise it is a variable.
411\n
412\n
413\a arg[1] & 2
414\n
415If this is zero, y_1 is a parameter. Otherwise it is a variable.
416\n
417\n
418\a arg[1] & 4
419\n
420If this is zero, y_2 is a parameter. Otherwise it is a variable.
421\n
422\n
423\a arg[1] & 8
424\n
425If this is zero, y_3 is a parameter. Otherwise it is a variable.
426\n
427\n
428\a arg[2 + j ] for j = 0, 1, 2, 3
429\n
430is the index corresponding to y_j.
431
432\param num_par
433is the total number of values in the vector \a parameter.
434
435\param parameter
436For j = 0, 1, 2, 3,
437if y_j is a parameter, \a parameter [ arg[2 + j] ] is its value.
438
439\param cap_order
440number of columns in the matrix containing the Taylor coefficients.
441
442\par Checked Assertions
443\li NumArg(CExpOp) == 6
444\li NumRes(CExpOp) == 1
445\li arg[0] < static_cast<size_t> ( CompareNe )
446\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
447\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
448<!-- end conditional_exp_op -->
449
450\param q
451is order of the Taylor coefficient of z that we are computing.
452
453\param r
454is the number of Taylor coefficient directions that we are computing.
455
456\par tpv
457We use the notation
458<code>tpv = (cap_order-1) * r + 1</code>
459which is the number of Taylor coefficients per variable
460
461\param taylor
462\b Input:
463For j = 0, 1, 2, 3, k = 1, ..., q,
464if y_j is a variable then
465<code>taylor [ arg[2+j] * tpv + 0 ]</code>
466is the zero order Taylor coefficient corresponding to y_j and
467<code>taylor [ arg[2+j] * tpv + (k-1)*r+1+ell</code> is its
468k-th order Taylor coefficient in the ell-th direction.
469\n
470\b Input:
471For j = 0, 1, 2, 3, k = 1, ..., q-1,
472<code>taylor [ i_z * tpv + 0 ]</code> 
473is the zero order Taylor coefficient corresponding to z and
474<code>taylor [ i_z * tpv + (k-1)*r+1+ell</code> is its
475k-th order Taylor coefficient in the ell-th direction.
476\n
477\b Output: <code>taylor [ i_z * tpv + (q-1)*r+1+ell ]</code>
478is the q-th order Taylor coefficient corresponding to z
479in the ell-th direction.
480*/
481template <class Base>
482inline void forward_cond_op_dir(
483        size_t         q           ,
484        size_t         r           ,
485        size_t         i_z         ,
486        const addr_t*  arg         , 
487        size_t         num_par     ,
488        const Base*    parameter   ,
489        size_t         cap_order   ,
490        Base*          taylor      )
491{       Base y_0, y_1, y_2, y_3;
492        Base zero(0);
493        size_t num_taylor_per_var = (cap_order-1) * r + 1;
494        Base* z = taylor + i_z * num_taylor_per_var;
495
496        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
497        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
498        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
499        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
500        CPPAD_ASSERT_UNKNOWN( 0 < q );
501        CPPAD_ASSERT_UNKNOWN( q < cap_order );
502
503        if( arg[1] & 1 )
504        {
505                y_0 = taylor[ arg[2] * num_taylor_per_var + 0 ];
506        }
507        else
508        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
509                y_0 = parameter[ arg[2] ];
510        }
511        if( arg[1] & 2 )
512        {
513                y_1 = taylor[ arg[3] * num_taylor_per_var + 0 ];
514        }
515        else
516        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
517                y_1 = parameter[ arg[3] ];
518        }
519        size_t m = (q-1) * r + 1;
520        for(size_t ell = 0; ell < r; ell++)
521        {       if( arg[1] & 4 )
522                {
523                        y_2 = taylor[ arg[4] * num_taylor_per_var + m + ell];
524                }
525                else    y_2 = zero;
526                if( arg[1] & 8 )
527                {
528                        y_3 = taylor[ arg[5] * num_taylor_per_var + m + ell];
529                }
530                else    y_3 = zero;
531                z[m+ell] = CondExpOp(
532                        CompareOp( arg[0] ),
533                        y_0,
534                        y_1,
535                        y_2,
536                        y_3
537                );
538        }
539        return;
540}
541
542/*!
543Compute zero order forward mode Taylor coefficients for op = CExpOp.
544
545<!-- replace conditional_exp_op -->
546The C++ source code coresponding to this operation is
547\verbatim
548        z = CondExpRel(y_0, y_1, y_2, y_3)
549\endverbatim
550where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
551
552\tparam Base
553base type for the operator; i.e., this operation was recorded
554using AD< \a Base > and computations by this routine are done using type
555\a Base.
556
557\param i_z
558is the AD variable index corresponding to the variable z.
559
560\param arg
561\n
562\a arg[0]
563is static cast to size_t from the enum type
564\verbatim
565        enum CompareOp {
566                CompareLt,
567                CompareLe,
568                CompareEq,
569                CompareGe,
570                CompareGt,
571                CompareNe
572        }
573\endverbatim
574for this operation.
575Note that arg[0] cannot be equal to CompareNe.
576\n
577\n
578\a arg[1] & 1
579\n
580If this is zero, y_0 is a parameter. Otherwise it is a variable.
581\n
582\n
583\a arg[1] & 2
584\n
585If this is zero, y_1 is a parameter. Otherwise it is a variable.
586\n
587\n
588\a arg[1] & 4
589\n
590If this is zero, y_2 is a parameter. Otherwise it is a variable.
591\n
592\n
593\a arg[1] & 8
594\n
595If this is zero, y_3 is a parameter. Otherwise it is a variable.
596\n
597\n
598\a arg[2 + j ] for j = 0, 1, 2, 3
599\n
600is the index corresponding to y_j.
601
602\param num_par
603is the total number of values in the vector \a parameter.
604
605\param parameter
606For j = 0, 1, 2, 3,
607if y_j is a parameter, \a parameter [ arg[2 + j] ] is its value.
608
609\param cap_order
610number of columns in the matrix containing the Taylor coefficients.
611
612\par Checked Assertions
613\li NumArg(CExpOp) == 6
614\li NumRes(CExpOp) == 1
615\li arg[0] < static_cast<size_t> ( CompareNe )
616\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
617\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
618<!-- end conditional_exp_op -->
619
620\param taylor
621\b Input:
622For j = 0, 1, 2, 3,
623if y_j is a variable then
624\a taylor [ \a arg[2+j] * cap_order + 0 ]
625is the zero order Taylor coefficient corresponding to y_j.
626\n
627\b Output: \a taylor [ \a i_z * \a cap_order + 0 ]
628is the zero order Taylor coefficient corresponding to z.
629*/
630template <class Base>
631inline void forward_cond_op_0(
632        size_t         i_z         ,
633        const addr_t*  arg         , 
634        size_t         num_par     ,
635        const Base*    parameter   ,
636        size_t         cap_order   ,
637        Base*          taylor      )
638{       Base y_0, y_1, y_2, y_3;
639        Base* z;
640
641        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
642        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
643        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
644        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
645
646        if( arg[1] & 1 )
647        {
648                y_0 = taylor[ arg[2] * cap_order + 0 ];
649        }
650        else
651        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
652                y_0 = parameter[ arg[2] ];
653        }
654        if( arg[1] & 2 )
655        {
656                y_1 = taylor[ arg[3] * cap_order + 0 ];
657        }
658        else
659        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
660                y_1 = parameter[ arg[3] ];
661        }
662        if( arg[1] & 4 )
663        {
664                y_2 = taylor[ arg[4] * cap_order + 0 ];
665        }
666        else
667        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[4]) < num_par );
668                y_2 = parameter[ arg[4] ];
669        }
670        if( arg[1] & 8 )
671        {
672                y_3 = taylor[ arg[5] * cap_order + 0 ];
673        }
674        else
675        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[5]) < num_par );
676                y_3 = parameter[ arg[5] ];
677        }
678        z = taylor + i_z * cap_order;
679        z[0] = CondExpOp(
680                CompareOp( arg[0] ),
681                y_0,
682                y_1,
683                y_2,
684                y_3
685        );
686        return;
687}
688
689/*!
690Compute reverse mode Taylor coefficients for op = CExpOp.
691
692This routine is given the partial derivatives of a function
693G( z , y , x , w , ... )
694and it uses them to compute the partial derivatives of
695\verbatim
696        H( y , x , w , u , ... ) = G[ z(y) , y , x , w , u , ... ]
697\endverbatim
698where y above represents y_0, y_1, y_2, y_3.
699
700<!-- replace conditional_exp_op -->
701The C++ source code coresponding to this operation is
702\verbatim
703        z = CondExpRel(y_0, y_1, y_2, y_3)
704\endverbatim
705where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
706
707\tparam Base
708base type for the operator; i.e., this operation was recorded
709using AD< \a Base > and computations by this routine are done using type
710\a Base.
711
712\param i_z
713is the AD variable index corresponding to the variable z.
714
715\param arg
716\n
717\a arg[0]
718is static cast to size_t from the enum type
719\verbatim
720        enum CompareOp {
721                CompareLt,
722                CompareLe,
723                CompareEq,
724                CompareGe,
725                CompareGt,
726                CompareNe
727        }
728\endverbatim
729for this operation.
730Note that arg[0] cannot be equal to CompareNe.
731\n
732\n
733\a arg[1] & 1
734\n
735If this is zero, y_0 is a parameter. Otherwise it is a variable.
736\n
737\n
738\a arg[1] & 2
739\n
740If this is zero, y_1 is a parameter. Otherwise it is a variable.
741\n
742\n
743\a arg[1] & 4
744\n
745If this is zero, y_2 is a parameter. Otherwise it is a variable.
746\n
747\n
748\a arg[1] & 8
749\n
750If this is zero, y_3 is a parameter. Otherwise it is a variable.
751\n
752\n
753\a arg[2 + j ] for j = 0, 1, 2, 3
754\n
755is the index corresponding to y_j.
756
757\param num_par
758is the total number of values in the vector \a parameter.
759
760\param parameter
761For j = 0, 1, 2, 3,
762if y_j is a parameter, \a parameter [ arg[2 + j] ] is its value.
763
764\param cap_order
765number of columns in the matrix containing the Taylor coefficients.
766
767\par Checked Assertions
768\li NumArg(CExpOp) == 6
769\li NumRes(CExpOp) == 1
770\li arg[0] < static_cast<size_t> ( CompareNe )
771\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
772\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
773<!-- end conditional_exp_op -->
774
775\param d
776is the order of the Taylor coefficient of z that we are  computing.
777
778\param taylor
779\b Input:
780For j = 0, 1, 2, 3 and k = 0 , ... , \a d,
781if y_j is a variable then
782\a taylor [ \a arg[2+j] * cap_order + k ]
783is the k-th order Taylor coefficient corresponding to y_j.
784\n
785\a taylor [ \a i_z * \a cap_order + k ]
786for k = 0 , ... , \a d
787is the k-th order Taylor coefficient corresponding to z.
788
789\param nc_partial
790number of columns in the matrix containing the Taylor coefficients.
791
792\param partial
793\b Input:
794For j = 0, 1, 2, 3 and k = 0 , ... , \a d,
795if y_j is a variable then
796\a partial [ \a arg[2+j] * nc_partial + k ]
797is the partial derivative of G( z , y , x , w , u , ... )
798with respect to the k-th order Taylor coefficient corresponding to y_j.
799\n
800\b Input: \a partial [ \a i_z * \a cap_order + k ]
801for k = 0 , ... , \a d
802is the partial derivative of G( z , y , x , w , u , ... )
803with respect to the k-th order Taylor coefficient corresponding to z.
804\n
805\b Output:
806For j = 0, 1, 2, 3 and k = 0 , ... , \a d,
807if y_j is a variable then
808\a partial [ \a arg[2+j] * nc_partial + k ]
809is the partial derivative of H( y , x , w , u , ... )
810with respect to the k-th order Taylor coefficient corresponding to y_j.
811
812*/
813template <class Base>
814inline void reverse_cond_op(
815        size_t         d           ,
816        size_t         i_z         ,
817        const addr_t*  arg         , 
818        size_t         num_par     ,
819        const Base*    parameter   ,
820        size_t         cap_order   ,
821        const Base*    taylor      ,
822        size_t         nc_partial  ,
823        Base*          partial     )
824{       Base y_0, y_1;
825        Base zero(0);
826        Base* pz;
827        Base* py_2;
828        Base* py_3;
829
830        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
831        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
832        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
833        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
834
835        pz = partial + i_z * nc_partial + 0;
836        if( arg[1] & 1 )
837        {
838                y_0 = taylor[ arg[2] * cap_order + 0 ];
839        }
840        else
841        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
842                y_0 = parameter[ arg[2] ];
843        }
844        if( arg[1] & 2 )
845        {
846                y_1 = taylor[ arg[3] * cap_order + 0 ];
847        }
848        else
849        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
850                y_1 = parameter[ arg[3] ];
851        }
852        if( arg[1] & 4 )
853        {
854                py_2 = partial + arg[4] * nc_partial;
855                size_t j = d + 1;
856                while(j--)
857                {       py_2[j] += CondExpOp(
858                                CompareOp( arg[0] ),
859                                y_0,
860                                y_1,
861                                pz[j],
862                                zero
863                        );
864                }
865        }
866        if( arg[1] & 8 )
867        {
868                py_3 = partial + arg[5] * nc_partial;
869                size_t j = d + 1;
870                while(j--)
871                {       py_3[j] += CondExpOp(
872                                CompareOp( arg[0] ),
873                                y_0,
874                                y_1,
875                                zero,
876                                pz[j]
877                        );
878                }
879        }
880        return;
881}
882
883/*!
884Compute forward Jacobian sparsity patterns for op = CExpOp.
885
886<!-- replace sparse_conditional_exp_op -->
887The C++ source code coresponding to this operation is
888\verbatim
889        z = CondExpRel(y_0, y_1, y_2, y_3)
890\endverbatim
891where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
892
893\tparam Vector_set
894is the type used for vectors of sets. It can be either
895\c sparse_pack, \c sparse_set, or \c sparse_list.
896
897\param i_z
898is the AD variable index corresponding to the variable z.
899
900\param arg
901\n
902\a arg[0]
903is static cast to size_t from the enum type
904\verbatim
905        enum CompareOp {
906                CompareLt,
907                CompareLe,
908                CompareEq,
909                CompareGe,
910                CompareGt,
911                CompareNe
912        }
913\endverbatim
914for this operation.
915Note that arg[0] cannot be equal to CompareNe.
916\n
917\n
918\a arg[1] & 1
919\n
920If this is zero, y_0 is a parameter. Otherwise it is a variable.
921\n
922\n
923\a arg[1] & 2
924\n
925If this is zero, y_1 is a parameter. Otherwise it is a variable.
926\n
927\n
928\a arg[1] & 4
929\n
930If this is zero, y_2 is a parameter. Otherwise it is a variable.
931\n
932\n
933\a arg[1] & 8
934\n
935If this is zero, y_3 is a parameter. Otherwise it is a variable.
936\n
937\n
938\a arg[2 + j ] for j = 0, 1, 2, 3
939\n
940is the index corresponding to y_j.
941
942\param num_par
943is the total number of values in the vector \a parameter.
944
945\par Checked Assertions
946\li NumArg(CExpOp) == 6
947\li NumRes(CExpOp) == 1
948\li arg[0] < static_cast<size_t> ( CompareNe )
949\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
950\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
951<!-- end sparse_conditional_exp_op -->
952
953\param sparsity
954\b Input:
955if y_2 is a variable, the set with index t is
956the sparsity pattern corresponding to y_2.
957This identifies which of the independent variables the variable y_2
958depends on.
959\n
960\b Input:
961if y_3 is a variable, the set with index t is
962the sparsity pattern corresponding to y_3.
963This identifies which of the independent variables the variable y_3
964depends on.
965\n
966\b Output:
967The set with index T is
968the sparsity pattern corresponding to z.
969This identifies which of the independent variables the variable z
970depends on.
971*/
972template <class Vector_set>
973inline void forward_sparse_jacobian_cond_op(
974        size_t             i_z           ,
975        const addr_t*      arg           , 
976        size_t             num_par       ,
977        Vector_set&        sparsity      )
978{
979        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
980        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
981        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
982        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
983
984# ifndef NDEBUG
985        if( arg[1] & 1 )
986        {
987        }
988        else
989        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
990        }
991        if( arg[1] & 2 )
992        {
993        }
994        else
995        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
996        }
997# endif
998        if( arg[1] & 4 )
999        {
1000                if( arg[1] & 8 )
1001                {
1002                        sparsity.binary_union(i_z, arg[4], arg[5], sparsity);
1003                }
1004                else
1005                {       CPPAD_ASSERT_UNKNOWN( size_t(arg[5]) < num_par );
1006                        sparsity.assignment(i_z, arg[4], sparsity);
1007                }
1008        }       
1009        else
1010        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[4]) < num_par );
1011                if( arg[1] & 8 )
1012                {
1013                        sparsity.assignment(i_z, arg[5], sparsity);
1014                }
1015                else
1016                {       CPPAD_ASSERT_UNKNOWN( size_t(arg[5]) < num_par );
1017                        sparsity.clear(i_z);
1018                }
1019        }
1020        return;
1021}
1022
1023/*!
1024Compute reverse Jacobian sparsity patterns for op = CExpOp.
1025
1026This routine is given the sparsity patterns
1027for a function G(z, y, x, ... )
1028and it uses them to compute the sparsity patterns for
1029\verbatim
1030        H( y, x, w , u , ... ) = G[ z(x,y) , y , x , w , u , ... ]
1031\endverbatim
1032where y represents the combination of y_0, y_1, y_2, and y_3.
1033
1034<!-- replace sparse_conditional_exp_op -->
1035The C++ source code coresponding to this operation is
1036\verbatim
1037        z = CondExpRel(y_0, y_1, y_2, y_3)
1038\endverbatim
1039where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
1040
1041\tparam Vector_set
1042is the type used for vectors of sets. It can be either
1043\c sparse_pack, \c sparse_set, or \c sparse_list.
1044
1045\param i_z
1046is the AD variable index corresponding to the variable z.
1047
1048\param arg
1049\n
1050\a arg[0]
1051is static cast to size_t from the enum type
1052\verbatim
1053        enum CompareOp {
1054                CompareLt,
1055                CompareLe,
1056                CompareEq,
1057                CompareGe,
1058                CompareGt,
1059                CompareNe
1060        }
1061\endverbatim
1062for this operation.
1063Note that arg[0] cannot be equal to CompareNe.
1064\n
1065\n
1066\a arg[1] & 1
1067\n
1068If this is zero, y_0 is a parameter. Otherwise it is a variable.
1069\n
1070\n
1071\a arg[1] & 2
1072\n
1073If this is zero, y_1 is a parameter. Otherwise it is a variable.
1074\n
1075\n
1076\a arg[1] & 4
1077\n
1078If this is zero, y_2 is a parameter. Otherwise it is a variable.
1079\n
1080\n
1081\a arg[1] & 8
1082\n
1083If this is zero, y_3 is a parameter. Otherwise it is a variable.
1084\n
1085\n
1086\a arg[2 + j ] for j = 0, 1, 2, 3
1087\n
1088is the index corresponding to y_j.
1089
1090\param num_par
1091is the total number of values in the vector \a parameter.
1092
1093\par Checked Assertions
1094\li NumArg(CExpOp) == 6
1095\li NumRes(CExpOp) == 1
1096\li arg[0] < static_cast<size_t> ( CompareNe )
1097\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
1098\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
1099<!-- end sparse_conditional_exp_op -->
1100
1101\param nz_compare
1102Are the derivatives with respect to left and right of the expression below
1103considered to be non-zero:
1104\code
1105        CondExpRel(left, right, if_true, if_false)
1106\endcode
1107This is used by the optimizer to obtain the correct dependency relations.
1108
1109
1110\param sparsity
1111if y_2 is a variable, the set with index t is
1112the sparsity pattern corresponding to y_2.
1113This identifies which of the dependent variables depend on the variable y_2.
1114On input, this pattern corresponds to the function G.
1115On ouput, it corresponds to the function H.
1116\n
1117\n
1118if y_3 is a variable, the set with index t is
1119the sparsity pattern corresponding to y_3.
1120This identifies which of the dependent variables depeond on the variable y_3.
1121On input, this pattern corresponds to the function G.
1122On ouput, it corresponds to the function H.
1123\n
1124\b Output:
1125The set with index T is
1126the sparsity pattern corresponding to z.
1127This identifies which of the dependent variables depend on the variable z.
1128On input and output, this pattern corresponds to the function G.
1129*/
1130template <class Vector_set>
1131inline void reverse_sparse_jacobian_cond_op(
1132        bool                nz_compare    ,
1133        size_t              i_z           ,
1134        const addr_t*       arg           , 
1135        size_t              num_par       ,
1136        Vector_set&         sparsity      )
1137{       
1138        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
1139        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
1140        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
1141        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
1142
1143# ifndef NDEBUG
1144        if( arg[1] & 1 )
1145        {
1146        }
1147        else
1148        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
1149        }
1150        if( arg[1] & 2 )
1151        {
1152        }
1153        else
1154        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
1155        }
1156        if( ! ( arg[1] & 4 ) )
1157        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[4]) < num_par );
1158        }
1159        if( ! ( arg[1] & 8 ) )
1160        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[5]) < num_par );
1161        }
1162# endif
1163        if( nz_compare )
1164        {       if( arg[1] & 1 )
1165                {
1166                        sparsity.binary_union(arg[2], arg[2], i_z, sparsity);
1167                }
1168                if( arg[1] & 2 )
1169                {
1170                        sparsity.binary_union(arg[3], arg[3], i_z, sparsity);
1171                }
1172        }
1173        // --------------------------------------------------------------------
1174        if( arg[1] & 4 )
1175        {
1176                sparsity.binary_union(arg[4], arg[4], i_z, sparsity);
1177        }
1178        if( arg[1] & 8 )
1179        {
1180                sparsity.binary_union(arg[5], arg[5], i_z, sparsity);
1181        }
1182        return;
1183}
1184
1185/*!
1186Compute reverse Hessian sparsity patterns for op = CExpOp.
1187
1188This routine is given the sparsity patterns
1189for a function G(z, y, x, ... )
1190and it uses them to compute the sparsity patterns for
1191\verbatim
1192        H( y, x, w , u , ... ) = G[ z(x,y) , y , x , w , u , ... ]
1193\endverbatim
1194where y represents the combination of y_0, y_1, y_2, and y_3.
1195
1196<!-- replace sparse_conditional_exp_op -->
1197The C++ source code coresponding to this operation is
1198\verbatim
1199        z = CondExpRel(y_0, y_1, y_2, y_3)
1200\endverbatim
1201where Rel is one of the following: Lt, Le, Eq, Ge, Gt.
1202
1203\tparam Vector_set
1204is the type used for vectors of sets. It can be either
1205\c sparse_pack, \c sparse_set, or \c sparse_list.
1206
1207\param i_z
1208is the AD variable index corresponding to the variable z.
1209
1210\param arg
1211\n
1212\a arg[0]
1213is static cast to size_t from the enum type
1214\verbatim
1215        enum CompareOp {
1216                CompareLt,
1217                CompareLe,
1218                CompareEq,
1219                CompareGe,
1220                CompareGt,
1221                CompareNe
1222        }
1223\endverbatim
1224for this operation.
1225Note that arg[0] cannot be equal to CompareNe.
1226\n
1227\n
1228\a arg[1] & 1
1229\n
1230If this is zero, y_0 is a parameter. Otherwise it is a variable.
1231\n
1232\n
1233\a arg[1] & 2
1234\n
1235If this is zero, y_1 is a parameter. Otherwise it is a variable.
1236\n
1237\n
1238\a arg[1] & 4
1239\n
1240If this is zero, y_2 is a parameter. Otherwise it is a variable.
1241\n
1242\n
1243\a arg[1] & 8
1244\n
1245If this is zero, y_3 is a parameter. Otherwise it is a variable.
1246\n
1247\n
1248\a arg[2 + j ] for j = 0, 1, 2, 3
1249\n
1250is the index corresponding to y_j.
1251
1252\param num_par
1253is the total number of values in the vector \a parameter.
1254
1255\par Checked Assertions
1256\li NumArg(CExpOp) == 6
1257\li NumRes(CExpOp) == 1
1258\li arg[0] < static_cast<size_t> ( CompareNe )
1259\li arg[1] != 0; i.e., not all of y_0, y_1, y_2, y_3 are parameters.
1260\li For j = 0, 1, 2, 3 if y_j is a parameter, arg[2+j] < num_par.
1261<!-- end sparse_conditional_exp_op -->
1262
1263
1264\param jac_reverse
1265\a jac_reverse[i_z]
1266is false (true) if the Jacobian of G with respect to z is always zero
1267(may be non-zero).
1268\n
1269\n
1270\a jac_reverse[ arg[4] ]
1271If y_2 is a variable,
1272\a jac_reverse[ arg[4] ]
1273is false (true) if the Jacobian with respect to y_2 is always zero
1274(may be non-zero).
1275On input, it corresponds to the function G,
1276and on output it corresponds to the function H.
1277\n
1278\n
1279\a jac_reverse[ arg[5] ]
1280If y_3 is a variable,
1281\a jac_reverse[ arg[5] ]
1282is false (true) if the Jacobian with respect to y_3 is always zero
1283(may be non-zero).
1284On input, it corresponds to the function G,
1285and on output it corresponds to the function H.
1286
1287\param hes_sparsity
1288The set with index \a i_z in \a hes_sparsity
1289is the Hessian sparsity pattern for the function G
1290where one of the partials is with respect to z.
1291\n
1292\n
1293If y_2 is a variable,
1294the set with index \a arg[4] in \a hes_sparsity
1295is the Hessian sparsity pattern
1296where one of the partials is with respect to y_2.
1297On input, this pattern corresponds to the function G.
1298On output, this pattern corresponds to the function H.
1299\n
1300\n
1301If y_3 is a variable,
1302the set with index \a arg[5] in \a hes_sparsity
1303is the Hessian sparsity pattern
1304where one of the partials is with respect to y_3.
1305On input, this pattern corresponds to the function G.
1306On output, this pattern corresponds to the function H.
1307*/
1308template <class Vector_set>
1309inline void reverse_sparse_hessian_cond_op(
1310        size_t               i_z           ,
1311        const addr_t*        arg           , 
1312        size_t               num_par       ,
1313        bool*                jac_reverse   ,
1314        Vector_set&          hes_sparsity  )
1315{       
1316
1317        CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < static_cast<size_t> (CompareNe) );
1318        CPPAD_ASSERT_UNKNOWN( NumArg(CExpOp) == 6 );
1319        CPPAD_ASSERT_UNKNOWN( NumRes(CExpOp) == 1 );
1320        CPPAD_ASSERT_UNKNOWN( arg[1] != 0 );
1321
1322# ifndef NDEBUG
1323        if( arg[1] & 1 )
1324        {
1325        }
1326        else
1327        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
1328        }
1329        if( arg[1] & 2 )
1330        {
1331        }
1332        else
1333        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[3]) < num_par );
1334        }
1335        if( ! ( arg[1] & 4 ) )
1336        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[4]) < num_par );
1337        }
1338        if( ! ( arg[1] & 8 ) )
1339        {       CPPAD_ASSERT_UNKNOWN( size_t(arg[5]) < num_par );
1340        }
1341# endif
1342        if( arg[1] & 4 )
1343        {
1344
1345                hes_sparsity.binary_union(arg[4], arg[4], i_z, hes_sparsity);
1346                jac_reverse[ arg[4] ] |= jac_reverse[i_z];
1347        }
1348        if( arg[1] & 8 )
1349        {
1350
1351                hes_sparsity.binary_union(arg[5], arg[5], i_z, hes_sparsity);
1352                jac_reverse[ arg[5] ] |= jac_reverse[i_z];
1353        }
1354        return;
1355}
1356
1357} // END_CPPAD_NAMESPACE
1358# endif
Note: See TracBrowser for help on using the repository browser.