|
Sacado Package Browser (Single Doxygen Collection) Version of the Day
|
|
Sacado Package Browser (Single Doxygen Collection) Version of the Day
|
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 #undef FAD_BINARYOP_MACRO 00515 00516 namespace Sacado { 00517 namespace Fad { 00518 00519 template <typename ExprT1, typename ExprT2> 00520 class MultiplicationOp {}; 00521 00522 template <typename T1, typename T2> 00523 class Expr< MultiplicationOp< Expr<T1>, Expr<T2> > > { 00524 00525 public: 00526 00527 typedef Expr<T1> ExprT1; 00528 typedef Expr<T2> ExprT2; 00529 typedef typename ExprT1::value_type value_type_1; 00530 typedef typename ExprT2::value_type value_type_2; 00531 typedef typename Sacado::Promote<value_type_1, 00532 value_type_2>::type value_type; 00533 00534 Expr(const ExprT1& expr1_, const ExprT2& expr2_) : 00535 expr1(expr1_), expr2(expr2_) {} 00536 00537 int size() const { 00538 int sz1 = expr1.size(), sz2 = expr2.size(); 00539 return sz1 > sz2 ? sz1 : sz2; 00540 } 00541 00542 bool hasFastAccess() const { 00543 return expr1.hasFastAccess() && expr2.hasFastAccess(); 00544 } 00545 00546 bool isPassive() const { 00547 return expr1.isPassive() && expr2.isPassive(); 00548 } 00549 00550 bool updateValue() const { 00551 return expr1.updateValue() && expr2.updateValue(); 00552 } 00553 00554 const value_type val() const { 00555 return expr1.val()*expr2.val(); 00556 } 00557 00558 const value_type dx(int i) const { 00559 if (expr1.size() > 0 && expr2.size() > 0) 00560 return expr1.val()*expr2.dx(i) + expr1.dx(i)*expr2.val(); 00561 else if (expr1.size() > 0) 00562 return expr1.dx(i)*expr2.val(); 00563 else 00564 return expr1.val()*expr2.dx(i); 00565 } 00566 00567 const value_type fastAccessDx(int i) const { 00568 return expr1.val()*expr2.fastAccessDx(i) + 00569 expr1.fastAccessDx(i)*expr2.val(); 00570 } 00571 00572 protected: 00573 00574 const ExprT1& expr1; 00575 const ExprT2& expr2; 00576 00577 }; 00578 00579 template <typename T1> 00580 class Expr< MultiplicationOp< Expr<T1>, typename Expr<T1>::value_type> > { 00581 00582 public: 00583 00584 typedef Expr<T1> ExprT1; 00585 typedef typename ExprT1::value_type value_type; 00586 typedef typename ExprT1::value_type ConstT; 00587 00588 Expr(const ExprT1& expr1_, const ConstT& c_) : 00589 expr1(expr1_), c(c_) {} 00590 00591 int size() const { 00592 return expr1.size(); 00593 } 00594 00595 bool hasFastAccess() const { 00596 return expr1.hasFastAccess(); 00597 } 00598 00599 bool isPassive() const { 00600 return expr1.isPassive(); 00601 } 00602 00603 bool updateValue() const { return expr1.updateValue(); } 00604 00605 const value_type val() const { 00606 return expr1.val()*c; 00607 } 00608 00609 const value_type dx(int i) const { 00610 return expr1.dx(i)*c; 00611 } 00612 00613 const value_type fastAccessDx(int i) const { 00614 return expr1.fastAccessDx(i)*c; 00615 } 00616 00617 protected: 00618 00619 const ExprT1& expr1; 00620 const ConstT& c; 00621 }; 00622 00623 template <typename T2> 00624 class Expr< MultiplicationOp< typename Expr<T2>::value_type, Expr<T2> > > { 00625 00626 public: 00627 00628 typedef Expr<T2> ExprT2; 00629 typedef typename ExprT2::value_type value_type; 00630 typedef typename ExprT2::value_type ConstT; 00631 00632 Expr(const ConstT& c_, const ExprT2& expr2_) : 00633 c(c_), expr2(expr2_) {} 00634 00635 int size() const { 00636 return expr2.size(); 00637 } 00638 00639 bool hasFastAccess() const { 00640 return expr2.hasFastAccess(); 00641 } 00642 00643 bool isPassive() const { 00644 return expr2.isPassive(); 00645 } 00646 00647 bool updateValue() const { return expr2.updateValue(); } 00648 00649 const value_type val() const { 00650 return c*expr2.val(); 00651 } 00652 00653 const value_type dx(int i) const { 00654 return c*expr2.dx(i); 00655 } 00656 00657 const value_type fastAccessDx(int i) const { 00658 return c*expr2.fastAccessDx(i); 00659 } 00660 00661 protected: 00662 00663 const ConstT& c; 00664 const ExprT2& expr2; 00665 }; 00666 00667 template <typename T1, typename T2> 00668 inline Expr< MultiplicationOp< Expr<T1>, Expr<T2> > > 00669 operator* (const Expr<T1>& expr1, const Expr<T2>& expr2) 00670 { 00671 typedef MultiplicationOp< Expr<T1>, Expr<T2> > expr_t; 00672 00673 return Expr<expr_t>(expr1, expr2); 00674 } 00675 00676 template <typename T> 00677 inline Expr< MultiplicationOp< Expr<T>, Expr<T> > > 00678 operator* (const Expr<T>& expr1, const Expr<T>& expr2) 00679 { 00680 typedef MultiplicationOp< Expr<T>, Expr<T> > expr_t; 00681 00682 return Expr<expr_t>(expr1, expr2); 00683 } 00684 00685 template <typename T> 00686 inline Expr< MultiplicationOp< typename Expr<T>::value_type, Expr<T> > > 00687 operator* (const typename Expr<T>::value_type& c, 00688 const Expr<T>& expr) 00689 { 00690 typedef typename Expr<T>::value_type ConstT; 00691 typedef MultiplicationOp< ConstT, Expr<T> > expr_t; 00692 00693 return Expr<expr_t>(c, expr); 00694 } 00695 00696 template <typename T> 00697 inline Expr< MultiplicationOp< Expr<T>, typename Expr<T>::value_type > > 00698 operator* (const Expr<T>& expr, 00699 const typename Expr<T>::value_type& c) 00700 { 00701 typedef typename Expr<T>::value_type ConstT; 00702 typedef MultiplicationOp< Expr<T>, ConstT > expr_t; 00703 00704 return Expr<expr_t>(expr, c); 00705 } 00706 } 00707 } 00708 00709 //-------------------------- Relational Operators ----------------------- 00710 00711 #define FAD_RELOP_MACRO(OP) \ 00712 namespace Sacado { \ 00713 namespace Fad { \ 00714 template <typename ExprT1, typename ExprT2> \ 00715 inline bool \ 00716 operator OP (const Expr<ExprT1>& expr1, \ 00717 const Expr<ExprT2>& expr2) \ 00718 { \ 00719 return expr1.val() OP expr2.val(); \ 00720 } \ 00721 \ 00722 template <typename ExprT2> \ 00723 inline bool \ 00724 operator OP (const typename Expr<ExprT2>::value_type& a, \ 00725 const Expr<ExprT2>& expr2) \ 00726 { \ 00727 return a OP expr2.val(); \ 00728 } \ 00729 \ 00730 template <typename ExprT1> \ 00731 inline bool \ 00732 operator OP (const Expr<ExprT1>& expr1, \ 00733 const typename Expr<ExprT1>::value_type& b) \ 00734 { \ 00735 return expr1.val() OP b; \ 00736 } \ 00737 } \ 00738 } 00739 00740 FAD_RELOP_MACRO(==) 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 00751 #undef FAD_RELOP_MACRO 00752 00753 namespace Sacado { 00754 00755 namespace Fad { 00756 00757 template <typename ExprT> 00758 inline bool operator ! (const Expr<ExprT>& expr) 00759 { 00760 return ! expr.val(); 00761 } 00762 00763 } // namespace Fad 00764 00765 } // namespace Sacado 00766 00767 //-------------------------- Boolean Operators ----------------------- 00768 namespace Sacado { 00769 00770 namespace Fad { 00771 00772 template <typename ExprT> 00773 bool toBool(const Expr<ExprT>& x) { 00774 bool is_zero = (x.val() == 0.0); 00775 for (int i=0; i<x.size(); i++) 00776 is_zero = is_zero && (x.dx(i) == 0.0); 00777 return !is_zero; 00778 } 00779 00780 } // namespace Fad 00781 00782 } // namespace Sacado 00783 00784 #define FAD_BOOL_MACRO(OP) \ 00785 namespace Sacado { \ 00786 namespace Fad { \ 00787 template <typename ExprT1, typename ExprT2> \ 00788 inline bool \ 00789 operator OP (const Expr<ExprT1>& expr1, \ 00790 const Expr<ExprT2>& expr2) \ 00791 { \ 00792 return toBool(expr1) OP toBool(expr2); \ 00793 } \ 00794 \ 00795 template <typename ExprT2> \ 00796 inline bool \ 00797 operator OP (const typename Expr<ExprT2>::value_type& a, \ 00798 const Expr<ExprT2>& expr2) \ 00799 { \ 00800 return a OP toBool(expr2); \ 00801 } \ 00802 \ 00803 template <typename ExprT1> \ 00804 inline bool \ 00805 operator OP (const Expr<ExprT1>& expr1, \ 00806 const typename Expr<ExprT1>::value_type& b) \ 00807 { \ 00808 return toBool(expr1) OP b; \ 00809 } \ 00810 } \ 00811 } 00812 00813 FAD_BOOL_MACRO(&&) 00814 FAD_BOOL_MACRO(||) 00815 00816 #undef FAD_BOOL_MACRO 00817 00818 //-------------------------- I/O Operators ----------------------- 00819 00820 namespace Sacado { 00821 00822 namespace Fad { 00823 00824 template <typename ExprT> 00825 std::ostream& operator << (std::ostream& os, const Expr<ExprT>& x) { 00826 os << x.val() << " ["; 00827 00828 for (int i=0; i< x.size(); i++) { 00829 os << " " << x.dx(i); 00830 } 00831 00832 os << " ]"; 00833 return os; 00834 } 00835 00836 } // namespace Fad 00837 00838 } // namespace Sacado 00839 00840 00841 #endif // SACADO_FAD_OPS_HPP
1.7.4
