CppAD: A C++ Algorithmic Differentiation Package  20130102
forward0sweep.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_FORWARD0SWEEP_INCLUDED
00003 # define CPPAD_FORWARD0SWEEP_INCLUDED
00004 
00005 /* --------------------------------------------------------------------------
00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-12 Bradley M. Bell
00007 
00008 CppAD is distributed under multiple licenses. This distribution is under
00009 the terms of the 
00010                     Eclipse Public License Version 1.0.
00011 
00012 A copy of this license is included in the COPYING file of this distribution.
00013 Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
00014 -------------------------------------------------------------------------- */
00015 
00016 CPPAD_BEGIN_NAMESPACE
00017 /*!
00018 \defgroup forward0sweep_hpp forward0sweep.hpp
00019 \{
00020 \file forward0sweep.hpp
00021 Compute zero order forward mode Taylor coefficients.
00022 */
00023 
00024 /*!
00025 \def CPPAD_FORWARD0SWEEP_TRACE
00026 This value is either zero or one. 
00027 Zero is the normal operational value.
00028 If it is one, a trace of every forward0sweep computation is printed.
00029 (Note that forward0sweep is not used if CPPAD_USE_FORWARD0SWEEP is zero).
00030 */
00031 # define CPPAD_FORWARD0SWEEP_TRACE 0
00032 
00033 /*!
00034 Compute zero order forward mode Taylor coefficients.
00035 
00036 \tparam Base
00037 The type used during the forward mode computations; i.e., the corresponding
00038 recording of operations used the type \c AD<Base>.
00039 
00040 \param s_out
00041 Is the stream where output corresponding to \c PriOp operations will
00042 be written.
00043 
00044 \param print
00045 If \a print is false,
00046 suppress the output that is otherwise generated by the c PriOp instructions.
00047 
00048 \param n
00049 is the number of independent variables on the tape.
00050 
00051 \param numvar
00052 is the total number of variables on the tape.
00053 This is also equal to the number of rows in the matrix \a Taylor; i.e.,
00054 \a Rec->num_rec_var().
00055 
00056 \param Rec
00057 The information stored in \a Rec
00058 is a recording of the operations corresponding to the function
00059 \f[
00060      F : {\bf R}^n \rightarrow {\bf R}^m
00061 \f]
00062 where \f$ n \f$ is the number of independent variables and
00063 \f$ m \f$ is the number of dependent variables.
00064 \n
00065 \n
00066 The object \a Rec is effectly constant.
00067 There are two exceptions to this.
00068 The first exception is that while palying back the tape
00069 the object \a Rec holds information about the current location
00070 with in the tape and this changes during palyback. 
00071 The second exception is the fact that the 
00072 zero order ( \a d = 0 ) versions of the VecAD operators LdpOp and LdvOp 
00073 modify the corresponding \a op_arg values returned by 
00074 \ref player::next_forward and \ref player::next_reverse; see the
00075 \link load_op.hpp LdpOp and LdvOp \endlink operations.
00076 
00077 \param J
00078 Is the number of columns in the coefficient matrix \a Taylor.
00079 This must be greater than or equal one.
00080 
00081 \param Taylor
00082 \b Input: For j = 1 , ... , \a n,
00083 \a Taylor [ j * J + 0 ]
00084 variable with index i on the tape 
00085 (independent variable with index (j-1) in the independent variable vector).
00086 \n
00087 \n
00088 \b Output: For i = \a n + 1, ... , \a numvar - 1,
00089 \a Taylor [ i * J + 0 ]
00090 is the zero order Taylor coefficient for the variable with 
00091 index i on the tape.
00092 
00093 \a return
00094 The return value is equal to the number of ComOp operations
00095 that have a different result from when the information in 
00096 \a Rec was recorded.
00097 (Note that if NDEBUG is true, there are no ComOp operations
00098 in Rec and hence this return value is always zero.)
00099 */
00100 
00101 template <class Base>
00102 size_t forward0sweep(
00103      std::ostream&         s_out,
00104      bool                  print,
00105      size_t                n,
00106      size_t                numvar,
00107      player<Base>         *Rec,
00108      size_t                J,
00109      Base                 *Taylor
00110 )
00111 {    CPPAD_ASSERT_UNKNOWN( J >= 1 );
00112 
00113      // op code for current instruction
00114      OpCode           op;
00115 
00116      // index for current instruction
00117      size_t         i_op;
00118 
00119      // next variables 
00120      size_t        i_var;
00121 
00122      // constant and non-constant version of the operation argument indices
00123      addr_t*         non_const_arg;
00124      const addr_t*   arg = 0;
00125 
00126      // initialize the comparision operator (ComOp) counter
00127      size_t compareCount = 0;
00128 
00129      // This is an order zero calculation, initialize vector indices
00130      pod_vector<size_t> VectorInd;  // address for each element
00131      pod_vector<bool>   VectorVar;  // is element a variable
00132      size_t  i = Rec->num_rec_vecad_ind();
00133      if( i > 0 )
00134      {    VectorInd.extend(i);
00135           VectorVar.extend(i);
00136           while(i--)
00137           {    VectorInd[i] = Rec->GetVecInd(i);
00138                VectorVar[i] = false;
00139           }
00140      }
00141 
00142      // work space used by UserOp.
00143      const size_t user_k = 0;     // order of this forward mode calculation
00144      vector<Base> user_tx;        // argument vector Taylor coefficients 
00145      vector<Base> user_ty;        // result vector Taylor coefficients 
00146      size_t user_index = 0;       // indentifier for this user_atomic operation
00147      size_t user_id    = 0;       // user identifier for this call to operator
00148      size_t user_i     = 0;       // index in result vector
00149      size_t user_j     = 0;       // index in argument vector
00150      size_t user_m     = 0;       // size of result vector
00151      size_t user_n     = 0;       // size of arugment vector
00152      // next expected operator in a UserOp sequence
00153      enum { user_start, user_arg, user_ret, user_end } user_state = user_start;
00154 
00155      // check numvar argument
00156      CPPAD_ASSERT_UNKNOWN( Rec->num_rec_var() == numvar );
00157 
00158      // length of the parameter vector (used by CppAD assert macros)
00159      const size_t num_par = Rec->num_rec_par();
00160 
00161         // length of the text vector (used by CppAD assert macros)
00162         const size_t num_text = Rec->num_rec_text();
00163 
00164      // pointer to the beginning of the parameter vector
00165      const Base* parameter = CPPAD_NULL;
00166      if( num_par > 0 )
00167           parameter = Rec->GetPar();
00168 
00169      // pointer to the beginning of the text vector
00170      const char* text = CPPAD_NULL;
00171      if( num_text > 0 )
00172           text = Rec->GetTxt(0);
00173 
00174      // skip the BeginOp at the beginning of the recording
00175      Rec->start_forward(op, arg, i_op, i_var);
00176      CPPAD_ASSERT_UNKNOWN( op == BeginOp );
00177 # if CPPAD_FORWARD0SWEEP_TRACE
00178      std::cout << std::endl;
00179 # endif
00180      while(op != EndOp)
00181      {
00182           // this op
00183           Rec->next_forward(op, arg, i_op, i_var);
00184 # ifndef NDEBUG
00185           if( i_op <= n )
00186           {    CPPAD_ASSERT_UNKNOWN((op == InvOp) | (op == BeginOp));
00187           }
00188           else CPPAD_ASSERT_UNKNOWN((op != InvOp) & (op != BeginOp));
00189 # endif
00190 
00191           // action to take depends on the case
00192           switch( op )
00193           {
00194                case AbsOp:
00195                forward_abs_op_0(i_var, arg[0], J, Taylor);
00196                break;
00197                // -------------------------------------------------
00198 
00199                case AddvvOp:
00200                forward_addvv_op_0(i_var, arg, parameter, J, Taylor);
00201                break;
00202                // -------------------------------------------------
00203 
00204                case AddpvOp:
00205                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
00206                forward_addpv_op_0(i_var, arg, parameter, J, Taylor);
00207                break;
00208                // -------------------------------------------------
00209 
00210                case AcosOp:
00211                // sqrt(1 - x * x), acos(x)
00212                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00213                forward_acos_op_0(i_var, arg[0], J, Taylor);
00214                break;
00215                // -------------------------------------------------
00216 
00217                case AsinOp:
00218                // sqrt(1 - x * x), asin(x)
00219                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00220                forward_asin_op_0(i_var, arg[0], J, Taylor);
00221                break;
00222                // -------------------------------------------------
00223 
00224                case AtanOp:
00225                // 1 + x * x, atan(x)
00226                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00227                forward_atan_op_0(i_var, arg[0], J, Taylor);
00228                break;
00229                // -------------------------------------------------
00230 
00231                case CSumOp:
00232                // CSumOp has a variable number of arguments and
00233                // next_forward thinks it one has one argument.
00234                // we must inform next_forward of this special case.
00235                Rec->forward_csum(op, arg, i_op, i_var);
00236                forward_csum_op(
00237                     0, i_var, arg, num_par, parameter, J, Taylor
00238                );
00239                break;
00240 
00241                // -------------------------------------------------
00242                case CExpOp:
00243                // Use the general case with d == 0 
00244                // (could create an optimzied verison for this case)
00245                forward_cond_op_0(
00246                     i_var, arg, num_par, parameter, J, Taylor
00247                );
00248                break;
00249                // ---------------------------------------------------
00250                case ComOp:
00251                forward_comp_op_0(
00252                compareCount, arg, num_par, parameter, J, Taylor
00253                );
00254                break;
00255                // ---------------------------------------------------
00256 
00257                case CosOp:
00258                // sin(x), cos(x)
00259                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00260                forward_cos_op_0(i_var, arg[0], J, Taylor);
00261                break;
00262                // ---------------------------------------------------
00263 
00264                case CoshOp:
00265                // sinh(x), cosh(x)
00266                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00267                forward_cosh_op_0(i_var, arg[0], J, Taylor);
00268                break;
00269                // -------------------------------------------------
00270 
00271                case DisOp:
00272                forward_dis_op_0(i_var, arg, J, Taylor);
00273                break;
00274                // -------------------------------------------------
00275 
00276                case DivvvOp:
00277                forward_divvv_op_0(i_var, arg, parameter, J, Taylor);
00278                break;
00279                // -------------------------------------------------
00280 
00281                case DivpvOp:
00282                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
00283                forward_divpv_op_0(i_var, arg, parameter, J, Taylor);
00284                break;
00285                // -------------------------------------------------
00286 
00287                case DivvpOp:
00288                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < num_par );
00289                forward_divvp_op_0(i_var, arg, parameter, J, Taylor);
00290                break;
00291                // -------------------------------------------------
00292 
00293                case EndOp:
00294                CPPAD_ASSERT_NARG_NRES(op, 0, 0);
00295                break;
00296                // -------------------------------------------------
00297 
00298                case ExpOp:
00299                forward_exp_op_0(i_var, arg[0], J, Taylor);
00300                break;
00301                // -------------------------------------------------
00302 
00303                case InvOp:
00304                break;
00305                // -------------------------------------------------
00306 
00307                case LdpOp:
00308                CPPAD_ASSERT_UNKNOWN( VectorInd.size() != 0 );
00309                CPPAD_ASSERT_UNKNOWN( VectorVar.size() != 0 );
00310                non_const_arg = Rec->forward_non_const_arg();
00311                forward_load_p_op_0(
00312                     i_var, 
00313                     non_const_arg, 
00314                     num_par, 
00315                     parameter, 
00316                     J, 
00317                     Taylor,
00318                     Rec->num_rec_vecad_ind(),
00319                     VectorVar.data(),
00320                     VectorInd.data()
00321                );
00322                break;
00323                // -------------------------------------------------
00324 
00325                case LdvOp:
00326                CPPAD_ASSERT_UNKNOWN( VectorInd.size() != 0 );
00327                CPPAD_ASSERT_UNKNOWN( VectorVar.size() != 0 );
00328                non_const_arg = Rec->forward_non_const_arg();
00329                forward_load_v_op_0(
00330                     i_var, 
00331                     non_const_arg, 
00332                     num_par, 
00333                     parameter, 
00334                     J, 
00335                     Taylor,
00336                     Rec->num_rec_vecad_ind(),
00337                     VectorVar.data(),
00338                     VectorInd.data()
00339                );
00340                break;
00341                // -------------------------------------------------
00342 
00343                case LogOp:
00344                forward_log_op_0(i_var, arg[0], J, Taylor);
00345                break;
00346                // -------------------------------------------------
00347 
00348                case MulvvOp:
00349                forward_mulvv_op_0(i_var, arg, parameter, J, Taylor);
00350                break;
00351                // -------------------------------------------------
00352 
00353                case MulpvOp:
00354                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
00355                forward_mulpv_op_0(i_var, arg, parameter, J, Taylor);
00356                break;
00357                // -------------------------------------------------
00358 
00359                case ParOp:
00360                forward_par_op_0(
00361                     i_var, arg, num_par, parameter, J, Taylor
00362                );
00363                break;
00364                // -------------------------------------------------
00365 
00366                case PowvpOp:
00367                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < num_par );
00368                forward_powvp_op_0(i_var, arg, parameter, J, Taylor);
00369                break;
00370                // -------------------------------------------------
00371 
00372                case PowpvOp:
00373                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
00374                forward_powpv_op_0(i_var, arg, parameter, J, Taylor);
00375                break;
00376                // -------------------------------------------------
00377 
00378                case PowvvOp:
00379                forward_powvv_op_0(i_var, arg, parameter, J, Taylor);
00380                break;
00381                // -------------------------------------------------
00382 
00383                case PriOp:
00384                if( print ) forward_pri_0(s_out,
00385                     i_var, arg, num_text, text, num_par, parameter, J, Taylor
00386                );
00387                break;
00388                // -------------------------------------------------
00389 
00390                case SignOp:
00391                // cos(x), sin(x)
00392                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00393                forward_sign_op_0(i_var, arg[0], J, Taylor);
00394                break;
00395                // -------------------------------------------------
00396 
00397                case SinOp:
00398                // cos(x), sin(x)
00399                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00400                forward_sin_op_0(i_var, arg[0], J, Taylor);
00401                break;
00402                // -------------------------------------------------
00403 
00404                case SinhOp:
00405                // cosh(x), sinh(x)
00406                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00407                forward_sinh_op_0(i_var, arg[0], J, Taylor);
00408                break;
00409                // -------------------------------------------------
00410 
00411                case SqrtOp:
00412                forward_sqrt_op_0(i_var, arg[0], J, Taylor);
00413                break;
00414                // -------------------------------------------------
00415 
00416                case StppOp:
00417                forward_store_pp_op_0(
00418                     i_var, 
00419                     arg, 
00420                     num_par, 
00421                     J, 
00422                     Taylor,
00423                     Rec->num_rec_vecad_ind(),
00424                     VectorVar.data(),
00425                     VectorInd.data()
00426                );
00427                break;
00428                // -------------------------------------------------
00429 
00430                case StpvOp:
00431                forward_store_pv_op_0(
00432                     i_var, 
00433                     arg, 
00434                     num_par, 
00435                     J, 
00436                     Taylor,
00437                     Rec->num_rec_vecad_ind(),
00438                     VectorVar.data(),
00439                     VectorInd.data()
00440                );
00441                break;
00442                // -------------------------------------------------
00443 
00444                case StvpOp:
00445                forward_store_vp_op_0(
00446                     i_var, 
00447                     arg, 
00448                     num_par, 
00449                     J, 
00450                     Taylor,
00451                     Rec->num_rec_vecad_ind(),
00452                     VectorVar.data(),
00453                     VectorInd.data()
00454                );
00455                break;
00456                // -------------------------------------------------
00457 
00458                case StvvOp:
00459                forward_store_vv_op_0(
00460                     i_var, 
00461                     arg, 
00462                     num_par, 
00463                     J, 
00464                     Taylor,
00465                     Rec->num_rec_vecad_ind(),
00466                     VectorVar.data(),
00467                     VectorInd.data()
00468                );
00469                break;
00470                // -------------------------------------------------
00471 
00472                case SubvvOp:
00473                forward_subvv_op_0(i_var, arg, parameter, J, Taylor);
00474                break;
00475                // -------------------------------------------------
00476 
00477                case SubpvOp:
00478                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
00479                forward_subpv_op_0(i_var, arg, parameter, J, Taylor);
00480                break;
00481                // -------------------------------------------------
00482 
00483                case SubvpOp:
00484                CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < num_par );
00485                forward_subvp_op_0(i_var, arg, parameter, J, Taylor);
00486                break;
00487                // -------------------------------------------------
00488 
00489                case TanOp:
00490                // tan(x)^2, tan(x)
00491                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00492                forward_tan_op_0(i_var, arg[0], J, Taylor);
00493                break;
00494                // -------------------------------------------------
00495 
00496                case TanhOp:
00497                // tanh(x)^2, tanh(x)
00498                CPPAD_ASSERT_UNKNOWN( i_var < numvar  );
00499                forward_tanh_op_0(i_var, arg[0], J, Taylor);
00500                break;
00501                // -------------------------------------------------
00502 
00503                case UserOp:
00504                // start or end an atomic operation sequence
00505                CPPAD_ASSERT_UNKNOWN( NumRes( UserOp ) == 0 );
00506                CPPAD_ASSERT_UNKNOWN( NumArg( UserOp ) == 4 );
00507                if( user_state == user_start )
00508                {    user_index = arg[0];
00509                     user_id    = arg[1];
00510                     user_n     = arg[2];
00511                     user_m     = arg[3];
00512                     if(user_tx.size() < user_n)
00513                          user_tx.resize(user_n);
00514                     if(user_ty.size() < user_m)
00515                          user_ty.resize(user_m);
00516                     user_j     = 0;
00517                     user_i     = 0;
00518                     user_state = user_arg;
00519                }
00520                else
00521                {    CPPAD_ASSERT_UNKNOWN( user_state == user_end );
00522                     CPPAD_ASSERT_UNKNOWN( user_index == size_t(arg[0]) );
00523                     CPPAD_ASSERT_UNKNOWN( user_id    == size_t(arg[1]) );
00524                     CPPAD_ASSERT_UNKNOWN( user_n     == size_t(arg[2]) );
00525                     CPPAD_ASSERT_UNKNOWN( user_m     == size_t(arg[3]) );
00526                     user_state = user_start;
00527                }
00528                break;
00529 
00530                case UsrapOp:
00531                // parameter argument in an atomic operation sequence
00532                CPPAD_ASSERT_UNKNOWN( user_state == user_arg );
00533                CPPAD_ASSERT_UNKNOWN( user_j < user_n );
00534                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_par );
00535                user_tx[user_j++] = parameter[ arg[0] ];
00536                if( user_j == user_n )
00537                {    // call users function for this operation
00538                     user_atomic<Base>::forward(user_index, user_id, 
00539                          user_k, user_n, user_m, user_tx, user_ty
00540                     );
00541                     user_state = user_ret;
00542                }
00543                break;
00544 
00545                case UsravOp:
00546                // variable argument in an atomic operation sequence
00547                CPPAD_ASSERT_UNKNOWN( user_state == user_arg );
00548                CPPAD_ASSERT_UNKNOWN( user_j < user_n );
00549                CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) <= i_var );
00550                user_tx[user_j++] = Taylor[ arg[0] * J + 0 ];
00551                if( user_j == user_n )
00552                {    // call users function for this operation
00553                     user_atomic<Base>::forward(user_index, user_id,
00554                          user_k, user_n, user_m, user_tx, user_ty
00555                     );
00556                     user_state = user_ret;
00557                }
00558                break;
00559 
00560                case UsrrpOp:
00561                // parameter result in an atomic operation sequence
00562                CPPAD_ASSERT_UNKNOWN( user_state == user_ret );
00563                CPPAD_ASSERT_UNKNOWN( user_i < user_m );
00564                user_i++;
00565                if( user_i == user_m )
00566                     user_state = user_end;
00567                break;
00568 
00569                case UsrrvOp:
00570                // variable result in an atomic operation sequence
00571                CPPAD_ASSERT_UNKNOWN( user_state == user_ret );
00572                CPPAD_ASSERT_UNKNOWN( user_i < user_m );
00573                Taylor[ i_var * J + 0 ] = user_ty[user_i++];
00574                if( user_i == user_m )
00575                     user_state = user_end;
00576                break;
00577                // -------------------------------------------------
00578 
00579                default:
00580                CPPAD_ASSERT_UNKNOWN(0);
00581           }
00582 # if CPPAD_FORWARD0SWEEP_TRACE
00583           size_t       d      = 0;
00584           size_t       i_tmp  = i_var;
00585           Base*        Z_tmp  = Taylor + i_var * J;
00586 
00587           printOp(
00588                std::cout, 
00589                Rec,
00590                i_tmp,
00591                op, 
00592                arg,
00593                d + 1, 
00594                Z_tmp, 
00595                0, 
00596                (Base *) CPPAD_NULL
00597           );
00598      }
00599      std::cout << std::endl;
00600 # else
00601      }
00602 # endif
00603      CPPAD_ASSERT_UNKNOWN( user_state == user_start );
00604      CPPAD_ASSERT_UNKNOWN( i_var + 1 == Rec->num_rec_var() );
00605 
00606      return compareCount;
00607 }
00608 
00609 /*! \} */
00610 CPPAD_END_NAMESPACE
00611 
00612 // preprocessor symbols that are local to this file
00613 # undef CPPAD_FORWARD0SWEEP_TRACE
00614 
00615 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines