Sacado Package Browser (Single Doxygen Collection) Version of the Day
Sacado_Fad_Ops.hpp
Go to the documentation of this file.
00001 // $Id$ 
00002 // $Source$ 
00003 // @HEADER
00004 // ***********************************************************************
00005 // 
00006 //                           Sacado Package
00007 //                 Copyright (2006) Sandia Corporation
00008 // 
00009 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00010 // the U.S. Government retains certain rights in this software.
00011 // 
00012 // This library is free software; you can redistribute it and/or modify
00013 // it under the terms of the GNU Lesser General Public License as
00014 // published by the Free Software Foundation; either version 2.1 of the
00015 // License, or (at your option) any later version.
00016 //  
00017 // This library is distributed in the hope that it will be useful, but
00018 // WITHOUT ANY WARRANTY; without even the implied warranty of
00019 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00020 // Lesser General Public License for more details.
00021 //  
00022 // You should have received a copy of the GNU Lesser General Public
00023 // License along with this library; if not, write to the Free Software
00024 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00025 // USA
00026 // Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
00027 // (etphipp@sandia.gov).
00028 // 
00029 // ***********************************************************************
00030 //
00031 // The forward-mode AD classes in Sacado are a derivative work of the
00032 // expression template classes in the Fad package by Nicolas Di Cesare.  
00033 // The following banner is included in the original Fad source code:
00034 //
00035 // ************ DO NOT REMOVE THIS BANNER ****************
00036 //
00037 //  Nicolas Di Cesare <Nicolas.Dicesare@ann.jussieu.fr>
00038 //  http://www.ann.jussieu.fr/~dicesare
00039 //
00040 //            CEMRACS 98 : C++ courses, 
00041 //         templates : new C++ techniques 
00042 //            for scientific computing 
00043 // 
00044 //********************************************************
00045 //
00046 //  A short implementation ( not all operators and 
00047 //  functions are overloaded ) of 1st order Automatic
00048 //  Differentiation in forward mode (FAD) using
00049 //  EXPRESSION TEMPLATES.
00050 //
00051 //********************************************************
00052 // @HEADER
00053 
00054 #ifndef SACADO_FAD_OPS_HPP
00055 #define SACADO_FAD_OPS_HPP
00056 
00057 #include "Sacado_Fad_Expression.hpp"
00058 #include "Sacado_cmath.hpp"
00059 #include <ostream>  // for std::ostream
00060 
00061 #define FAD_UNARYOP_MACRO(OPNAME,OP,VALUE,DX,FASTACCESSDX)    \
00062 namespace Sacado {              \
00063   namespace Fad {             \
00064                   \
00065     template <typename ExprT>           \
00066     class OP {};              \
00067                   \
00068     template <typename ExprT>           \
00069     class Expr< OP<ExprT> > {           \
00070     public:               \
00071                   \
00072       typedef typename ExprT::value_type value_type;      \
00073                   \
00074       Expr(const ExprT& expr_) : expr(expr_)  {}      \
00075                   \
00076       int size() const { return expr.size(); }        \
00077                   \
00078       bool hasFastAccess() const { return expr.hasFastAccess(); } \
00079                   \
00080       bool isPassive() const { return expr.isPassive();}    \
00081                   \
00082       bool updateValue() const { return expr.updateValue(); }   \
00083                         \
00084       value_type val() const {            \
00085   return VALUE;             \
00086       }                 \
00087                   \
00088       value_type dx(int i) const {          \
00089   return DX;              \
00090       }                 \
00091                   \
00092       value_type fastAccessDx(int i) const {        \
00093   return FASTACCESSDX;            \
00094       }                 \
00095                   \
00096     protected:                \
00097                   \
00098       const ExprT& expr;            \
00099     };                  \
00100                   \
00101     template <typename T>           \
00102     inline Expr< OP< Expr<T> > >          \
00103     OPNAME (const Expr<T>& expr)          \
00104     {                 \
00105       typedef OP< Expr<T> > expr_t;         \
00106                         \
00107       return Expr<expr_t>(expr);          \
00108     }                 \
00109   }                 \
00110 }
00111 
00112 FAD_UNARYOP_MACRO(operator+,
00113       UnaryPlusOp, 
00114       expr.val(),
00115       expr.dx(i),
00116       expr.fastAccessDx(i))
00117 FAD_UNARYOP_MACRO(operator-,
00118       UnaryMinusOp, 
00119       -expr.val(),
00120       -expr.dx(i),
00121       -expr.fastAccessDx(i))
00122 FAD_UNARYOP_MACRO(exp,
00123       ExpOp, 
00124       std::exp(expr.val()),
00125       std::exp(expr.val())*expr.dx(i),
00126       std::exp(expr.val())*expr.fastAccessDx(i))
00127 FAD_UNARYOP_MACRO(log,
00128       LogOp, 
00129       std::log(expr.val()),
00130       expr.dx(i)/expr.val(),
00131       expr.fastAccessDx(i)/expr.val())
00132 FAD_UNARYOP_MACRO(log10,
00133       Log10Op, 
00134       std::log10(expr.val()),
00135       expr.dx(i)/( std::log(value_type(10))*expr.val()),
00136       expr.fastAccessDx(i) / ( std::log(value_type(10))*expr.val()))
00137 FAD_UNARYOP_MACRO(sqrt,
00138       SqrtOp, 
00139       std::sqrt(expr.val()),
00140       expr.dx(i)/(value_type(2)* std::sqrt(expr.val())),
00141       expr.fastAccessDx(i)/(value_type(2)* std::sqrt(expr.val())))
00142 FAD_UNARYOP_MACRO(cos,
00143       CosOp, 
00144       std::cos(expr.val()),
00145       -expr.dx(i)* std::sin(expr.val()),
00146       -expr.fastAccessDx(i)* std::sin(expr.val()))
00147 FAD_UNARYOP_MACRO(sin,
00148       SinOp, 
00149       std::sin(expr.val()),
00150       expr.dx(i)* std::cos(expr.val()),
00151       expr.fastAccessDx(i)* std::cos(expr.val()))
00152 FAD_UNARYOP_MACRO(tan,
00153       TanOp, 
00154       std::tan(expr.val()),
00155       expr.dx(i)*
00156         (value_type(1)+ std::tan(expr.val())* std::tan(expr.val())),
00157       expr.fastAccessDx(i)*
00158         (value_type(1)+ std::tan(expr.val())* std::tan(expr.val())))
00159 FAD_UNARYOP_MACRO(acos,
00160       ACosOp, 
00161       std::acos(expr.val()),
00162       -expr.dx(i)/ std::sqrt(value_type(1)-expr.val()*expr.val()),
00163       -expr.fastAccessDx(i) /
00164         std::sqrt(value_type(1)-expr.val()*expr.val()))
00165 FAD_UNARYOP_MACRO(asin,
00166       ASinOp, 
00167       std::asin(expr.val()),
00168       expr.dx(i)/ std::sqrt(value_type(1)-expr.val()*expr.val()),
00169       expr.fastAccessDx(i) /
00170         std::sqrt(value_type(1)-expr.val()*expr.val()))
00171 FAD_UNARYOP_MACRO(atan,
00172       ATanOp, 
00173       std::atan(expr.val()),
00174       expr.dx(i)/(value_type(1)+expr.val()*expr.val()),
00175       expr.fastAccessDx(i)/(value_type(1)+expr.val()*expr.val()))
00176 FAD_UNARYOP_MACRO(cosh,
00177       CoshOp, 
00178       std::cosh(expr.val()),
00179       expr.dx(i)* std::sinh(expr.val()),
00180       expr.fastAccessDx(i)* std::sinh(expr.val()))
00181 FAD_UNARYOP_MACRO(sinh,
00182       SinhOp, 
00183       std::sinh(expr.val()),
00184       expr.dx(i)* std::cosh(expr.val()),
00185       expr.fastAccessDx(i)* std::cosh(expr.val()))
00186 FAD_UNARYOP_MACRO(tanh,
00187       TanhOp, 
00188       std::tanh(expr.val()),
00189       expr.dx(i)/( std::cosh(expr.val())* std::cosh(expr.val())),
00190       expr.fastAccessDx(i) / 
00191         ( std::cosh(expr.val())* std::cosh(expr.val())))
00192 FAD_UNARYOP_MACRO(acosh,
00193       ACoshOp, 
00194       std::acosh(expr.val()),
00195       expr.dx(i)/ std::sqrt((expr.val()-value_type(1)) * 
00196                (expr.val()+value_type(1))),
00197       expr.fastAccessDx(i)/ std::sqrt((expr.val()-value_type(1)) * 
00198              (expr.val()+value_type(1))))
00199 FAD_UNARYOP_MACRO(asinh,
00200       ASinhOp, 
00201       std::asinh(expr.val()),
00202       expr.dx(i)/ std::sqrt(value_type(1)+expr.val()*expr.val()),
00203       expr.fastAccessDx(i)/ std::sqrt(value_type(1)+
00204              expr.val()*expr.val()))
00205 FAD_UNARYOP_MACRO(atanh,
00206       ATanhOp, 
00207       std::atanh(expr.val()),
00208       expr.dx(i)/(value_type(1)-expr.val()*expr.val()),
00209       expr.fastAccessDx(i)/(value_type(1)-
00210              expr.val()*expr.val()))
00211 FAD_UNARYOP_MACRO(abs,
00212       AbsOp, 
00213       std::abs(expr.val()),
00214       expr.val() >= 0 ? value_type(+expr.dx(i)) : 
00215         value_type(-expr.dx(i)),
00216       expr.val() >= 0 ? value_type(+expr.fastAccessDx(i)) : 
00217         value_type(-expr.fastAccessDx(i)))
00218 FAD_UNARYOP_MACRO(fabs,
00219       FAbsOp, 
00220       std::fabs(expr.val()),
00221       expr.val() >= 0 ? value_type(+expr.dx(i)) : 
00222         value_type(-expr.dx(i)),
00223       expr.val() >= 0 ? value_type(+expr.fastAccessDx(i)) : 
00224         value_type(-expr.fastAccessDx(i)))
00225 
00226 #undef FAD_UNARYOP_MACRO
00227 
00228 #define FAD_BINARYOP_MACRO(OPNAME,OP,VALUE,DX,FASTACCESSDX,VAL_CONST_DX_1,VAL_CONST_DX_2,CONST_DX_1,CONST_DX_2,CONST_FASTACCESSDX_1,CONST_FASTACCESSDX_2) \
00229 namespace Sacado {              \
00230   namespace Fad {             \
00231                   \
00232     template <typename ExprT1, typename ExprT2>       \
00233     class OP {};              \
00234                   \
00235     template <typename T1, typename T2>         \
00236     class Expr< OP< Expr<T1>, Expr<T2> > > {        \
00237                   \
00238     public:               \
00239                   \
00240       typedef Expr<T1> ExprT1;            \
00241       typedef Expr<T2> ExprT2;            \
00242       typedef typename ExprT1::value_type value_type_1;     \
00243       typedef typename ExprT2::value_type value_type_2;     \
00244       typedef typename Sacado::Promote<value_type_1,      \
00245                value_type_2>::type value_type;  \
00246                   \
00247       Expr(const ExprT1& expr1_, const ExprT2& expr2_) :    \
00248   expr1(expr1_), expr2(expr2_) {}         \
00249                   \
00250       int size() const {            \
00251   int sz1 = expr1.size(), sz2 = expr2.size();     \
00252   return sz1 > sz2 ? sz1 : sz2;         \
00253       }                 \
00254                   \
00255       bool hasFastAccess() const {          \
00256   return expr1.hasFastAccess() && expr2.hasFastAccess();    \
00257       }                 \
00258                   \
00259       bool isPassive() const {            \
00260   return expr1.isPassive() && expr2.isPassive();      \
00261       }                 \
00262                   \
00263       bool updateValue() const {          \
00264   return expr1.updateValue() && expr2.updateValue();    \
00265       }                 \
00266                   \
00267       const value_type val() const {          \
00268   return VALUE;             \
00269       }                 \
00270                   \
00271       const value_type dx(int i) const {        \
00272   return DX;              \
00273       }                 \
00274                   \
00275       const value_type fastAccessDx(int i) const {      \
00276   return FASTACCESSDX;            \
00277       }                 \
00278                         \
00279     protected:                \
00280                   \
00281       const ExprT1& expr1;            \
00282       const ExprT2& expr2;            \
00283                   \
00284     };                  \
00285                   \
00286     template <typename T1>            \
00287     class Expr< OP< Expr<T1>, typename Expr<T1>::value_type> > {  \
00288                   \
00289     public:               \
00290                   \
00291       typedef Expr<T1> ExprT1;            \
00292       typedef typename ExprT1::value_type value_type;     \
00293       typedef typename ExprT1::value_type ConstT;     \
00294                   \
00295       Expr(const ExprT1& expr1_, const ConstT& c_) :      \
00296   expr1(expr1_), c(c_) {}           \
00297                   \
00298       int size() const {            \
00299   return expr1.size();            \
00300       }                 \
00301                   \
00302       bool hasFastAccess() const {          \
00303   return expr1.hasFastAccess();         \
00304       }                 \
00305                   \
00306       bool isPassive() const {            \
00307   return expr1.isPassive();         \
00308       }                 \
00309                   \
00310       bool updateValue() const { return expr1.updateValue(); }    \
00311                   \
00312       const value_type val() const {          \
00313   return VAL_CONST_DX_2;            \
00314       }                 \
00315                   \
00316       const value_type dx(int i) const {        \
00317   return CONST_DX_2;            \
00318       }                 \
00319                   \
00320       const value_type fastAccessDx(int i) const {      \
00321   return CONST_FASTACCESSDX_2;          \
00322       }                 \
00323                   \
00324     protected:                \
00325                   \
00326       const ExprT1& expr1;            \
00327       const ConstT& c;              \
00328     };                  \
00329                   \
00330     template <typename T2>            \
00331     class Expr< OP< typename Expr<T2>::value_type, Expr<T2> > > { \
00332                   \
00333     public:               \
00334                   \
00335       typedef Expr<T2> ExprT2;            \
00336       typedef typename ExprT2::value_type value_type;     \
00337       typedef typename ExprT2::value_type ConstT;     \
00338                   \
00339       Expr(const ConstT& c_, const ExprT2& expr2_) :      \
00340   c(c_), expr2(expr2_) {}           \
00341                   \
00342       int size() const {            \
00343   return expr2.size();            \
00344       }                 \
00345                   \
00346       bool hasFastAccess() const {          \
00347   return expr2.hasFastAccess();         \
00348       }                 \
00349                   \
00350       bool isPassive() const {            \
00351   return expr2.isPassive();         \
00352       }                 \
00353                   \
00354       bool updateValue() const { return expr2.updateValue(); }    \
00355                   \
00356       const value_type val() const {          \
00357   return VAL_CONST_DX_1;            \
00358       }                 \
00359                   \
00360       const value_type dx(int i) const {        \
00361   return CONST_DX_1;            \
00362       }                 \
00363                   \
00364       const value_type fastAccessDx(int i) const {      \
00365   return CONST_FASTACCESSDX_1;          \
00366       }                 \
00367                         \
00368     protected:                \
00369                   \
00370       const ConstT& c;              \
00371       const ExprT2& expr2;            \
00372     };                  \
00373                   \
00374     template <typename T1, typename T2>         \
00375     inline Expr< OP< Expr<T1>, Expr<T2> > >       \
00376     OPNAME (const Expr<T1>& expr1, const Expr<T2>& expr2)   \
00377     {                 \
00378       typedef OP< Expr<T1>, Expr<T2> > expr_t;        \
00379                       \
00380       return Expr<expr_t>(expr1, expr2);        \
00381     }                 \
00382                   \
00383     template <typename T>           \
00384     inline Expr< OP< Expr<T>, Expr<T> > >       \
00385     OPNAME (const Expr<T>& expr1, const Expr<T>& expr2)     \
00386     {                 \
00387       typedef OP< Expr<T>, Expr<T> > expr_t;        \
00388                       \
00389       return Expr<expr_t>(expr1, expr2);        \
00390     }                 \
00391                   \
00392     template <typename T>           \
00393     inline Expr< OP< typename Expr<T>::value_type, Expr<T> > >    \
00394     OPNAME (const typename Expr<T>::value_type& c,      \
00395       const Expr<T>& expr)          \
00396     {                 \
00397       typedef typename Expr<T>::value_type ConstT;      \
00398       typedef OP< ConstT, Expr<T> > expr_t;       \
00399                   \
00400       return Expr<expr_t>(c, expr);         \
00401     }                 \
00402                   \
00403     template <typename T>           \
00404     inline Expr< OP< Expr<T>, typename Expr<T>::value_type > >    \
00405     OPNAME (const Expr<T>& expr,          \
00406       const typename Expr<T>::value_type& c)      \
00407     {                 \
00408       typedef typename Expr<T>::value_type ConstT;      \
00409       typedef OP< Expr<T>, ConstT > expr_t;       \
00410                   \
00411       return Expr<expr_t>(expr, c);         \
00412     }                 \
00413   }                 \
00414 }
00415 
00416 
00417 FAD_BINARYOP_MACRO(operator+,
00418        AdditionOp, 
00419        expr1.val() + expr2.val(),
00420        expr1.dx(i) + expr2.dx(i),
00421        expr1.fastAccessDx(i) + expr2.fastAccessDx(i),
00422        c + expr2.val(),
00423        expr1.val() + c,
00424        expr2.dx(i),
00425        expr1.dx(i),
00426        expr2.fastAccessDx(i),
00427        expr1.fastAccessDx(i))
00428 FAD_BINARYOP_MACRO(operator-,
00429        SubtractionOp, 
00430        expr1.val() - expr2.val(),
00431        expr1.dx(i) - expr2.dx(i),
00432        expr1.fastAccessDx(i) - expr2.fastAccessDx(i),
00433        c - expr2.val(),
00434        expr1.val() - c,
00435        -expr2.dx(i),
00436        expr1.dx(i),
00437        -expr2.fastAccessDx(i),
00438        expr1.fastAccessDx(i))
00439 // FAD_BINARYOP_MACRO(operator*,
00440 //       MultiplicationOp, 
00441 //       expr1.val() * expr2.val(),
00442 //       expr1.val()*expr2.dx(i) + expr1.dx(i)*expr2.val(),
00443 //       expr1.val()*expr2.fastAccessDx(i) + 
00444 //         expr1.fastAccessDx(i)*expr2.val(),
00445 //       c * expr2.val(),
00446 //       expr1.val() * c,
00447 //       c*expr2.dx(i),
00448 //       expr1.dx(i)*c,
00449 //       c*expr2.fastAccessDx(i),
00450 //       expr1.fastAccessDx(i)*c)
00451 FAD_BINARYOP_MACRO(operator/,
00452        DivisionOp, 
00453        expr1.val() / expr2.val(),
00454        (expr1.dx(i)*expr2.val() - expr2.dx(i)*expr1.val()) /
00455          (expr2.val()*expr2.val()),
00456        (expr1.fastAccessDx(i)*expr2.val() - 
00457           expr2.fastAccessDx(i)*expr1.val()) /
00458           (expr2.val()*expr2.val()),
00459        c / expr2.val(),
00460        expr1.val() / c,
00461        -expr2.dx(i)*c / (expr2.val()*expr2.val()),
00462        expr1.dx(i)/c,
00463        -expr2.fastAccessDx(i)*c / (expr2.val()*expr2.val()),
00464        expr1.fastAccessDx(i)/c)
00465 FAD_BINARYOP_MACRO(atan2,
00466        Atan2Op,
00467        std::atan2(expr1.val(), expr2.val()),
00468        (expr2.val()*expr1.dx(i) - expr1.val()*expr2.dx(i))/
00469       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00470        (expr2.val()*expr1.fastAccessDx(i) - expr1.val()*expr2.fastAccessDx(i))/
00471       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00472        std::atan2(c, expr2.val()),
00473        std::atan2(expr1.val(), c),
00474        (-c*expr2.dx(i)) / (c*c + expr2.val()*expr2.val()),
00475        (c*expr1.dx(i))/ (expr1.val()*expr1.val() + c*c),
00476        (-c*expr2.fastAccessDx(i))/ (c*c + expr2.val()*expr2.val()),
00477        (c*expr1.fastAccessDx(i))/ (expr1.val()*expr1.val() + c*c))
00478 FAD_BINARYOP_MACRO(pow,
00479        PowerOp,
00480        std::pow(expr1.val(), expr2.val()),
00481        expr1.val() == value_type(0) ? value_type(0) : value_type((expr2.dx(i)*std::log(expr1.val())+expr2.val()*expr1.dx(i)/expr1.val())*std::pow(expr1.val(),expr2.val())),
00482        expr1.val() == value_type(0) ? value_type(0.0) : value_type((expr2.fastAccessDx(i)*std::log(expr1.val())+expr2.val()*expr1.fastAccessDx(i)/expr1.val())*std::pow(expr1.val(),expr2.val())),
00483        std::pow(c, expr2.val()),
00484        std::pow(expr1.val(), c),
00485        c == value_type(0) ? value_type(0) : value_type(expr2.dx(i)*std::log(c)*std::pow(c,expr2.val())),
00486        expr1.val() == value_type(0) ? value_type(0.0) : value_type(c*expr1.dx(i)/expr1.val()*std::pow(expr1.val(),c)),
00487        c == value_type(0) ? value_type(0) : value_type(expr2.fastAccessDx(i)*std::log(c)*std::pow(c,expr2.val())),
00488        expr1.val() == value_type(0) ? value_type(0.0) : value_type(c*expr1.fastAccessDx(i)/expr1.val()*std::pow(expr1.val(),c)))
00489 FAD_BINARYOP_MACRO(max,
00490                    MaxOp,
00491                    std::max(expr1.val(), expr2.val()),
00492                    expr1.val() >= expr2.val() ? expr1.dx(i) : expr2.dx(i),
00493                    expr1.val() >= expr2.val() ? expr1.fastAccessDx(i) : 
00494                                                 expr2.fastAccessDx(i),
00495        std::max(c, expr2.val()),
00496        std::max(expr1.val(), c),
00497                    c >= expr2.val() ? value_type(0) : expr2.dx(i),
00498                    expr1.val() >= c ? expr1.dx(i) : value_type(0),
00499                    c >= expr2.val() ? value_type(0) : expr2.fastAccessDx(i),
00500                    expr1.val() >= c ? expr1.fastAccessDx(i) : value_type(0))
00501 FAD_BINARYOP_MACRO(min,
00502                    MinOp,
00503                    std::min(expr1.val(), expr2.val()),
00504                    expr1.val() <= expr2.val() ? expr1.dx(i) : expr2.dx(i),
00505                    expr1.val() <= expr2.val() ? expr1.fastAccessDx(i) : 
00506                                                 expr2.fastAccessDx(i),
00507        std::min(c, expr2.val()),
00508        std::min(expr1.val(), c),
00509                    c <= expr2.val() ? value_type(0) : expr2.dx(i),
00510                    expr1.val() <= c ? expr1.dx(i) : value_type(0),
00511                    c <= expr2.val() ? value_type(0) : expr2.fastAccessDx(i),
00512                    expr1.val() <= c ? expr1.fastAccessDx(i) : value_type(0))
00513 
00514 
00515 #undef FAD_BINARYOP_MACRO
00516 
00517 namespace Sacado {              
00518   namespace Fad {             
00519                   
00520     template <typename ExprT1, typename ExprT2>       
00521     class MultiplicationOp {};              
00522                   
00523     template <typename T1, typename T2>         
00524     class Expr< MultiplicationOp< Expr<T1>, Expr<T2> > > {        
00525                   
00526     public:               
00527                   
00528       typedef Expr<T1> ExprT1;            
00529       typedef Expr<T2> ExprT2;            
00530       typedef typename ExprT1::value_type value_type_1;     
00531       typedef typename ExprT2::value_type value_type_2;     
00532       typedef typename Sacado::Promote<value_type_1,      
00533                value_type_2>::type value_type;  
00534                   
00535       Expr(const ExprT1& expr1_, const ExprT2& expr2_) :    
00536   expr1(expr1_), expr2(expr2_) {}         
00537                   
00538       int size() const {            
00539   int sz1 = expr1.size(), sz2 = expr2.size();     
00540   return sz1 > sz2 ? sz1 : sz2;         
00541       }                 
00542                   
00543       bool hasFastAccess() const {          
00544   return expr1.hasFastAccess() && expr2.hasFastAccess();    
00545       }                 
00546                   
00547       bool isPassive() const {            
00548   return expr1.isPassive() && expr2.isPassive();      
00549       }
00550 
00551       bool updateValue() const { 
00552   return expr1.updateValue() && expr2.updateValue(); 
00553       }
00554 
00555       const value_type val() const {          
00556   return expr1.val()*expr2.val();
00557       }                 
00558                   
00559       const value_type dx(int i) const {
00560   if (expr1.size() > 0 && expr2.size() > 0)
00561     return expr1.val()*expr2.dx(i) + expr1.dx(i)*expr2.val();
00562   else if (expr1.size() > 0)
00563     return expr1.dx(i)*expr2.val();
00564   else
00565     return expr1.val()*expr2.dx(i);
00566       }                 
00567                   
00568       const value_type fastAccessDx(int i) const {      
00569   return expr1.val()*expr2.fastAccessDx(i) + 
00570     expr1.fastAccessDx(i)*expr2.val();
00571       }                 
00572                         
00573     protected:                
00574                   
00575       const ExprT1& expr1;            
00576       const ExprT2& expr2;            
00577                   
00578     };                  
00579                   
00580     template <typename T1>            
00581     class Expr< MultiplicationOp< Expr<T1>, typename Expr<T1>::value_type> > {  
00582                   
00583     public:               
00584                   
00585       typedef Expr<T1> ExprT1;            
00586       typedef typename ExprT1::value_type value_type;     
00587       typedef typename ExprT1::value_type ConstT;     
00588                   
00589       Expr(const ExprT1& expr1_, const ConstT& c_) :      
00590   expr1(expr1_), c(c_) {}           
00591                   
00592       int size() const {            
00593   return expr1.size();            
00594       }                 
00595                   
00596       bool hasFastAccess() const {          
00597   return expr1.hasFastAccess();         
00598       }                 
00599                   
00600       bool isPassive() const {            
00601   return expr1.isPassive();         
00602       }
00603 
00604       bool updateValue() const { return expr1.updateValue(); }
00605                   
00606       const value_type val() const {          
00607   return expr1.val()*c;           
00608       }                 
00609                   
00610       const value_type dx(int i) const {        
00611   return expr1.dx(i)*c;           
00612       }                 
00613                   
00614       const value_type fastAccessDx(int i) const {      
00615   return expr1.fastAccessDx(i)*c;         
00616       }                 
00617                   
00618     protected:                
00619                   
00620       const ExprT1& expr1;            
00621       const ConstT& c;              
00622     };                  
00623                   
00624     template <typename T2>            
00625     class Expr< MultiplicationOp< typename Expr<T2>::value_type, Expr<T2> > > { 
00626                   
00627     public:               
00628                   
00629       typedef Expr<T2> ExprT2;            
00630       typedef typename ExprT2::value_type value_type;     
00631       typedef typename ExprT2::value_type ConstT;     
00632                   
00633       Expr(const ConstT& c_, const ExprT2& expr2_) :      
00634   c(c_), expr2(expr2_) {}           
00635                   
00636       int size() const {            
00637   return expr2.size();            
00638       }                 
00639                   
00640       bool hasFastAccess() const {          
00641   return expr2.hasFastAccess();         
00642       }                 
00643                   
00644       bool isPassive() const {            
00645   return expr2.isPassive();         
00646       }
00647 
00648       bool updateValue() const { return expr2.updateValue(); }
00649                   
00650       const value_type val() const {          
00651   return c*expr2.val();           
00652       }                 
00653                   
00654       const value_type dx(int i) const {        
00655   return c*expr2.dx(i);           
00656       }                 
00657                   
00658       const value_type fastAccessDx(int i) const {      
00659   return c*expr2.fastAccessDx(i);         
00660       }                 
00661                         
00662     protected:                
00663                   
00664       const ConstT& c;              
00665       const ExprT2& expr2;            
00666     };                  
00667                   
00668     template <typename T1, typename T2>         
00669     inline Expr< MultiplicationOp< Expr<T1>, Expr<T2> > >       
00670     operator* (const Expr<T1>& expr1, const Expr<T2>& expr2)    
00671     {                 
00672       typedef MultiplicationOp< Expr<T1>, Expr<T2> > expr_t;        
00673                       
00674       return Expr<expr_t>(expr1, expr2);        
00675     }                 
00676                   
00677     template <typename T>           
00678     inline Expr< MultiplicationOp< Expr<T>, Expr<T> > >       
00679     operator* (const Expr<T>& expr1, const Expr<T>& expr2)      
00680     {                 
00681       typedef MultiplicationOp< Expr<T>, Expr<T> > expr_t;        
00682                       
00683       return Expr<expr_t>(expr1, expr2);        
00684     }                 
00685                   
00686     template <typename T>           
00687     inline Expr< MultiplicationOp< typename Expr<T>::value_type, Expr<T> > >    
00688     operator* (const typename Expr<T>::value_type& c,     
00689       const Expr<T>& expr)          
00690     {                 
00691       typedef typename Expr<T>::value_type ConstT;      
00692       typedef MultiplicationOp< ConstT, Expr<T> > expr_t;       
00693                   
00694       return Expr<expr_t>(c, expr);         
00695     }                 
00696                   
00697     template <typename T>           
00698     inline Expr< MultiplicationOp< Expr<T>, typename Expr<T>::value_type > >    
00699     operator* (const Expr<T>& expr,         
00700       const typename Expr<T>::value_type& c)      
00701     {                 
00702       typedef typename Expr<T>::value_type ConstT;      
00703       typedef MultiplicationOp< Expr<T>, ConstT > expr_t;       
00704                   
00705       return Expr<expr_t>(expr, c);         
00706     }                 
00707   }                 
00708 }
00709 
00710 //-------------------------- Relational Operators -----------------------
00711 
00712 #define FAD_RELOP_MACRO(OP)           \
00713 namespace Sacado {              \
00714   namespace Fad {             \
00715     template <typename ExprT1, typename ExprT2>       \
00716     inline bool               \
00717     operator OP (const Expr<ExprT1>& expr1,       \
00718      const Expr<ExprT2>& expr2)       \
00719     {                 \
00720       return expr1.val() OP expr2.val();        \
00721     }                 \
00722                   \
00723     template <typename ExprT2>            \
00724     inline bool               \
00725     operator OP (const typename Expr<ExprT2>::value_type& a,    \
00726      const Expr<ExprT2>& expr2)       \
00727     {                 \
00728       return a OP expr2.val();            \
00729     }                 \
00730                   \
00731     template <typename ExprT1>            \
00732     inline bool               \
00733     operator OP (const Expr<ExprT1>& expr1,       \
00734      const typename Expr<ExprT1>::value_type& b)    \
00735     {                 \
00736       return expr1.val() OP b;            \
00737     }                 \
00738   }                 \
00739 }
00740 
00741 FAD_RELOP_MACRO(==)
00742 FAD_RELOP_MACRO(!=)
00743 FAD_RELOP_MACRO(<)
00744 FAD_RELOP_MACRO(>)
00745 FAD_RELOP_MACRO(<=)
00746 FAD_RELOP_MACRO(>=)
00747 FAD_RELOP_MACRO(<<=)
00748 FAD_RELOP_MACRO(>>=)
00749 FAD_RELOP_MACRO(&)
00750 FAD_RELOP_MACRO(|)
00751 
00752 #undef FAD_RELOP_MACRO
00753 
00754 namespace Sacado {
00755 
00756   namespace Fad {
00757 
00758     template <typename ExprT>
00759     inline bool operator ! (const Expr<ExprT>& expr) 
00760     {
00761       return ! expr.val();
00762     }
00763 
00764   } // namespace Fad
00765 
00766 } // namespace Sacado
00767 
00768 //-------------------------- Boolean Operators -----------------------
00769 namespace Sacado {
00770 
00771   namespace Fad {
00772 
00773     template <typename ExprT>
00774     bool toBool(const Expr<ExprT>& x) {
00775       bool is_zero = (x.val() == 0.0);
00776       for (int i=0; i<x.size(); i++)
00777   is_zero = is_zero && (x.dx(i) == 0.0);
00778       return !is_zero;
00779     }
00780 
00781   } // namespace Fad
00782 
00783 } // namespace Sacado
00784 
00785 #define FAD_BOOL_MACRO(OP)            \
00786 namespace Sacado {              \
00787   namespace Fad {             \
00788     template <typename ExprT1, typename ExprT2>       \
00789     inline bool               \
00790     operator OP (const Expr<ExprT1>& expr1,       \
00791      const Expr<ExprT2>& expr2)       \
00792     {                 \
00793       return toBool(expr1) OP toBool(expr2);        \
00794     }                 \
00795                   \
00796     template <typename ExprT2>            \
00797     inline bool               \
00798     operator OP (const typename Expr<ExprT2>::value_type& a,    \
00799      const Expr<ExprT2>& expr2)       \
00800     {                 \
00801       return a OP toBool(expr2);          \
00802     }                 \
00803                   \
00804     template <typename ExprT1>            \
00805     inline bool               \
00806     operator OP (const Expr<ExprT1>& expr1,       \
00807      const typename Expr<ExprT1>::value_type& b)    \
00808     {                 \
00809       return toBool(expr1) OP b;          \
00810     }                 \
00811   }                 \
00812 }
00813 
00814 FAD_BOOL_MACRO(&&)
00815 FAD_BOOL_MACRO(||)
00816 
00817 #undef FAD_BOOL_MACRO
00818 
00819 //-------------------------- I/O Operators -----------------------
00820 
00821 namespace Sacado {
00822 
00823   namespace Fad {
00824 
00825     template <typename ExprT>
00826     std::ostream& operator << (std::ostream& os, const Expr<ExprT>& x) {
00827       os << x.val() << " [";
00828       
00829       for (int i=0; i< x.size(); i++) {
00830         os << " " << x.dx(i);
00831       }
00832 
00833       os << " ]";
00834       return os;
00835     }
00836 
00837   } // namespace Fad
00838 
00839 } // namespace Sacado
00840 
00841 
00842 #endif // SACADO_FAD_OPS_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines