Sacado Package Browser (Single Doxygen Collection) Version of the Day
Sacado_rad2.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 //
00004 //                           Sacado Package
00005 //                 Copyright (2007) Sandia Corporation
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
00025 // (etphipp@sandia.gov).
00026 //
00027 // ***********************************************************************
00028 // @HEADER
00029 
00030 // Extension of the RAD package (Reverse Automatic Differentiation) --
00031 // a package specialized for function and gradient evaluations -- to
00032 // Hessian-vector products.
00033 // This variant is for Hessian-vector products of "double" variables, and
00034 // Sacado::Rad2d::ADvar should be equivalent to Sacado::Rad2::ADvar<double>,
00035 // but this nontemplated code may easier to understand.  It relies on ops
00036 // implemented in Sacado_radops2.cpp.
00037 // Written in 2007 by David M. Gay at Sandia National Labs, Albuquerque, NM.
00038 
00039 #ifndef SACADO_RAD2_H
00040 #define SACADO_RAD2_H
00041 
00042 #include <stddef.h>
00043 #include <cmath>
00044 #include <math.h>
00045 
00046 #include "Sacado_ConfigDefs.h"
00047 
00048 #if defined(RAD_DEBUG_BLOCKKEEP) && !defined(HAVE_SACADO_UNINIT)
00049 #undef RAD_DEBUG_BLOCKKEEP
00050 #endif
00051 
00052 #ifndef SACADO_NO_NAMESPACE
00053 namespace Sacado {
00054 namespace Rad2d { // "2" for 2nd derivatives, "d" for "double"
00055 #endif
00056 
00057 // -DNO_USING_STDCC is needed, e.g., with Sun CC 5.7
00058 #ifndef RAD_NO_USING_STDCC
00059   // Bring math functions into scope
00060   using std::exp;
00061   using std::log;
00062   using std::log10;
00063   using std::sqrt;
00064   using std::cos;
00065   using std::sin;
00066   using std::tan;
00067   using std::acos;
00068   using std::asin;
00069   using std::atan;
00070   using std::cosh;
00071   using std::sinh;
00072   using std::tanh;
00073   using std::abs;
00074   using std::fabs;
00075   using std::atan2;
00076   using std::pow;
00077 #endif //NO_USING_STDCC
00078 
00079  class ADvar;
00080  class ADvari;
00081  class ADvar1;
00082  class ADvar2;
00083  class ConstADvar;
00084  class Derp;
00085  class IndepADvar;
00086 
00087  struct
00088 ADmemblock {  // We get memory in ADmemblock chunks and never give it back,
00089     // but reuse it once computations start anew after call(s) on
00090     // ADcontext::Gradcomp() or ADcontext::Weighted_Gradcomp().
00091   ADmemblock *next;
00092   double memblk[2000];
00093   };
00094 
00095  struct
00096 ADvari_block {
00097   enum { Gulp = 1021 };
00098   ADvari_block *next, *prev;
00099   ADvari **limit;
00100   ADvari *pADvari[Gulp];
00101   };
00102 
00103  class
00104 ADcontext { // A singleton class: one instance in radops.c
00105   ADmemblock *Busy, *Free;
00106   char *Mbase;
00107   size_t Mleft;
00108   ADvari **Ailimit, **Ainext;
00109   ADvari_block *Aibusy, *Aifree;
00110   ADmemblock First;
00111   ADvari_block AiFirst;
00112   void *new_ADmemblock(size_t);
00113   void new_ADvari_block();
00114  public:
00115   ADcontext();
00116   void *Memalloc(size_t len);
00117   static void Gradcomp(int);
00118   static inline void Gradcomp() { Gradcomp(1); }
00119   static void Hvprod(int, ADvar**, double*, double*);
00120   static void Weighted_Gradcomp(int, ADvar**, double*);
00121   inline void ADvari_record(ADvari *x) {
00122     if (Ainext >= Ailimit)
00123       new_ADvari_block();
00124     *Ainext++ = x;
00125     }
00126   };
00127 
00128 inline void *ADcontext::Memalloc(size_t len) {
00129     if (Mleft >= len)
00130       return Mbase + (Mleft -= len);
00131     return new_ADmemblock(len);
00132     }
00133 
00134  class
00135 CADcontext: public ADcontext {
00136  protected:
00137   bool fpval_implies_const;
00138  public:
00139   friend class ADvar;
00140   CADcontext(): ADcontext() { fpval_implies_const = false; }
00141   static const double One, negOne;
00142   };
00143 
00144  class
00145 Derp {    // one derivative-propagation operation
00146  public:
00147   friend class ADvarn;
00148   static Derp *LastDerp;
00149   Derp *next;
00150   const double *a;
00151   const ADvari *b;
00152   mutable ADvari *c;
00153   void *operator new(size_t);
00154   void operator delete(void*) {} /*Should never be called.*/
00155   inline Derp(){};
00156   Derp(const ADvari *);
00157   Derp(const double *, const ADvari *);
00158   Derp(const double *, const ADvari *, const ADvari *);
00159   /* c->aval += a * b->aval; */
00160   };
00161 
00162 inline Derp::Derp(const ADvari *c1): c((ADvari*)c1) {
00163     next = LastDerp;
00164     LastDerp = this;
00165     }
00166 
00167 inline Derp::Derp(const double *a1, const ADvari *c1):
00168     a(a1), c((ADvari*)c1) {
00169     next = LastDerp;
00170     LastDerp = this;
00171     }
00172 
00173 inline Derp::Derp(const double *a1, const ADvari *b1, const ADvari *c1):
00174     a(a1), b(b1), c((ADvari*)c1) {
00175     next = LastDerp;
00176     LastDerp = this;
00177     }
00178 
00179  enum Advari_Opclass {
00180   Hv_const,
00181   Hv_copy,
00182   Hv_binary,
00183   Hv_unary,
00184   Hv_negate,
00185   Hv_plusLR,
00186   Hv_minusLR,
00187   Hv_timesL,
00188   Hv_timesLR,
00189   Hv_quotLR,
00190   Hv_nary
00191   };
00192 
00193  extern ADvari& ADf1(double f, double g, double h, const ADvari &x);
00194  extern ADvari& ADf2(double f, double gx, double gy, double hxx,
00195       double hxy, double hyy, const ADvari &x, const ADvari &y);
00196  extern ADvari& ADfn(double f, int n, const ADvar *x, const double *g, const double *h);
00197 
00198  class
00199 ADvari {  // implementation of an ADvar
00200  public:
00201   static ADcontext adc;
00202   Advari_Opclass opclass;
00203   double Val;   // result of this operation
00204   mutable double aval;  // adjoint -- partial of final result w.r.t. this Val
00205   mutable double dO;  // deriv of op w.r.t. t in x + t*p
00206   mutable double aO;  // adjoint (in Hv computation) of op
00207   mutable double adO; // adjoint (in Hv computation) of dO
00208   void *operator new(size_t len) { return ADvari::adc.Memalloc(len); }
00209   void operator delete(void*) {} /*Should never be called.*/
00210   inline ADvari(Advari_Opclass oc, double t):
00211     opclass(oc), Val(t), aval(0.), dO(0.)
00212     { if (oc != Hv_const) ADvari::adc.ADvari_record(this); }
00213   inline ADvari(Advari_Opclass oc, double t, double ta):
00214     opclass(oc), Val(t), aval(ta), dO(0.)
00215     { if (oc != Hv_const) ADvari::adc.ADvari_record(this); }
00216  private:
00217   inline ADvari(): Val(0.), aval(0.), dO(0.) {} // prevent construction without value (?)
00218  public:
00219   friend class ConstADvari;
00220 #ifdef RAD_AUTO_AD_Const
00221   friend class ADcontext;
00222   friend class ADvar1;
00223   friend class ADvar;
00224   friend class ConstADvar;
00225   friend class IndepADvar;
00226  private:
00227   ADvari *Next;
00228   static ADvari *First_ADvari, **Last_ADvari;
00229  protected:
00230   IndepADvar *padv;
00231  public:
00232   ADvari(const IndepADvar *, double);
00233 #endif
00234 #define F friend
00235 #define R ADvari&
00236 #define Ai const ADvari&
00237 #define T1(r,f) F r f(Ai);
00238 #define T2(r,f) \
00239 F r f(Ai,Ai); \
00240 F r f(double,Ai); \
00241 F r f(Ai,double);
00242   T1(R,operator+)
00243   T2(R,operator+)
00244   T1(R,operator-)
00245   T2(R,operator-)
00246   T2(R,operator*)
00247   T2(R,operator/)
00248   T1(R,abs)
00249   T1(R,acos)
00250   T1(R,acosh)
00251   T1(R,asin)
00252   T1(R,asinh)
00253   T1(R,atan)
00254   T1(R,atanh)
00255   T2(R,atan2)
00256   T2(R,max)
00257   T2(R,min)
00258   T1(R,cos)
00259   T1(R,cosh)
00260   T1(R,exp)
00261   T1(R,log)
00262   T1(R,log10)
00263   T2(R,pow)
00264   T1(R,sin)
00265   T1(R,sinh)
00266   T1(R,sqrt)
00267   T1(R,tan)
00268   T1(R,tanh)
00269   T1(R,fabs)
00270   T1(R,copy)
00271   T2(int,operator<)
00272   T2(int,operator<=)
00273   T2(int,operator==)
00274   T2(int,operator!=)
00275   T2(int,operator>=)
00276   T2(int,operator>)
00277 #undef T2
00278 #undef T1
00279 #undef Ai
00280 #undef R
00281 #undef F
00282 
00283   friend ADvari& ADf1(double f, double g, double h, const ADvari &x);
00284   friend ADvari& ADf2(double f, double gx, double gy, double hxx,
00285       double hxy, double hyy, const ADvari &x, const ADvari &y);
00286   friend ADvari& ADfn(double f, int n, const ADvar *x, const double *g, const double *h);
00287   };
00288 
00289  inline void* Derp::operator new(size_t len) { return ADvari::adc.Memalloc(len); }
00290 
00291 
00292  class
00293 ADvar1: public ADvari { // simplest unary ops
00294  public:
00295   Derp d;
00296   ADvar1(Advari_Opclass oc, double val1, const double *a1, const ADvari *c1):
00297     ADvari(oc,val1), d(a1,this,c1) {}
00298 #ifdef RAD_AUTO_AD_Const
00299   ADvar1(const IndepADvar *, const IndepADvar &);
00300   ADvar1(const IndepADvar *, const ADvari &);
00301 #endif
00302   };
00303 
00304  class
00305 ConstADvari: public ADvari {
00306  private:
00307   ConstADvari *prevcad;
00308   ConstADvari() {}; // prevent construction without value (?)
00309   static ConstADvari *lastcad;
00310  public:
00311   static CADcontext cadc;
00312   inline void *operator new(size_t len) { return ConstADvari::cadc.Memalloc(len); }
00313   inline ConstADvari(double t): ADvari(Hv_copy, t) { prevcad = lastcad; lastcad = this; }
00314   static void aval_reset(void);
00315   };
00316 
00317  class
00318 IndepADvar
00319 {
00320  private:
00321   inline IndepADvar& operator=(const IndepADvar &x) {
00322     /* private to prevent assignment */
00323 #ifdef RAD_AUTO_AD_Const
00324     if (cv)
00325       cv->padv = 0;
00326     cv = new ADvar1(this,x);
00327     return *this;
00328 #else
00329 #ifdef RAD_EQ_ALIAS
00330     cv = x.cv;
00331     return *this;
00332 #else
00333     return ADvar_operatoreq(this,*x.cv);
00334 #endif
00335 #endif /* RAD_AUTO_AD_Const */
00336     }
00337  protected:
00338   static void AD_Const(const IndepADvar&);
00339   ADvari *cv;
00340  public:
00341   typedef double value_type;
00342   friend class ADvar;
00343   friend class ADvar1;
00344   friend class ADvarn;
00345   friend class ADcontext;
00346   IndepADvar(double);
00347   IndepADvar(int);
00348   IndepADvar(long);
00349   IndepADvar& operator=(double);
00350 #ifdef RAD_AUTO_AD_Const
00351   inline IndepADvar(const IndepADvar &x) { cv = x.cv ? new ADvar1(this, x) : 0; };
00352   inline IndepADvar() { cv = 0; }
00353   inline ~IndepADvar() {
00354       if (cv)
00355         cv->padv = 0;
00356     }
00357 #else
00358   inline IndepADvar() {
00359 #ifndef RAD_EQ_ALIAS
00360     cv = 0;
00361 #endif
00362     }
00363   inline ~IndepADvar() {}
00364   friend IndepADvar& ADvar_operatoreq(IndepADvar*, const ADvari&);
00365 #endif
00366 
00367   friend void AD_Const(const IndepADvar&);
00368 
00369   inline operator ADvari&() { return *cv; }
00370   inline operator ADvari&() const { return *(ADvari*)cv; }
00371 
00372   inline double val() const { return cv->Val; }
00373   inline double adj() const { return cv->aval; }
00374   static inline void Gradcomp(int wantgrad)
00375         { ADcontext::Gradcomp(wantgrad); }
00376   static inline void Gradcomp()
00377         { ADcontext::Gradcomp(1); }
00378   static inline void Hvprod(int n, ADvar **vp, double *v, double *hv)
00379         { ADcontext::Hvprod(n, vp, v, hv); }
00380   static inline void aval_reset() { ConstADvari::aval_reset(); }
00381   static inline void Weighted_Gradcomp(int n, ADvar **v, double *w)
00382         { ADcontext::Weighted_Gradcomp(n, v, w); }
00383 
00384 
00385 #define Ai const ADvari&
00386 #define AI const IndepADvar&
00387 #define T2(r,f) \
00388  r f(AI,AI);\
00389  r f(Ai,AI);\
00390  r f(AI,Ai);\
00391  r f(double,AI);\
00392  r f(AI,double);
00393 #define T1(f) friend ADvari& f(AI);
00394 
00395 #define F friend ADvari&
00396 T2(F, operator+)
00397 T2(F, operator-)
00398 T2(F, operator*)
00399 T2(F, operator/)
00400 T2(F, atan2)
00401 T2(F, max)
00402 T2(F, min)
00403 T2(F, pow)
00404 #undef F
00405 #define F friend int
00406 T2(F, operator<)
00407 T2(F, operator<=)
00408 T2(F, operator==)
00409 T2(F, operator!=)
00410 T2(F, operator>=)
00411 T2(F, operator>)
00412 
00413 T1(operator+)
00414 T1(operator-)
00415 T1(abs)
00416 T1(acos)
00417 T1(acosh)
00418 T1(asin)
00419 T1(asinh)
00420 T1(atan)
00421 T1(atanh)
00422 T1(cos)
00423 T1(cosh)
00424 T1(exp)
00425 T1(log)
00426 T1(log10)
00427 T1(sin)
00428 T1(sinh)
00429 T1(sqrt)
00430 T1(tan)
00431 T1(tanh)
00432 T1(fabs)
00433 T1(copy)
00434 
00435 #undef T1
00436 #undef T2
00437 #undef F
00438 #undef Ai
00439 #undef AI
00440 
00441   };
00442 
00443  class
00444 ADvar: public IndepADvar {    // an "active" variable
00445   void ADvar_ctr(double d);
00446  public:
00447   inline ADvar() { /* cv = 0; */ }
00448   inline ADvar(double d) { ADvar_ctr(d); }
00449   inline ADvar(int i) { ADvar_ctr((double)i); }
00450   inline ADvar(long i)  { ADvar_ctr((double)i); }
00451   inline ~ADvar() {}
00452 #ifdef RAD_AUTO_AD_Const
00453   friend class ADvar1;
00454   inline ADvar(const IndepADvar &x) { cv = x.cv ? new ADvar1(this, x) : 0; }
00455   inline ADvar(ADvari &x) { cv = &x; x.padv = this; }
00456   inline ADvar& operator=(const IndepADvar &x) {
00457     if (cv)
00458       cv->padv = 0;
00459     cv = new ADvar1(this,x);
00460     return *this;
00461     }
00462   inline ADvar& operator=(const ADvari &x) {
00463     if (cv)
00464       cv->padv = 0;
00465     cv = new ADvar1(this, x);
00466     return *this;
00467     }
00468 #else
00469   friend ADvar& ADvar_operatoreq(ADvar*, const ADvari&);
00470 #ifdef RAD_EQ_ALIAS
00471   /* allow aliasing v and w after "v = w;" */
00472   inline ADvar(const IndepADvar &x) { cv = x.cv; }
00473   inline ADvar(const ADvari &x) { cv = (ADvari*)&x; }
00474   inline ADvar& operator=(const ADvari &x) { cv = (ADvari*)&x;   return *this; }
00475   inline ADvar& operator=(const IndepADvar &x)  { cv = (ADvari*)x.cv; return *this; }
00476 #else
00477   ADvar(const IndepADvar &x) { cv = x.cv ?
00478     new ADvar1(Hv_copy, x.cv->Val, &CADcontext::One, x.cv) : 0; }
00479   ADvar(const ADvari &x) { cv = new ADvar1(Hv_copy, x.Val, &CADcontext::One, &x); }
00480   inline ADvar& operator=(const ADvari &x) { return ADvar_operatoreq(this,x); }
00481   inline ADvar& operator=(const IndepADvar &x)    { return ADvar_operatoreq(this,*x.cv); }
00482 #endif /* RAD_EQ_ALIAS */
00483 #endif /* RAD_AUTO_AD_Const */
00484   ADvar& operator=(double);
00485   ADvar& operator+=(const ADvari&);
00486   ADvar& operator+=(double);
00487   ADvar& operator-=(const ADvari&);
00488   ADvar& operator-=(double);
00489   ADvar& operator*=(const ADvari&);
00490   ADvar& operator*=(double);
00491   ADvar& operator/=(const ADvari&);
00492   ADvar& operator/=(double);
00493   inline static bool get_fpval_implies_const(void)
00494     { return ConstADvari::cadc.fpval_implies_const; }
00495   inline static void set_fpval_implies_const(bool newval)
00496     { ConstADvari::cadc.fpval_implies_const = newval; }
00497   inline static bool setget_fpval_implies_const(bool newval) {
00498     bool oldval = ConstADvari::cadc.fpval_implies_const;
00499     ConstADvari::cadc.fpval_implies_const = newval;
00500     return oldval;
00501     }
00502   static inline void Gradcomp(int wantgrad)
00503         { ADcontext::Gradcomp(wantgrad); }
00504   static inline void Gradcomp()
00505         { ADcontext::Gradcomp(1); }
00506   static inline void Hvprod(int n, ADvar **vp, double *v, double *hv)
00507         { ADcontext::Hvprod(n, vp, v, hv); }
00508   static inline void aval_reset() { ConstADvari::aval_reset(); }
00509   static inline void Weighted_Gradcomp(int n, ADvar **v, double *w)
00510         { ADcontext::Weighted_Gradcomp(n, v, w); }
00511   };
00512 
00513  inline void AD_Const(const IndepADvar&v) { IndepADvar::AD_Const(v); }
00514 
00515  class
00516 ConstADvar: public ADvar {
00517  private: // disable op=
00518   ConstADvar& operator+=(const ADvari&);
00519   ConstADvar& operator+=(double);
00520   ConstADvar& operator-=(const ADvari&);
00521   ConstADvar& operator-=(double);
00522   ConstADvar& operator*=(const ADvari&);
00523   ConstADvar& operator*=(double);
00524   ConstADvar& operator/=(const ADvari&);
00525   ConstADvar& operator/=(double);
00526   void ConstADvar_ctr(double);
00527  public:
00528   inline ConstADvar(double d) { ConstADvar_ctr(d); }
00529   inline ConstADvar(int i)  { ConstADvar_ctr((double)i); }
00530   inline ConstADvar(long i) { ConstADvar_ctr((double)i); }
00531   ConstADvar(const ADvari &x);
00532 #ifdef RAD_AUTO_AD_Const
00533   ConstADvar(const IndepADvar &x) { cv = new ADvar1(this,x); }
00534 #endif
00535   inline ~ConstADvar(){}
00536 #ifdef RAD_NO_CONST_UPDATE
00537  private:
00538 #endif
00539   ConstADvar();
00540   inline ConstADvar& operator=(double d) { cv->Val = d; return *this; }
00541   inline ConstADvar& operator=(const IndepADvar& d) { cv->Val = d.val(); return *this; }
00542  };
00543 
00544  class
00545 ADvar1s: public ADvar1 { // unary ops with partials
00546  public:
00547   double pL;  // deriv of op w.r.t. left operand L
00548   ADvar1s(double val1, double d1, const ADvari *c1):
00549     ADvar1(Hv_timesL,val1,&pL,c1), pL(d1) {}
00550   };
00551 
00552  class
00553 ADvar1g: public ADvar1 { // unary ops with partials
00554  public:
00555   double pL;  // deriv of op w.r.t. left operand L
00556   double pL2; // partial of op w.r.t. L,L
00557   ADvar1g(double val1, double d1, double d2, const ADvari *c1):
00558     ADvar1(Hv_unary,val1,&pL,c1), pL(d1), pL2(d2) {}
00559   };
00560 
00561  class
00562 ADvar2: public ADvari { // basic binary ops
00563  public:
00564   Derp dL, dR;
00565   ADvar2(Advari_Opclass oc, double val1, const ADvari *Lcv, const double *Lc,
00566       const ADvari *Rcv, const double *Rc): ADvari(oc,val1) {
00567     dR.next = Derp::LastDerp;
00568     dL.next = &dR;
00569     Derp::LastDerp = &dL;
00570     dL.a = Lc;
00571     dL.c = (ADvari*)Lcv;
00572     dR.a = Rc;
00573     dR.c = (ADvari*)Rcv;
00574     dL.b = dR.b = this;
00575     }
00576   };
00577 
00578  class
00579 ADvar2q: public ADvar2 { // binary ops with partials
00580  public:
00581   double pL;  // deriv of op w.r.t. left operand L
00582   double pR;  // deriv of op w.r.t. right operand R
00583   double pLR; // second partial w.r.t. L,R
00584   double pR2; // second partial w.r.t. R,R
00585   ADvar2q(double val1, double Lp, double Rp, double LR, double R2,
00586     const ADvari *Lcv, const ADvari *Rcv);
00587   };
00588 
00589  class
00590 ADvar2g: public ADvar2 { // general binary ops with partials
00591  public:
00592   double pL;  // deriv of op w.r.t. left operand L
00593   double pR;  // deriv of op w.r.t. right operand R
00594   double pL2; // second partial w.r.t. L,L
00595   double pLR; // second partial w.r.t. L,R
00596   double pR2; // second partial w.r.t. R,R
00597   ADvar2g(double val1, double Lp, double Rp, double L2, double LR, double R2,
00598     const ADvari *Lcv, const ADvari *Rcv);
00599   };
00600 
00601  class
00602 ADvarn: public ADvari { // n-ary ops with partials g and 2nd partials h (lower triangle, rowwise)
00603  public:
00604   int n;
00605   double *G, *H;
00606   Derp *D;
00607   ADvarn(double val1, int n1, const ADvar *x, const double *g, const double *h);
00608   };
00609 
00610 inline ADvari &operator+(ADvari &T) { return T; }
00611 inline ADvari &operator+(const ADvari &T) { return (ADvari&) T; }
00612 
00613 inline int operator<(const ADvari &L, const ADvari &R) { return L.Val < R.Val; }
00614 inline int operator<(const ADvari &L, double R) { return L.Val < R; }
00615 inline int operator<(double L, const ADvari &R) { return L < R.Val; }
00616 
00617 inline int operator<=(const ADvari &L, const ADvari &R) { return L.Val <= R.Val; }
00618 inline int operator<=(const ADvari &L, double R) { return L.Val <= R; }
00619 inline int operator<=(double L, const ADvari &R) { return L <= R.Val; }
00620 
00621 inline int operator==(const ADvari &L, const ADvari &R) { return L.Val == R.Val; }
00622 inline int operator==(const ADvari &L, double R) { return L.Val == R; }
00623 inline int operator==(double L, const ADvari &R) { return L == R.Val; }
00624 
00625 inline int operator!=(const ADvari &L, const ADvari &R) { return L.Val != R.Val; }
00626 inline int operator!=(const ADvari &L, double R) { return L.Val != R; }
00627 inline int operator!=(double L, const ADvari &R) { return L != R.Val; }
00628 
00629 inline int operator>=(const ADvari &L, const ADvari &R) { return L.Val >= R.Val; }
00630 inline int operator>=(const ADvari &L, double R) { return L.Val >= R; }
00631 inline int operator>=(double L, const ADvari &R) { return L >= R.Val; }
00632 
00633 inline int operator>(const ADvari &L, const ADvari &R) { return L.Val > R.Val; }
00634 inline int operator>(const ADvari &L, double R) { return L.Val > R; }
00635 inline int operator>(double L, const ADvari &R) { return L > R.Val; }
00636 
00637 inline ADvari& copy(const IndepADvar &x)
00638 { return *(new ADvar1(Hv_copy, x.cv->Val, &CADcontext::One, x.cv)); }
00639 
00640 inline ADvari& copy(const ADvari &x)
00641 { return *(new ADvar1(Hv_copy, x.Val, &CADcontext::One, &x)); }
00642 
00643 inline ADvari& abs(const ADvari &x)
00644 { return fabs(x); }
00645 
00646 #define A (ADvari*)
00647 #define T inline
00648 #define F ADvari&
00649 #define Ai const ADvari&
00650 #define AI const IndepADvar&
00651 #define D double
00652 #define T2(r,f) \
00653  T r f(Ai L, AI R) { return f(L, *A R.cv); }\
00654  T r f(AI L, Ai R) { return f(*A L.cv, R); }\
00655  T r f(AI L, AI R) { return f(*A L.cv, *A R.cv); }\
00656  T r f(AI L, D R) { return f(*A L.cv, R); }\
00657  T r f(D L, AI R) { return f(L, *A R.cv); }
00658 
00659 T2(F, operator+)
00660 T2(F, operator-)
00661 T2(F, operator*)
00662 T2(F, operator/)
00663 T2(F, atan2)
00664 T2(F, pow)
00665 T2(F, max)
00666 T2(F, min)
00667 T2(int, operator<)
00668 T2(int, operator<=)
00669 T2(int, operator==)
00670 T2(int, operator!=)
00671 T2(int, operator>=)
00672 T2(int, operator>)
00673 
00674 #undef T2
00675 #undef D
00676 
00677 #define T1(f)\
00678  T F f(AI x) { return f(*A x.cv); }
00679 
00680 T1(operator+)
00681 T1(operator-)
00682 T1(abs)
00683 T1(acos)
00684 T1(acosh)
00685 T1(asin)
00686 T1(asinh)
00687 T1(atan)
00688 T1(atanh)
00689 T1(cos)
00690 T1(cosh)
00691 T1(exp)
00692 T1(log)
00693 T1(log10)
00694 T1(sin)
00695 T1(sinh)
00696 T1(sqrt)
00697 T1(tan)
00698 T1(tanh)
00699 T1(fabs)
00700 
00701 #undef T1
00702 #undef AI
00703 #undef Ai
00704 #undef F
00705 #undef T
00706 #undef A
00707 
00708 #ifndef SACADO_NO_NAMESPACE
00709 } // namespace Rad2d
00710 } // namespace Sacado
00711 #endif // SACADO_NAMESPACE
00712 #endif // SACADO_RAD2_H
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines