CppAD: A C++ Algorithmic Differentiation Package  20130918
forward.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_FORWARD_INCLUDED
00003 # define CPPAD_FORWARD_INCLUDED
00004 
00005 /* --------------------------------------------------------------------------
00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 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 // documened after Forward but included here so easy to see
00017 # include <cppad/local/capacity_order.hpp>
00018 # include <cppad/local/num_skip.hpp>
00019 
00020 namespace CppAD { // BEGIN_CPPAD_NAMESPACE
00021 /*!
00022 \file forward.hpp
00023 User interface to forward mode computations
00024 */
00025 
00026 /*!
00027 Arbitrary order, one direction, forward mode Taylor coefficieints.
00028 
00029 \tparam Base
00030 The type used during the forward mode computations; i.e., the corresponding
00031 recording of operations used the type AD<Base>.
00032 
00033 \tparam Vector
00034 is a Simple Vector class with eleements of type Base.
00035 
00036 \param q
00037 is the hightest order for this forward mode computation; i.e., 
00038 after this calculation there will be <code>q+1</code>
00039 Taylor coefficients per variable.
00040 
00041 \param xq
00042 contains Taylor coefficients for the independent variables.
00043 The size of xq must either be n or <code>(q+1)*n</code>,
00044 We define <code>p = q + 1 - xq.size()/n</code>. 
00045 For <code>j = 0 , ... , n-1</code>,
00046 <code>k = p, ... , q</code>, are
00047 <code>xq[ (q+1-p)*j + k - p ]</code>
00048 is the k-th order coefficient for the j-th independent variable.
00049 
00050 \param s
00051 Is the stream where output corresponding to PriOp operations will written.
00052 
00053 \return
00054 contains Taylor coefficients for the independent variables.
00055 The size of the return value y is <code>m*(q+1-p)</code>.
00056 For <code>i = 0, ... , m-1</code>,
00057 <code>k = p, ..., q</code>,
00058 <code>y[(q+1-p)*i + (k-p)]</code> 
00059 is the k-th order coefficient for the i-th dependent variable.
00060 
00061 \par taylor_
00062 The Taylor coefficients up to order p-1 are inputs
00063 and the coefficents from order p through q are outputs.
00064 Let <code>N = num_var_tape_</code>, and
00065 <code>C = cap_order_taylor_</code>.
00066 Note that for
00067 <code>i = 1 , ..., N-1</code>,
00068 <code>k = 0 , ..., q</code>,
00069 <code>taylor_[ C*i + k ]</code>
00070 is the k-th order cofficent,
00071 for the i-th varaible on the tape.
00072 (The first independent variable has index one on the tape 
00073 and there is no variable with index zero.)
00074 */
00075 
00076 template <typename Base>
00077 template <typename Vector>
00078 Vector ADFun<Base>::Forward(
00079      size_t q                    , 
00080      const Vector& xq            , 
00081      std::ostream& s             )
00082 {    // temporary indices
00083      size_t i, j, k;
00084 
00085      // number of independent variables
00086      size_t n = ind_taddr_.size();
00087 
00088      // number of dependent variables
00089      size_t m = dep_taddr_.size();
00090 
00091      // check Vector is Simple Vector class with Base type elements
00092      CheckSimpleVector<Base, Vector>();
00093 
00094 
00095      CPPAD_ASSERT_KNOWN(
00096           size_t(xq.size()) == n || size_t(xq.size()) == n*(q+1),
00097           "Forward(q, xq): xq.size() is not equal n or n*(q+1)"
00098      );
00099 
00100      // lowest order we are computing
00101      size_t p = q + 1 - size_t(xq.size()) / n;
00102      CPPAD_ASSERT_UNKNOWN( p == 0 || p == q );
00103      CPPAD_ASSERT_KNOWN(
00104           q <= num_order_taylor_ || p == 0,
00105           "Forward(q, xq): Number of Taylor coefficient orders stored in this"
00106           " ADFun\nis less than q and xq.size() != n*(q+1)."
00107      );  
00108      // does taylor_ need more orders or different number of directions
00109      if( cap_order_taylor_ <= q )
00110      {    if( p == 0 )
00111           {    // no need to copy old values during capacity_order
00112                num_order_taylor_ = 0;
00113           }
00114           size_t c = std::max(q + 1, cap_order_taylor_);
00115           capacity_order(c);
00116      }
00117      CPPAD_ASSERT_UNKNOWN( cap_order_taylor_ > q );
00118 
00119      // short hand notation for order capacity
00120      size_t C = cap_order_taylor_;
00121 
00122      // set Taylor coefficients for independent variables
00123      for(j = 0; j < n; j++)
00124      {    CPPAD_ASSERT_UNKNOWN( ind_taddr_[j] < num_var_tape_  );
00125 
00126           // ind_taddr_[j] is operator taddr for j-th independent variable
00127           CPPAD_ASSERT_UNKNOWN( play_.GetOp( ind_taddr_[j] ) == InvOp );
00128 
00129           if( p ==  q )
00130                taylor_[ C * ind_taddr_[j] + q] = xq[j];
00131           else
00132           {    for(k = 0; k <= q; k++)
00133                     taylor_[ C * ind_taddr_[j] + k] = xq[ (q+1)*j + k];
00134           }
00135      }
00136 
00137      // evaluate the derivatives
00138      CPPAD_ASSERT_UNKNOWN( cskip_op_.size() == play_.num_op_rec() );
00139      CPPAD_ASSERT_UNKNOWN( load_op_.size()  == play_.num_load_op_rec() );
00140      if( q == 0 )
00141      {
00142           compare_change_ = forward0sweep(s, true,
00143                n, num_var_tape_, &play_, C, 
00144                taylor_.data(), cskip_op_.data(), load_op_
00145           );
00146      }
00147      else if( p == 0 )
00148      {    compare_change_ = forward1sweep(s, true, p, q, 
00149                n, num_var_tape_, &play_, C, 
00150                taylor_.data(), cskip_op_.data(), load_op_
00151           );
00152      }
00153      else
00154      {    forward1sweep(s, true, p, q, 
00155                n, num_var_tape_, &play_, C, 
00156                taylor_.data(), cskip_op_.data(), load_op_
00157           );
00158      }
00159 
00160      // return Taylor coefficients for dependent variables
00161      Vector yq;
00162      if( p == q )
00163      {    yq.resize(m);
00164           for(i = 0; i < m; i++)
00165           {    CPPAD_ASSERT_UNKNOWN( dep_taddr_[i] < num_var_tape_  );
00166                yq[i] = taylor_[ C * dep_taddr_[i] + q];
00167           }
00168      }
00169      else
00170      {    yq.resize(m * (q+1) );
00171           for(i = 0; i < m; i++)   
00172           {    for(k = 0; k <= q; k++)
00173                     yq[ (q+1) * i + k] = 
00174                          taylor_[ C * dep_taddr_[i] + k ]; 
00175           }
00176      }
00177 # ifndef NDEBUG
00178      if( check_for_nan_ )
00179      {    bool ok = true;
00180           if( p == 0 )
00181           {    for(i = 0; i < m; i++)
00182                {    // Visual Studio 2012, CppAD required in front of isnan ?
00183                     ok &= ! CppAD::isnan( yq[ (q+1) * i + 0 ] );
00184                }
00185           } 
00186           CPPAD_ASSERT_KNOWN(ok,
00187                "yq = f.Forward(q, xq): has a zero order Taylor coefficient "
00188                "with the value nan."
00189           );  
00190           if( 0 < q )
00191           {    for(i = 0; i < m; i++)
00192                {    for(k = p; k <= q; k++)
00193                     {    // Studio 2012, CppAD required in front of isnan ?
00194                          ok &= ! CppAD::isnan( yq[ (q+1-p)*i + k-p ] );
00195                     }
00196                }
00197           }
00198           CPPAD_ASSERT_KNOWN(ok,
00199           "yq = f.Forward(q, xq): has a non-zero order Taylor coefficient\n"
00200           "with the value nan (but zero order coefficients are not nan)."
00201           );
00202      }
00203 # endif
00204 
00205      // now we have q + 1  taylor_ coefficient orders per variable
00206      num_order_taylor_ = q + 1;
00207 
00208      return yq;
00209 }
00210 
00211 
00212 } // END_CPPAD_NAMESPACE
00213 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines