Sacado_Fad_Ops.hpp

Go to the documentation of this file.
00001 // $Id: Sacado_Fad_Ops.hpp,v 1.13.2.3 2009/03/05 19:41:29 etphipp Exp $ 
00002 // $Source: /space/CVS/Trilinos/packages/sacado/src/Sacado_Fad_Ops.hpp,v $ 
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 <cmath>
00059 #include <algorithm>  // for std::min and std::max
00060 #include <ostream>  // for std::ostream
00061 
00062 #define FAD_UNARYOP_MACRO(OPNAME,OP,VALUE,DX,FASTACCESSDX)    \
00063 namespace Sacado {              \
00064   namespace Fad {             \
00065                   \
00066     template <typename ExprT>           \
00067     class OP {                \
00068     public:               \
00069                   \
00070       typedef typename ExprT::value_type value_type;      \
00071                   \
00072       OP() {}                           \
00073                   \
00074       static value_type computeValue(const ExprT& expr) {   \
00075   return VALUE;             \
00076       }                 \
00077                   \
00078       static value_type computeDx(int i, const ExprT& expr) {   \
00079   return DX;              \
00080       }                 \
00081                   \
00082       static value_type computeFastAccessDx(int i, const ExprT& expr) { \
00083   return FASTACCESSDX;            \
00084       }                 \
00085     };                  \
00086                   \
00087     template <typename T>           \
00088     inline Expr< UnaryExpr< Expr<T>, OP > >       \
00089     OPNAME (const Expr<T>& expr)          \
00090     {                 \
00091       typedef UnaryExpr< Expr<T>, OP > expr_t;        \
00092                         \
00093       return Expr<expr_t>(expr_t(expr));        \
00094     }                 \
00095   }                 \
00096 }
00097 
00098 FAD_UNARYOP_MACRO(operator+,
00099       UnaryPlusOp, 
00100       expr.val(),
00101       expr.dx(i),
00102       expr.fastAccessDx(i))
00103 FAD_UNARYOP_MACRO(operator-,
00104       UnaryMinusOp, 
00105       -expr.val(),
00106       -expr.dx(i),
00107       -expr.fastAccessDx(i))
00108 FAD_UNARYOP_MACRO(exp,
00109       ExpOp, 
00110       std::exp(expr.val()),
00111       std::exp(expr.val())*expr.dx(i),
00112       std::exp(expr.val())*expr.fastAccessDx(i))
00113 FAD_UNARYOP_MACRO(log,
00114       LogOp, 
00115       std::log(expr.val()),
00116       expr.dx(i)/expr.val(),
00117       expr.fastAccessDx(i)/expr.val())
00118 FAD_UNARYOP_MACRO(log10,
00119       Log10Op, 
00120       std::log10(expr.val()),
00121       expr.dx(i)/( std::log(value_type(10))*expr.val()),
00122       expr.fastAccessDx(i) / ( std::log(value_type(10))*expr.val()))
00123 FAD_UNARYOP_MACRO(sqrt,
00124       SqrtOp, 
00125       std::sqrt(expr.val()),
00126       expr.dx(i)/(value_type(2)* std::sqrt(expr.val())),
00127       expr.fastAccessDx(i)/(value_type(2)* std::sqrt(expr.val())))
00128 FAD_UNARYOP_MACRO(cos,
00129       CosOp, 
00130       std::cos(expr.val()),
00131       -expr.dx(i)* std::sin(expr.val()),
00132       -expr.fastAccessDx(i)* std::sin(expr.val()))
00133 FAD_UNARYOP_MACRO(sin,
00134       SinOp, 
00135       std::sin(expr.val()),
00136       expr.dx(i)* std::cos(expr.val()),
00137       expr.fastAccessDx(i)* std::cos(expr.val()))
00138 FAD_UNARYOP_MACRO(tan,
00139       TanOp, 
00140       std::tan(expr.val()),
00141       expr.dx(i)*
00142         (value_type(1)+ std::tan(expr.val())* std::tan(expr.val())),
00143       expr.fastAccessDx(i)*
00144         (value_type(1)+ std::tan(expr.val())* std::tan(expr.val())))
00145 FAD_UNARYOP_MACRO(acos,
00146       ACosOp, 
00147       std::acos(expr.val()),
00148       -expr.dx(i)/ std::sqrt(value_type(1)-expr.val()*expr.val()),
00149       -expr.fastAccessDx(i) /
00150         std::sqrt(value_type(1)-expr.val()*expr.val()))
00151 FAD_UNARYOP_MACRO(asin,
00152       ASinOp, 
00153       std::asin(expr.val()),
00154       expr.dx(i)/ std::sqrt(value_type(1)-expr.val()*expr.val()),
00155       expr.fastAccessDx(i) /
00156         std::sqrt(value_type(1)-expr.val()*expr.val()))
00157 FAD_UNARYOP_MACRO(atan,
00158       ATanOp, 
00159       std::atan(expr.val()),
00160       expr.dx(i)/(value_type(1)+expr.val()*expr.val()),
00161       expr.fastAccessDx(i)/(value_type(1)+expr.val()*expr.val()))
00162 FAD_UNARYOP_MACRO(cosh,
00163       CoshOp, 
00164       std::cosh(expr.val()),
00165       expr.dx(i)* std::sinh(expr.val()),
00166       expr.fastAccessDx(i)* std::sinh(expr.val()))
00167 FAD_UNARYOP_MACRO(sinh,
00168       SinhOp, 
00169       std::sinh(expr.val()),
00170       expr.dx(i)* std::cosh(expr.val()),
00171       expr.fastAccessDx(i)* std::cosh(expr.val()))
00172 FAD_UNARYOP_MACRO(tanh,
00173       TanhOp, 
00174       std::tanh(expr.val()),
00175       expr.dx(i)/( std::cosh(expr.val())* std::cosh(expr.val())),
00176       expr.fastAccessDx(i) / 
00177         ( std::cosh(expr.val())* std::cosh(expr.val())))
00178 FAD_UNARYOP_MACRO(acosh,
00179       ACoshOp, 
00180       acosh(expr.val()),
00181       expr.dx(i)/ std::sqrt((expr.val()-value_type(1)) * 
00182                (expr.val()+value_type(1))),
00183       expr.fastAccessDx(i)/ std::sqrt((expr.val()-value_type(1)) * 
00184              (expr.val()+value_type(1))))
00185 FAD_UNARYOP_MACRO(asinh,
00186       ASinhOp, 
00187       asinh(expr.val()),
00188       expr.dx(i)/ std::sqrt(value_type(1)+expr.val()*expr.val()),
00189       expr.fastAccessDx(i)/ std::sqrt(value_type(1)+
00190              expr.val()*expr.val()))
00191 FAD_UNARYOP_MACRO(atanh,
00192       ATanhOp, 
00193       atanh(expr.val()),
00194       expr.dx(i)/(value_type(1)-expr.val()*expr.val()),
00195       expr.fastAccessDx(i)/(value_type(1)-
00196              expr.val()*expr.val()))
00197 FAD_UNARYOP_MACRO(abs,
00198       AbsOp, 
00199       std::abs(expr.val()),
00200       expr.val() >= 0 ? value_type(+expr.dx(i)) : 
00201         value_type(-expr.dx(i)),
00202       expr.val() >= 0 ? value_type(+expr.fastAccessDx(i)) : 
00203         value_type(-expr.fastAccessDx(i)))
00204 FAD_UNARYOP_MACRO(fabs,
00205       FAbsOp, 
00206       std::fabs(expr.val()),
00207       expr.val() >= 0 ? value_type(+expr.dx(i)) : 
00208         value_type(-expr.dx(i)),
00209       expr.val() >= 0 ? value_type(+expr.fastAccessDx(i)) : 
00210         value_type(-expr.fastAccessDx(i)))
00211 
00212 #undef FAD_UNARYOP_MACRO
00213 
00214 #define FAD_BINARYOP_MACRO(OPNAME,OP,VALUE,DX,FASTACCESSDX,CONST_DX_1,CONST_DX_2,CONST_FASTACCESSDX_1,CONST_FASTACCESSDX_2) \
00215 namespace Sacado {              \
00216   namespace Fad {             \
00217                   \
00218     template <typename ExprT1, typename ExprT2>       \
00219     class OP {                \
00220                   \
00221     public:               \
00222                   \
00223       typedef typename ExprT1::value_type value_type_1;     \
00224       typedef typename ExprT2::value_type value_type_2;     \
00225       typedef typename Sacado::Promote<value_type_1,      \
00226                value_type_2>::type value_type;  \
00227                   \
00228       OP() {}                                             \
00229                   \
00230       static value_type             \
00231       computeValue(const ExprT1& expr1, const ExprT2& expr2) {          \
00232   return VALUE;             \
00233       }                 \
00234                   \
00235       static value_type             \
00236       computeDx(int i, const ExprT1& expr1,       \
00237     const ExprT2& expr2) {                \
00238   return DX;              \
00239       }                 \
00240                   \
00241       static value_type             \
00242       computeFastAccessDx(int i, const ExprT1& expr1,     \
00243         const ExprT2& expr2) {      \
00244   return FASTACCESSDX;            \
00245       }                 \
00246     };                  \
00247                   \
00248     template <typename ExprT1>            \
00249     class OP<ExprT1, ConstExpr<typename ExprT1::value_type> > {   \
00250                   \
00251     public:               \
00252                   \
00253       typedef typename ExprT1::value_type value_type;     \
00254       typedef ConstExpr<typename ExprT1::value_type> ExprT2;    \
00255                   \
00256       OP() {}                                             \
00257                   \
00258       static value_type             \
00259       computeValue(const ExprT1& expr1, const ExprT2& expr2) {          \
00260   return VALUE;             \
00261       }                 \
00262                   \
00263       static value_type             \
00264       computeDx(int i, const ExprT1& expr1,       \
00265     const ExprT2& expr2) {                \
00266   return CONST_DX_2;            \
00267       }                 \
00268                   \
00269       static value_type             \
00270       computeFastAccessDx(int i, const ExprT1& expr1,     \
00271         const ExprT2& expr2) {      \
00272   return CONST_FASTACCESSDX_2;          \
00273       }                 \
00274     };                  \
00275                   \
00276     template <typename ExprT2>            \
00277     class OP<ConstExpr<typename ExprT2::value_type>, ExprT2 > {   \
00278                   \
00279     public:               \
00280                   \
00281       typedef typename ExprT2::value_type value_type;     \
00282       typedef ConstExpr<typename ExprT2::value_type> ExprT1;    \
00283                   \
00284       OP() {}                                             \
00285                   \
00286       static value_type             \
00287       computeValue(const ExprT1& expr1, const ExprT2& expr2){         \
00288   return VALUE;             \
00289       }                 \
00290                   \
00291       static value_type             \
00292       computeDx(int i, const ExprT1& expr1,       \
00293     const ExprT2& expr2) {                \
00294   return CONST_DX_1;            \
00295       }                 \
00296                   \
00297       static value_type             \
00298       computeFastAccessDx(int i, const ExprT1& expr1,     \
00299         const ExprT2& expr2) {      \
00300   return CONST_FASTACCESSDX_1;          \
00301       }                 \
00302     };                  \
00303                   \
00304     template <typename T1, typename T2>         \
00305     inline Expr< BinaryExpr< Expr<T1>, Expr<T2>, OP > >     \
00306     OPNAME (const Expr<T1>& expr1, const Expr<T2>& expr2)   \
00307     {                 \
00308       typedef BinaryExpr< Expr<T1>, Expr<T2>, OP > expr_t;    \
00309                       \
00310       return Expr<expr_t>(expr_t(expr1, expr2));      \
00311     }                 \
00312                   \
00313     template <typename T>           \
00314     inline Expr< BinaryExpr< Expr<T>, Expr<T>, OP > >     \
00315     OPNAME (const Expr<T>& expr1, const Expr<T>& expr2)     \
00316     {                 \
00317       typedef BinaryExpr< Expr<T>, Expr<T>, OP > expr_t;    \
00318                       \
00319       return Expr<expr_t>(expr_t(expr1, expr2));      \
00320     }                 \
00321                   \
00322     template <typename T>           \
00323     inline Expr< BinaryExpr< ConstExpr<typename Expr<T>::value_type>, \
00324            Expr<T>, OP > >        \
00325     OPNAME (const typename Expr<T>::value_type& c,      \
00326       const Expr<T>& expr)          \
00327     {                 \
00328       typedef ConstExpr<typename Expr<T>::value_type> ConstT;   \
00329       typedef BinaryExpr< ConstT, Expr<T>, OP > expr_t;     \
00330                   \
00331       return Expr<expr_t>(expr_t(ConstT(c), expr));     \
00332     }                 \
00333                   \
00334     template <typename T>           \
00335     inline Expr< BinaryExpr< Expr<T>,         \
00336            ConstExpr<typename Expr<T>::value_type>, \
00337            OP > >         \
00338     OPNAME (const Expr<T>& expr,          \
00339       const typename Expr<T>::value_type& c)      \
00340     {                 \
00341       typedef ConstExpr<typename Expr<T>::value_type> ConstT;   \
00342       typedef BinaryExpr< Expr<T>, ConstT, OP > expr_t;     \
00343                   \
00344       return Expr<expr_t>(expr_t(expr, ConstT(c)));     \
00345     }                 \
00346   }                 \
00347 }
00348 
00349 
00350 FAD_BINARYOP_MACRO(operator+,
00351        AdditionOp, 
00352        expr1.val() + expr2.val(),
00353        expr1.dx(i) + expr2.dx(i),
00354        expr1.fastAccessDx(i) + expr2.fastAccessDx(i),
00355        expr2.dx(i),
00356        expr1.dx(i),
00357        expr2.fastAccessDx(i),
00358        expr1.fastAccessDx(i))
00359 FAD_BINARYOP_MACRO(operator-,
00360        SubtractionOp, 
00361        expr1.val() - expr2.val(),
00362        expr1.dx(i) - expr2.dx(i),
00363        expr1.fastAccessDx(i) - expr2.fastAccessDx(i),
00364        -expr2.dx(i),
00365        expr1.dx(i),
00366        -expr2.fastAccessDx(i),
00367        expr1.fastAccessDx(i))
00368 FAD_BINARYOP_MACRO(operator*,
00369        MultiplicationOp, 
00370        expr1.val() * expr2.val(),
00371        expr1.val()*expr2.dx(i) + expr1.dx(i)*expr2.val(),
00372        expr1.val()*expr2.fastAccessDx(i) + 
00373          expr1.fastAccessDx(i)*expr2.val(),
00374        expr1.val()*expr2.dx(i),
00375        expr1.dx(i)*expr2.val(),
00376        expr1.val()*expr2.fastAccessDx(i),
00377        expr1.fastAccessDx(i)*expr2.val())
00378 FAD_BINARYOP_MACRO(operator/,
00379        DivisionOp, 
00380        expr1.val() / expr2.val(),
00381        (expr1.dx(i)*expr2.val() - expr2.dx(i)*expr1.val()) /
00382          (expr2.val()*expr2.val()),
00383        (expr1.fastAccessDx(i)*expr2.val() - 
00384           expr2.fastAccessDx(i)*expr1.val()) /
00385           (expr2.val()*expr2.val()),
00386        -expr2.dx(i)*expr1.val() / (expr2.val()*expr2.val()),
00387        expr1.dx(i)/expr2.val(),
00388        -expr2.fastAccessDx(i)*expr1.val() / 
00389          (expr2.val()*expr2.val()),
00390        expr1.fastAccessDx(i)/expr2.val())
00391 FAD_BINARYOP_MACRO(atan2,
00392        Atan2Op,
00393        std::atan2(expr1.val(), expr2.val()),
00394        (expr2.val()*expr1.dx(i) - expr1.val()*expr2.dx(i))/
00395       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00396        (expr2.val()*expr1.fastAccessDx(i) - expr1.val()*expr2.fastAccessDx(i))/
00397       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00398        (-expr1.val()*expr2.dx(i))/
00399       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00400        (expr2.val()*expr1.dx(i))/
00401       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00402        (-expr1.val()*expr2.fastAccessDx(i))/
00403       (expr1.val()*expr1.val() + expr2.val()*expr2.val()),
00404        (expr2.val()*expr1.fastAccessDx(i))/
00405       (expr1.val()*expr1.val() + expr2.val()*expr2.val()))
00406 FAD_BINARYOP_MACRO(pow,
00407        PowerOp,
00408        std::pow(expr1.val(), expr2.val()),
00409        (expr2.dx(i)* std::log(expr1.val())+expr2.val()*expr1.dx(i)/
00410         expr1.val())* std::pow(expr1.val(),expr2.val()),
00411        (expr2.fastAccessDx(i)* std::log(expr1.val())+
00412         expr2.val()*expr1.fastAccessDx(i)/
00413         expr1.val())* std::pow(expr1.val(),expr2.val()),
00414        expr2.dx(i)* std::log(expr1.val()) *
00415          std::pow(expr1.val(),expr2.val()),
00416        expr2.val()*expr1.dx(i)/
00417        expr1.val()* std::pow(expr1.val(),expr2.val()),
00418        expr2.fastAccessDx(i)* std::log(expr1.val()) *
00419          std::pow(expr1.val(),expr2.val()),
00420        expr2.val()*expr1.fastAccessDx(i)/
00421          expr1.val()* std::pow(expr1.val(),expr2.val()))
00422 FAD_BINARYOP_MACRO(max,
00423                    MaxOp,
00424                    std::max(expr1.val(), expr2.val()),
00425                    expr1.val() >= expr2.val() ? expr1.dx(i) : expr2.dx(i),
00426                    expr1.val() >= expr2.val() ? expr1.fastAccessDx(i) : 
00427                                                 expr2.fastAccessDx(i),
00428                    expr1.val() >= expr2.val() ? value_type(0) : expr2.dx(i),
00429                    expr1.val() >= expr2.val() ? expr1.dx(i) : value_type(0),
00430                    expr1.val() >= expr2.val() ? value_type(0) : 
00431                                                 expr2.fastAccessDx(i),
00432                    expr1.val() >= expr2.val() ? expr1.fastAccessDx(i) : 
00433                                                 value_type(0))
00434 FAD_BINARYOP_MACRO(min,
00435                    MinOp,
00436                    std::min(expr1.val(), expr2.val()),
00437                    expr1.val() <= expr2.val() ? expr1.dx(i) : expr2.dx(i),
00438                    expr1.val() <= expr2.val() ? expr1.fastAccessDx(i) : 
00439                                                 expr2.fastAccessDx(i),
00440                    expr1.val() <= expr2.val() ? value_type(0) : expr2.dx(i),
00441                    expr1.val() <= expr2.val() ? expr1.dx(i) : value_type(0),
00442                    expr1.val() <= expr2.val() ? value_type(0) : 
00443                                                 expr2.fastAccessDx(i),
00444                    expr1.val() <= expr2.val() ? expr1.fastAccessDx(i) : 
00445                                                 value_type(0))
00446 
00447 #undef FAD_BINARYOP_MACRO
00448 
00449   // The general definition of max/min works for Fad variables too, except
00450   // we need to add a case when the argument types are different.  This 
00451   // can't conflict with the general definition, so we need to use
00452   // Substitution Failure Is Not An Error
00453 #include "Sacado_mpl_disable_if.hpp"
00454 #include "Sacado_mpl_is_same.hpp"
00455 
00456 #define FAD_SFINAE_BINARYOP_MACRO(OPNAME,OP,VALUE,DX,FASTACCESSDX,CONST_DX_1,CONST_DX_2,CONST_FASTACCESSDX_1,CONST_FASTACCESSDX_2) \
00457 namespace Sacado {              \
00458   namespace Fad {             \
00459                   \
00460     template <typename ExprT1, typename ExprT2>       \
00461     class OP {                \
00462                   \
00463     public:               \
00464                   \
00465       typedef typename ExprT1::value_type value_type_1;     \
00466       typedef typename ExprT2::value_type value_type_2;     \
00467       typedef typename Sacado::Promote<value_type_1,      \
00468                value_type_2>::type value_type;  \
00469                   \
00470       OP() {}                                             \
00471                   \
00472       static value_type             \
00473       computeValue(const ExprT1& expr1, const ExprT2& expr2) {          \
00474   return VALUE;             \
00475       }                 \
00476                   \
00477       static value_type             \
00478       computeDx(int i, const ExprT1& expr1,       \
00479     const ExprT2& expr2) {                \
00480   return DX;              \
00481       }                 \
00482                   \
00483       static value_type             \
00484       computeFastAccessDx(int i, const ExprT1& expr1,     \
00485         const ExprT2& expr2) {      \
00486   return FASTACCESSDX;            \
00487       }                 \
00488     };                  \
00489                   \
00490     template <typename ExprT1>            \
00491     class OP<ExprT1, ConstExpr<typename ExprT1::value_type> > {   \
00492                   \
00493     public:               \
00494                   \
00495       typedef typename ExprT1::value_type value_type;     \
00496       typedef ConstExpr<typename ExprT1::value_type> ExprT2;    \
00497                   \
00498       OP() {}                                             \
00499                   \
00500       static value_type             \
00501       computeValue(const ExprT1& expr1, const ExprT2& expr2) {          \
00502   return VALUE;             \
00503       }                 \
00504                   \
00505       static value_type             \
00506       computeDx(int i, const ExprT1& expr1,       \
00507     const ExprT2& expr2) {                \
00508   return CONST_DX_2;            \
00509       }                 \
00510                   \
00511       static value_type             \
00512       computeFastAccessDx(int i, const ExprT1& expr1,     \
00513         const ExprT2& expr2) {      \
00514   return CONST_FASTACCESSDX_2;          \
00515       }                 \
00516     };                  \
00517                   \
00518     template <typename ExprT2>            \
00519     class OP<ConstExpr<typename ExprT2::value_type>, ExprT2 > {   \
00520                   \
00521     public:               \
00522                   \
00523       typedef typename ExprT2::value_type value_type;     \
00524       typedef ConstExpr<typename ExprT2::value_type> ExprT1;    \
00525                   \
00526       OP() {}                                             \
00527                   \
00528       static value_type             \
00529       computeValue(const ExprT1& expr1, const ExprT2& expr2){         \
00530   return VALUE;             \
00531       }                 \
00532                   \
00533       static value_type             \
00534       computeDx(int i, const ExprT1& expr1,       \
00535     const ExprT2& expr2) {                \
00536   return CONST_DX_1;            \
00537       }                 \
00538                   \
00539       static value_type             \
00540       computeFastAccessDx(int i, const ExprT1& expr1,     \
00541         const ExprT2& expr2) {      \
00542   return CONST_FASTACCESSDX_1;          \
00543       }                 \
00544     };                  \
00545                   \
00546     template <typename T1, typename T2>         \
00547     inline                                                              \
00548     typename                                                            \
00549     mpl::disable_if< mpl::is_same<T1,T2>,                               \
00550                      Expr<BinaryExpr<Expr<T1>, Expr<T2>, OP> > >::type  \
00551     OPNAME (const Expr<T1>& expr1, const Expr<T2>& expr2)   \
00552     {                 \
00553       typedef BinaryExpr< Expr<T1>, Expr<T2>, OP > expr_t;    \
00554                       \
00555       return Expr<expr_t>(expr_t(expr1, expr2));      \
00556     }                 \
00557                   \
00558     template <typename T>           \
00559     inline Expr< BinaryExpr< ConstExpr<typename Expr<T>::value_type>, \
00560            Expr<T>, OP > >        \
00561     OPNAME (const typename Expr<T>::value_type& c,      \
00562       const Expr<T>& expr)          \
00563     {                 \
00564       typedef ConstExpr<typename Expr<T>::value_type> ConstT;   \
00565       typedef BinaryExpr< ConstT, Expr<T>, OP > expr_t;     \
00566                   \
00567       return Expr<expr_t>(expr_t(ConstT(c), expr));     \
00568     }                 \
00569                   \
00570     template <typename T>           \
00571     inline Expr< BinaryExpr< Expr<T>,         \
00572            ConstExpr<typename Expr<T>::value_type>, \
00573            OP > >         \
00574     OPNAME (const Expr<T>& expr,          \
00575       const typename Expr<T>::value_type& c)      \
00576     {                 \
00577       typedef ConstExpr<typename Expr<T>::value_type> ConstT;   \
00578       typedef BinaryExpr< Expr<T>, ConstT, OP > expr_t;     \
00579                   \
00580       return Expr<expr_t>(expr_t(expr, ConstT(c)));     \
00581     }                 \
00582   }                 \
00583 }
00584 
00585 // FAD_SFINAE_BINARYOP_MACRO(max,
00586 //                    MaxOp,
00587 //                    std::max(expr1.val(), expr2.val()),
00588 //                    expr1.val() >= expr2.val() ? expr1.dx(i) : expr2.dx(i),
00589 //                    expr1.val() >= expr2.val() ? expr1.fastAccessDx(i) : 
00590 //                                                 expr2.fastAccessDx(i),
00591 //                    expr1.val() >= expr2.val() ? value_type(0) : expr2.dx(i),
00592 //                    expr1.val() >= expr2.val() ? expr1.dx(i) : value_type(0),
00593 //                    expr1.val() >= expr2.val() ? value_type(0) : 
00594 //                                                 expr2.fastAccessDx(i),
00595 //                    expr1.val() >= expr2.val() ? expr1.fastAccessDx(i) : 
00596 //                                                 value_type(0))
00597 // FAD_SFINAE_BINARYOP_MACRO(min,
00598 //                    MinOp,
00599 //                    std::min(expr1.val(), expr2.val()),
00600 //                    expr1.val() <= expr2.val() ? expr1.dx(i) : expr2.dx(i),
00601 //                    expr1.val() <= expr2.val() ? expr1.fastAccessDx(i) : 
00602 //                                                 expr2.fastAccessDx(i),
00603 //                    expr1.val() <= expr2.val() ? value_type(0) : expr2.dx(i),
00604 //                    expr1.val() <= expr2.val() ? expr1.dx(i) : value_type(0),
00605 //                    expr1.val() <= expr2.val() ? value_type(0) : 
00606 //                                                 expr2.fastAccessDx(i),
00607 //                    expr1.val() <= expr2.val() ? expr1.fastAccessDx(i) : 
00608 //                                                 value_type(0))
00609 
00610 #undef FAD_SFINAE_BINARYOP_MACRO
00611 
00612 //-------------------------- Relational Operators -----------------------
00613 
00614 #define FAD_RELOP_MACRO(OP)           \
00615 namespace Sacado {              \
00616   namespace Fad {             \
00617     template <typename ExprT1, typename ExprT2>       \
00618     inline bool               \
00619     operator OP (const Expr<ExprT1>& expr1,       \
00620      const Expr<ExprT2>& expr2)       \
00621     {                 \
00622       return expr1.val() OP expr2.val();        \
00623     }                 \
00624                   \
00625     template <typename ExprT2>            \
00626     inline bool               \
00627     operator OP (const typename Expr<ExprT2>::value_type& a,    \
00628      const Expr<ExprT2>& expr2)       \
00629     {                 \
00630       return a OP expr2.val();            \
00631     }                 \
00632                   \
00633     template <typename ExprT1>            \
00634     inline bool               \
00635     operator OP (const Expr<ExprT1>& expr1,       \
00636      const typename Expr<ExprT1>::value_type& b)    \
00637     {                 \
00638       return expr1.val() OP b;            \
00639     }                 \
00640   }                 \
00641 }
00642 
00643 FAD_RELOP_MACRO(==)
00644 FAD_RELOP_MACRO(!=)
00645 FAD_RELOP_MACRO(<)
00646 FAD_RELOP_MACRO(>)
00647 FAD_RELOP_MACRO(<=)
00648 FAD_RELOP_MACRO(>=)
00649 FAD_RELOP_MACRO(<<=)
00650 FAD_RELOP_MACRO(>>=)
00651 FAD_RELOP_MACRO(&)
00652 FAD_RELOP_MACRO(|)
00653 
00654 #undef FAD_RELOP_MACRO
00655 
00656 namespace Sacado {
00657 
00658   namespace Fad {
00659 
00660     template <typename ExprT>
00661     inline bool operator ! (const Expr<ExprT>& expr) 
00662     {
00663       return ! expr.val();
00664     }
00665 
00666   } // namespace Fad
00667 
00668 } // namespace Sacado
00669 
00670 //-------------------------- I/O Operators -----------------------
00671 
00672 namespace Sacado {
00673 
00674   namespace Fad {
00675 
00676     template <typename ExprT>
00677     std::ostream& operator << (std::ostream& os, const Expr<ExprT>& x) {
00678       os << x.val() << "  [";
00679       
00680       for (int i=0; i< x.size(); i++) {
00681   os << " " << x.dx(i);
00682       }
00683 
00684       os << " ]\n";
00685       return os;
00686     }
00687 
00688   } // namespace Fad
00689 
00690 } // namespace Sacado
00691 
00692 
00693 #endif // SACADO_FAD_OPS_HPP

Generated on Tue Oct 20 12:55:05 2009 for Sacado Package Browser (Single Doxygen Collection) by doxygen 1.4.7