44#if defined(RAD_DEBUG_BLOCKKEEP) && !defined(HAVE_SACADO_UNINIT)
45#undef RAD_DEBUG_BLOCKKEEP
49#ifndef RAD_AUTO_AD_Const
50#define RAD_AUTO_AD_Const
64#ifndef RAD_AUTO_AD_Const
65#ifdef RAD_DEBUG_BLOCKKEEP
74#ifndef RAD_NO_USING_STDCC
95#ifdef RAD_AUTO_AD_Const
96#undef RAD_DEBUG_BLOCKKEEP
98#ifdef RAD_DEBUG_BLOCKKEEP
99#if !(RAD_DEBUG_BLOCKKEEP > 0)
100#undef RAD_DEBUG_BLOCKKEEP
102extern "C" void _uninit_f2c(
void *
x,
int type,
long len);
108struct UninitType<float> {
109 static const int utype = 4;
113struct UninitType<double> {
114 static const int utype = 5;
118struct UninitType<
std::complex<T> > {
119 static const int utype = UninitType<T>::utype + 2;
128 template<
typename T>
class
141#define Dtype typename DoubleAvoid<Double>::dtype
142#define Ttype typename DoubleAvoid<Double>::ttype
144 template<
typename Double>
class IndepADvar;
145 template<
typename Double>
class ConstADvar;
146 template<
typename Double>
class ConstADvari;
147 template<
typename Double>
class ADvar;
148 template<
typename Double>
class ADvar1;
149 template<
typename Double>
class ADvar1g;
150 template<
typename Double>
class ADvar1s;
151 template<
typename Double>
class ADvar2;
152 template<
typename Double>
class ADvar2g;
153 template<
typename Double>
class ADvar2q;
154 template<
typename Double>
class ADvari;
155 template<
typename Double>
class ADvari_block;
156 template<
typename Double>
class ADvarn;
157 template<
typename Double>
class Derp;
159 template<
typename Double>
struct
167 template<
typename Double>
class
171 enum { Gulp = 1021 };
177 template<
typename Double>
class
196 double First0[(
sizeof(
ADMemblock) +
sizeof(
double) - 1) /
sizeof(
double)];
197 double First1[(
sizeof(
ADVari_block) +
sizeof(
double) - 1) /
sizeof(
double)];
198 void *new_ADmemblock(
size_t);
199 void new_ADvari_block();
204#ifdef RAD_DEBUG_BLOCKKEEP
211 void *Memalloc(
size_t len);
212 static void Gradcomp(
int wantgrad);
213 static void Hvprod(
int,
ADVar**, Double*, Double*);
215 static void Weighted_Gradcomp(
int,
ADVar**, Double*);
218 if (Ainext >= Ailimit)
224 template<
typename Double>
class
227 bool fpval_implies_const;
229 friend class ADvar<Double>;
233 template<
typename Double>
class
236 friend class ADvarn<Double>;
244 Derp(
const ADVari *);
245 Derp(
const Double *,
const ADVari *);
246 Derp(
const Double *,
const ADVari *,
const ADVari *);
253#define Ai const ADvari<Double>&
254#define AI const IndepADvar<Double>&
255#define T template<typename Double>
282#define F ADvari<Double>&
331#define SNS Sacado::Rad2
389 template<
typename Double> ADvari<Double>&
390ADf1(Double f, Double g, Double h,
const ADvari<Double> &
x);
392 template<
typename Double> ADvari<Double>&
393ADf2(Double f, Double gx, Double gy, Double hxx,
394 Double hxy, Double hyy,
const ADvari<Double> &
x,
const ADvari<Double> &
y);
396 template<
typename Double> ADvari<Double>&
397ADfn(Double f,
int n,
const IndepADvar<Double> *
x,
const Double *g,
const Double *h);
399 template<
typename Double>
class
405#ifdef RAD_AUTO_AD_Const
415 static ADvari *First_ADvari, **Last_ADvari;
417 inline void ADvari_padv(
const IndepADVar *v) {
419 Last_ADvari = &this->Next;
426 static int gcgen_cur, last_opno, zap_gcgen, zap_gcgen1, zap_opno;
427 static FILE *debug_file;
436 void *
operator new(
size_t len) {
439 rv->gcgen = gcgen_cur;
440 rv->opno = ++last_opno;
441 if (last_opno == zap_opno && gcgen_cur == zap_gcgen)
448 void operator delete(
void*) {}
450 opclass(oc), Val(t), aval(0.), dO(0.)
453 opclass(oc), Val(t), aval(ta), dO(0.)
456 inline ADvari(): Val(0.), aval(0.), dO(0.) {}
459#ifdef RAD_AUTO_AD_Const
461 friend class ADvar<Double>;
462 friend class ADvar1<Double>;
464 friend class ADvar2<Double>;
467 ADvari(
const IndepADVar *, Double);
468 ADvari(
const IndepADVar *, Double, Double);
469 ADvari(
const IndepADVar *, Double, Double,
int);
473#define Ai const ADvari&
474#define T1(r,f) F r f <>(Ai);
479F r f <>(double,Ai); \
480F r f <>(Ai,double); \
526 friend ADvari& ADf1<>(Double f, Double g, Double h,
const ADvari &
x);
527 friend ADvari& ADf2<>(Double f, Double gx, Double gy, Double hxx,
528 Double hxy, Double hyy,
const ADvari &
x,
const ADvari &
y);
529 friend ADvari& ADfn<>(Double f,
int n,
const IndepADVar *
x,
530 const Double *g,
const Double *h);
532 inline operator Double() {
return this->Val; }
533 inline operator Double()
const {
return this->Val; }
536 template<
typename Double>
class
542 ADVari(oc,val1), d(a1,this,c1) {}
543#ifdef RAD_AUTO_AD_Const
544 typedef typename ADVari::IndepADVar IndepADVar;
546 ADvar1(
const IndepADVar*,
const IndepADVar&);
547 ADvar1(
const IndepADVar*,
const ADVari&);
549 const ADVari *c1,
const ADVar *v):
550 ADVari(oc,val1), d(a1,this,c1) {
552 this->ADvari_padv(v);
558 template<
typename Double>
class
570 static void aval_reset(
void);
574 template<
typename Double>
class
579#ifdef RAD_AUTO_AD_Const
583#elif defined(RAD_EQ_ALIAS)
595 friend class ADvar<Double>;
597 friend class ADvar1<Double>;
598 friend class ADvarn<Double>;
606#ifdef RAD_AUTO_AD_Const
610 inline IndepADvar() { cv = 0; }
611 inline ~IndepADvar() {
626 inline operator ADVari&()
const {
627 ADVari *tcv = this->cv;
632 inline operator ADVari*()
const {
633 ADVari *tcv = this->cv;
639 inline operator ADVari&()
const {
return *this->cv; }
640 inline operator ADVari*()
const {
return this->cv; }
646 friend void AD_Const1<>(Double*,
const IndepADvar&);
657 static inline void Hvprod(
int n,
ADVar **vp, Double *v, Double *hv)
667#define Ai const ADVari&
668#define AI const IndepADvar&
682#define T1(f) friend ADVari& f<> (AI);
684#define F friend ADVari&
733 template<
typename Double>
class
743 void ADvar_ctr(Double d) {
745 if (ConstADVari::cadc.fpval_implies_const)
748#ifdef RAD_AUTO_AD_Const
750 x->ADvari_padv(
this);
758 friend class ADvar1<Double>;
766#ifdef RAD_AUTO_AD_Const
768 this->cv =
x.cv ?
new ADVar1(
this,
x) : 0;
770 ADvar(ADVari &
x) :IndepADVar(&
x) {
x.ADvari_padv((IndepADVar*)
this);}
771 inline ADvar& operator=(IndepADVar &
x) {
774 this->cv =
new ADVar1(
this,
x);
777 inline ADvar& operator=(ADVari &
x) {
780 this->cv =
new ADVar1(
this,
x);
784 friend ADvar& ADvar_operatoreq<>(ADvar*,
const ADVari&);
787 inline ADvar(
const IndepADVar &
x) { this->cv = (ADVari*)
x.cv; }
788 inline ADvar(
const ADVari &
x) { this->cv = (ADVari*)&
x; }
789 inline ADvar& operator=(IndepADVar &
x) { this->cv = (ADVari*)
x.cv;
return *
this; }
790 inline ADvar& operator=(
const ADVari &
x) { this->cv = (ADVari*)&
x;
return *
this; }
792 friend ADvar& ADvar_operatoreq<>(ADvar*,
const ADVari&);
801 ADvar& operator=(Double);
802 ADvar& operator+=(
const ADVari&);
803 ADvar& operator+=(Double);
804 ADvar& operator-=(
const ADVari&);
805 ADvar& operator-=(Double);
806 ADvar& operator*=(
const ADVari&);
807 ADvar& operator*=(Double);
808 ADvar& operator/=(
const ADVari&);
809 ADvar& operator/=(Double);
811 {
return ConstADVari::cadc.fpval_implies_const; }
813 { ConstADVari::cadc.fpval_implies_const = newval; }
815 bool oldval = ConstADVari::cadc.fpval_implies_const;
816 ConstADVari::cadc.fpval_implies_const = newval;
823 static inline void aval_reset() { ConstADVari::aval_reset(); }
828template<
typename Double>
832template<
typename Double>
835 template<
typename Double>
class
852 void ConstADvar_ctr(Double);
862#ifdef RAD_NO_CONST_UPDATE
870 template<
typename Double>
class
878#ifdef RAD_AUTO_AD_Const
880 typedef typename ADVar1::ADVar
ADVar;
881 ADvar1s(Double val1, Double a1,
const ADVari *c1,
const ADVar *v):
886 template<
typename Double>
class
897 template<
typename Double>
class
904 const ADVari *Rcv,
const Double *Rc):
906 dR.
next = DErp::LastDerp;
908 DErp::LastDerp = &dL;
915#ifdef RAD_AUTO_AD_Const
918 const ADVari *Rcv,
const Double *Rc,
ADVar *v):
920 dR.
next = DErp::LastDerp;
922 DErp::LastDerp = &dL;
929 this->ADvari_padv(v);
934 template<
typename Double>
class
944 ADvar2q(Double val1, Double Lp, Double Rp, Double LR, Double R2,
947 pL(Lp), pR(Rp), pLR(LR), pR2(R2) {}
948#ifdef RAD_AUTO_AD_Const
949 typedef typename ADVar2::ADVar
ADVar;
950 ADvar2q(Double val1, Double Lp, Double Rp, Double LR, Double R2,
951 const ADVari *Lcv,
const ADVari *Rcv,
const ADVar *v):
952 ADVar2(
Hv_quotLR,val1,Lcv,&pL,Rcv,&pR,v),
953 pL(Lp), pR(Rp), pLR(LR), pR2(R2) {}
957 template<
typename Double>
class
967 ADvar2g(Double val1, Double Lp, Double Rp, Double L2, Double LR, Double R2,
970 pL(Lp), pR(Rp), pL2(L2), pLR(LR), pR2(R2) { }
973 template<
typename Double>
class
990 a1 = G = (Double*)ADVari::adc.Memalloc(n1*
sizeof(*g));
991 d1 =
D = (
DErp*)ADVari::adc.Memalloc(n1*
sizeof(
DErp));
992 dlast = DErp::LastDerp;
993 for(
i = 0;
i < n1;
i++, d1++) {
1001 DErp::LastDerp = dlast;
1002 nh = (n1*(n1+1)) >> 1;
1003 a1 = H = (
double*)ADVari::adc.Memalloc(nh *
sizeof(*h));
1004 for(
i = 0;
i < nh;
i++)
1009template<
typename Double>
1012template<
typename Double>
1014template<
typename Double>
1016template<
typename Double>
1019template<
typename Double>
1021template<
typename Double>
1023template<
typename Double>
1026template<
typename Double>
1028template<
typename Double>
1030template<
typename Double>
1033template<
typename Double>
1035template<
typename Double>
1037template<
typename Double>
1040template<
typename Double>
1042template<
typename Double>
1044template<
typename Double>
1047template<
typename Double>
1049template<
typename Double>
1051template<
typename Double>
1054template<
typename Double>
1057 return Mbase + (Mleft -= len);
1058 return new_ADmemblock(len);
1061template<
typename Double>
1067template<
typename Double>
1073template<
typename Double>
1089#ifdef RAD_AUTO_AD_Const
1095#ifndef RAD_DEBUG_gcgen1
1096#define RAD_DEBUG_gcgen1 -1
1107 template<
typename Double>
1116 Mbase = (
char*)First->memblk;
1117 Mleft =
sizeof(First->memblk);
1122 fb->
limit = Ailimit = fb->
pADvari + ADVari_block::Gulp;
1123 rad_need_reinit = 0;
1124#ifdef RAD_DEBUG_BLOCKKEEP
1125 rad_busy_blocks = 0;
1131template<
typename Double>
void*
1136#ifdef RAD_AUTO_AD_Const
1139#ifdef RAD_Const_WARN
1145 if ((rad_need_reinit & 1) &&
this == &ADVari::adc) {
1146 rad_need_reinit &= ~1;
1148 Aibusy = b = AiFirst;
1151 Ailimit = b->
limit = (Ainext = b->
pADvari) + ADVari_block::Gulp;
1152#ifdef RAD_DEBUG_BLOCKKEEP
1153 Mleft = rad_mleft_save;
1154 if (Mleft <
sizeof(First->memblk))
1156 UninitType<Double>::utype,
1157 (
sizeof(First->memblk) - Mleft)
1159 if ((mb = Busy->
next)) {
1160 if (!(mb0 = rad_Oldcurmb))
1162 for(;; mb = mb->
next) {
1164 UninitType<Double>::utype,
1165 sizeof(First->memblk)
1171 rad_Oldcurmb = Busy;
1172 if (rad_busy_blocks >= RAD_DEBUG_BLOCKKEEP) {
1173 rad_busy_blocks = 0;
1177 for(mb = Busy; mb != mb0; mb = mb1) {
1184 Mbase = (
char*)First->
memblk;
1185 Mleft =
sizeof(First->memblk);
1192 for(mb = Busy; mb != mb0; mb = mb1) {
1199 Mbase = (
char*)First->
memblk;
1200 Mleft =
sizeof(First->memblk);
1201#ifdef RAD_AUTO_AD_Const
1202 if (ADVari::adc.rad_need_reinit & 2) {
1203 ADVari::adc.rad_need_reinit &= ~2;
1204 *ADVari::Last_ADvari = 0;
1205 ADVari::Last_ADvari = &ADVari::First_ADvari;
1206 anext = ADVari::First_ADvari;
1208 while((
a = anext)) {
1211#ifdef RAD_Const_WARN
1212 if ((
i =
a->opno) > 0)
1228 return Mbase + (Mleft -= len);
1235#ifdef RAD_DEBUG_BLOCKKEEP
1240 return (Mbase = (
char*)
x->memblk) +
1241 (Mleft =
sizeof(First->memblk) - len);
1244template<
typename Double>
void
1249 ob->
limit = Ailimit;
1254 Aibusy = Aibusy->
next = nb;
1255 nb->
limit = Ailimit = (Ainext = nb->
pADvari) + ADVari_block::Gulp;
1261template<
typename Double>
void
1266 if (ADVari::adc.rad_need_reinit) {
1267 for(d = DErp::LastDerp; d; d = d->
next)
1270 if (!(ADVari::adc.rad_need_reinit & 1)) {
1271 ADVari::adc.rad_need_reinit = 1;
1272 ADVari::adc.rad_mleft_save = ADVari::adc.Mleft;
1273 ADVari::adc.Mleft = 0;
1276 if (ADVari::gcgen_cur == ADVari::zap_gcgen1) {
1278 if (!(fname = getenv(
"RAD_DEBUG_FILE")))
1279 fname =
"rad_debug.out";
1283 ADVari::debug_file = fopen(fname,
"w");
1284 ADVari::zap_gcgen1 = -1;
1287 if ((d = DErp::LastDerp) != 0) {
1288#ifdef RAD_AUTO_AD_Const
1289 ADVari::adc.rad_need_reinit |= 2;
1294 if (ADVari::debug_file)
1296 fprintf(ADVari::debug_file,
"%d\t%d\t%g + %g * %g",
1297 d->
c->opno, d->
b->opno, d->
c->
aval,
1300 fprintf(ADVari::debug_file,
" = %g\n", d->
c->
aval);
1301 fflush(ADVari::debug_file);
1302 }
while((d = d->
next));
1306 while((d = d->
next));
1311 if (ADVari::debug_file) {
1312 fclose(ADVari::debug_file);
1313 ADVari::debug_file = 0;
1317 ADVari::gcgen_cur++;
1318 ADVari::last_opno = 0;
1322template<
typename Double>
void
1327#ifdef RAD_Const_WARN
1331#ifdef RAD_AUTO_AD_Const
1336 if (ADVari::adc.rad_need_reinit) {
1337 for(d = DErp::LastDerp; d; d = d->
next)
1340 if (!(ADVari::adc.rad_need_reinit & 1)) {
1341 ADVari::adc.rad_need_reinit = 1;
1342 ADVari::adc.rad_mleft_save = ADVari::adc.Mleft;
1343 ADVari::adc.Mleft = 0;
1346 if (ADVari::gcgen_cur == ADVari::zap_gcgen1) {
1348 if (!(fname = getenv(
"RAD_DEBUG_FILE")))
1349 fname =
"rad_debug.out";
1353 ADVari::debug_file = fopen(fname,
"w");
1354 ADVari::zap_gcgen1 = -1;
1357 if ((d = DErp::LastDerp) != 0) {
1358 for(
i = 0;
i < n;
i++)
1361 if (ADVari::debug_file)
1363 fprintf(ADVari::debug_file,
"%d\t%d\t%g + %g * %g",
1366 fprintf(ADVari::debug_file,
" = %g\n", d->
c->
aval);
1367 fflush(ADVari::debug_file);
1368 }
while((d = d->
next));
1372 while((d = d->
next));
1375 if (ADVari::debug_file) {
1376 fclose(ADVari::debug_file);
1377 ADVari::debug_file = 0;
1380#ifdef RAD_AUTO_AD_Const
1381 *ADVari::Last_ADvari = 0;
1382 ADVari::Last_ADvari = &ADVari::First_ADvari;
1383 if ((anext = ADVari::First_ADvari) && !(ADVari::adc.rad_need_reinit & 2)) {
1384 ADVari::adc.rad_need_reinit = 3;
1385 while((
a = anext)) {
1388#ifdef RAD_Const_WARN
1389 if ((
i =
a->opno) > 0)
1404 ADVari::gcgen_cur++;
1405 ADVari::last_opno = 0;
1409 template<
typename Double>
1417 template<
typename Double>
1425 template<
typename Double>
1433 template<
typename Double>
1441 template<
typename Double>
1448 template<
typename Double>
void
1455 template<
typename Double>
1463 template<
typename Double>
1471 template<
typename Double>
1479 template<
typename Double>
1485 ConstADVari *ncv =
new ConstADVari(v.
val());
1486#ifdef RAD_AUTO_AD_Const
1493 template<
typename Double>
1504#ifdef RAD_AUTO_AD_Const
1506 template<
typename Double>
1510 this->ADvari_padv(
x);
1513 template<
typename Double>
1517 this->ADvari_padv(
x);
1520 template<
typename Double>
1522 ADVari(
Hv_copy,
y.cv->Val), d((const Double*)&ADcontext<Double>::One, (ADVari*)this,
y.cv)
1524 this->ADvari_padv(
x);
1527 template<
typename Double>
1529 ADVari(
Hv_copy,
y.Val), d((const Double*)&ADcontext<Double>::One, this, &
y)
1531 this->ADvari_padv(
x);
1536 template<
typename Double>
1542 template<
typename Double>
1547 template<
typename Double>
1551#ifdef RAD_AUTO_AD_Const
1562 template<
typename Double>
1566#ifdef RAD_AUTO_AD_Const
1572 this->cv = ConstADVari::cadc.fpval_implies_const
1579 template<
typename Double>
1585 template<
typename Double>
1591#ifdef RAD_AUTO_AD_Const
1592#define RAD_ACA ,this
1597 template<
typename Double>
1606 template<
typename Double>
1612 template<
typename Double>
1620 template<
typename Double>
1626 template<
typename Double>
1632 template<
typename Double>
1641 template<
typename Double>
1647 template<
typename Double>
1655 template<
typename Double>
1661 template<
typename Double>
1667 template<
typename Double>
1676 template<
typename Double>
1682 template<
typename Double>
1690 template<
typename Double>
1696 template<
typename Double>
1699 Double Lv = L.
Val, Rv =
R.Val, pL = 1. / Rv, q = Lv/Rv, qpL = q*pL;
1703 template<
typename Double>
1707 Double Lv = Lcv->
Val, Rv =
R.Val, pL = 1. / Rv, q = Lv/Rv, qpL = q*pL;
1712 template<
typename Double>
1718 template<
typename Double>
1721 Double recip = 1. /
R.Val;
1722 Double q = L * recip;
1723 Double d1 = -q*recip;
1727 template<
typename Double>
1735 template<
typename Double>
1738 Double t = v.
Val, t1 = 1. - t*t, d1 = -1./std::sqrt(t1);
1742 template<
typename Double>
1745 Double d1, t, t1, t2;
1747 t1 = std::sqrt(t2 = t*t - 1.);
1752 template<
typename Double>
1757 d1 = 1. / std::sqrt(t1 = 1. - t*t);
1761 template<
typename Double>
1764 Double d1, t, t1, t2, td;
1766 t1 = std::sqrt(t2 = t*t + 1.);
1771 return *(
new ADvar1g<Double>(td*std::log(t*td + t1), d1, -(t/t2)*d1, &v));
1774 template<
typename Double>
1777 Double t = v.
Val, d1 = 1./(1. + t*t);
1781 template<
typename Double>
1784 Double t = v.
Val, d1 = 1./(1. - t*t);
1785 return *(
new ADvar1g<Double>(0.5*std::log((1.+t)/(1.-t)), d1, (t+t)*d1*d1, &v));
1788 template<
typename Double>
1791 Double R2, t, t2,
x, x2,
y, y2;
1802 template<
typename Double>
1805 Double t, x2,
y, y2;
1813 template<
typename Double>
1816 Double t,
x, x2, y2;
1824 template<
typename Double>
1831 template<
typename Double>
1839 template<
typename Double>
1847 template<
typename Double>
1854 template<
typename Double>
1862 template<
typename Double>
1870 template<
typename Double>
1873 Double t = std::cos(v.
Val);
1877 template<
typename Double>
1880 Double t = std::cosh(v.
Val);
1884 template<
typename Double>
1887 Double t = std::exp(v.
Val);
1891 template<
typename Double>
1894 Double
x = v.
Val, d1 = 1. /
x;
1898 template<
typename Double>
1901 static double num = 1. / std::log(10.);
1909 template<
typename Double>
1912 Double
dx, dy, t,
x, xlog, xym1,
y;
1923 template<
typename Double>
1926 Double dy, t, xlog,
y;
1934 template<
typename Double>
1944 template<
typename Double>
1947 Double t = std::sin(v.
Val);
1951 template<
typename Double>
1954 Double t = std::sinh(v.
Val);
1958 template<
typename Double>
1961 Double t = std::sqrt(v.
Val);
1962 Double d1 = 0.5 / t;
1966 template<
typename Double>
1970 rv = std::tan(v.
Val);
1971 t = 1. / std::cos(v.
Val);
1976 template<
typename Double>
1980 rv = std::tanh(v.
Val);
1981 t = 1. / std::cosh(v.
Val);
1986 template<
typename Double>
1991 if ((t = v.
Val) < 0) {
1998 template<
typename Double>
2005 if ((t = v.
Val) < 0) {
2012 template<
typename Double>
2018 template<
typename Double>
2024 template<
typename Double>
2026ADf2(Double f, Double gx, Double gy, Double hxx, Double hxy, Double hyy,
2031 template<
typename Double>
2033ADf2(Double f, Double gx, Double gy, Double hxx, Double hxy, Double hyy,
2038 template<
typename Double>
2040ADf2(Double f, Double gx, Double gy, Double hxx, Double hxy, Double hyy,
2045 template<
typename Double>
2047ADf2(Double f, Double gx, Double gy, Double hxx, Double hxy, Double hyy,
2052 template<
typename Double>
2058 template<
typename Double>
2064 template<
typename Double>
2068 ADVari *
a, *aL, *aR, **ap, **ape;
2071 Double aO, adO, *g, *h, *h0, t, tL, tR;
2073 for(
i = 0;
i < n;
i++) {
2076 a->aO =
a->adO = 0.;
2078 ADVari::adc.Aibusy->limit = ADVari::adc.Ainext;
2080 for(b0 = 0, b = ADVari::adc.AiFirst; b; b0 = b, b = b->
next) {
2085 a->
aO =
a->adO = 0.;
2086 switch(
a->opclass) {
2122 for(
i = 0;
i < m;
i++)
2123 t += g[
i] * d[
i].
c->dO;
2133 for(b = b0; b; b = b->
prev) {
2140 switch(
a->opclass) {
2190 aL->
adO += adO * tL;
2195 aL->
aO += aO * (tL = aR->
Val) + adO*aR->
dO;
2196 aR->
aO += aO * (tR = aL->
Val) + adO*aL->
dO;
2197 aL->
adO += adO * tL;
2198 aR->
adO += adO * tR;
2218 for(
i = 0;
i < m;
i++) {
2220 aL->
adO += adO * (t = g[
i]);
2223 for(h = h0, j = 0; j <=
i; j++)
2224 d[j].
c->aO += t * *h++;
2226 for(k = j; j < m; j++)
2227 d[j].
c->aO += t * *(h += k++);
2234 for(
i = 0;
i < n;
i++) {
2241 template<
typename Double>
2248#define A (ADvari<Double>*)
2249#ifdef RAD_Const_WARN
2250#define C(x) (((x)->opno < 0) ? RAD_Const_Warn(x) : 0, *A x)
2254#define T template<typename Double> inline
2255#define F ADvari<Double>&
2256#define Ai const ADvari<Double>&
2257#define AI const IndepADvar<Double>&
2260 T r f(Ai L, AI R) { return f(L, C(R.cv)); }\
2261 T r f(AI L, Ai R) { return f(C(L.cv), R); }\
2262 T r f(AI L, AI R) { return f(C(L.cv), C(R.cv)); }\
2263 T r f(AI L, D R) { return f(C(L.cv), R); }\
2264 T r f(Ai L, Dtype R) { return f(L, (D)R); }\
2265 T r f(AI L, Dtype R) { return f(C(L.cv), (D)R); }\
2266 T r f(Ai L, long R) { return f(L, (D)R); }\
2267 T r f(AI L, long R) { return f(C(L.cv), (D)R); }\
2268 T r f(Ai L, int R) { return f(L, (D)R); }\
2269 T r f(AI L, int R) { return f(C(L.cv), (D)R); }\
2270 T r f(D L, AI R) { return f(L, C(R.cv)); }\
2271 T r f(Dtype L, Ai R) { return f((D)L, R); }\
2272 T r f(Dtype L, AI R) { return f((D)L, C(R.cv)); }\
2273 T r f(long L, Ai R) { return f((D)L, R); }\
2274 T r f(long L, AI R) { return f((D)L, C(R.cv)); }\
2275 T r f(int L, Ai R) { return f((D)L, R); }\
2276 T r f(int L, AI R) { return f((D)L, C(R.cv)); }
2297 T F f(AI x) { return f(C(x.cv)); }
int RAD_Const_Warn(void *v)
expr expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c *expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 c
static void aval_reset(void)
ADvari_block< Double > ADVari_block
ADvar1g< Double > ADVar1g
ADvar2g< Double > ADVar2g
ADvar2q< Double > ADVar2q
static void Weighted_Gradcomp(int, ADVar **, Double *)
static void Hvprod(int, ADVar **, Double *, Double *)
ADvar1s< Double > ADVar1s
void * Memalloc(size_t len)
ADmemblock< Double > ADMemblock
static const Double negOne
void * new_ADmemblock(size_t)
void ADvari_record(ADVari *x)
ADvar1(Advari_Opclass oc, Double val1, const Double *a1, const ADVari *c1)
ADvar1g(Double val1, Double d1, Double d2, const ADVari *c1)
ADvar1s(Double val1, Double a1, const ADVari *c1)
ADvar2(Advari_Opclass oc, Double val1, const ADVari *Lcv, const Double *Lc, const ADVari *Rcv, const Double *Rc)
ADvar2g(Double val1, Double Lp, Double Rp, Double L2, Double LR, Double R2, const ADVari *Lcv, const ADVari *Rcv)
ADvar2q(Double val1, Double Lp, Double Rp, Double LR, Double R2, const ADVari *Lcv, const ADVari *Rcv)
ConstADvari< Double > ConstADVari
static void Weighted_Gradcomp(int n, ADvar **v, Double *w)
ADvar & operator/=(const ADVari &)
ADvar & operator=(const ADVari &x)
IndepADvar< Double > IndepADVar
ADvar & operator*=(const ADVari &)
static void Gradcomp(int wantgrad)
ADvar(const IndepADVar &x)
static void set_fpval_implies_const(bool newval)
ADvar & operator-=(const ADVari &)
ADvar & operator+=(const ADVari &)
static bool get_fpval_implies_const(void)
static bool setget_fpval_implies_const(bool newval)
ADvar & operator=(IndepADVar &x)
ScalarType< value_type >::type scalar_type
IndepADvar< Double > IndepADVar
ADvari(Advari_Opclass oc, Double t)
static ADcontext< Double > adc
ADvari(Advari_Opclass oc, Double t, Double ta)
ADVari::IndepADVar IndepADVar
ADvarn(Double val1, int n1, const IndepADVar *x, const Double *g, const Double *h)
ConstADvar & operator=(ADVari &d)
ADVar::IndepADVar IndepADVar
void ConstADvar_ctr(Double)
ConstADvar & operator*=(Double)
ConstADvar & operator+=(ADVari &)
ConstADvar & operator-=(Double)
ConstADvar & operator*=(ADVari &)
ConstADvar & operator/=(ADVari &)
ADVar::ConstADVari ConstADVari
ConstADvar & operator-=(ADVari &)
ConstADvar & operator+=(Double)
ConstADvar & operator/=(Double)
ConstADvar & operator=(Double d)
static ConstADvari * lastcad
static void aval_reset(void)
static CADcontext< Double > cadc
IndepADvar & operator=(IndepADvar &x)
static void AD_Const(const IndepADvar &)
static void Hvprod(int n, ADVar **vp, Double *v, Double *hv)
static void Weighted_Gradcomp(int n, ADVar **v, Double *w)
static void Gradcomp(int wantgrad)
ADvari< Double > & operator+(const ADvari< Double > &T)
ADvari< Double > & atan(const ADvari< Double > &v)
ADvari< Double > & abs(const ADvari< Double > &v)
void AD_Const1(Double *, const IndepADvar< Double > &)
ADvari< Double > & exp(const ADvari< Double > &v)
void AD_Const(const IndepADvar< Double > &)
ADvari< Double > & atanh(const ADvari< Double > &v)
ADvari< Double > & atan2(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & asinh(const ADvari< Double > &v)
ADvari< Double > & max(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & acosh(const ADvari< Double > &v)
ADvari< Double > & pow(const ADvari< Double > &L, const ADvari< Double > &R)
int operator<=(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & operator*(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & tan(const ADvari< Double > &v)
ADvari< Double > & tanh(const ADvari< Double > &v)
ADvari< Double > & log10(const ADvari< Double > &v)
ADvari< Double > & ADf2(Double f, Double gx, Double gy, const IndepADvar< Double > &x, const IndepADvar< Double > &y)
ADvari< Double > & min(const ADvari< Double > &L, const ADvari< Double > &R)
int operator>(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & log(const ADvari< Double > &v)
ADvari< Double > & sqrt(const ADvari< Double > &v)
int operator==(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & fabs(const ADvari< Double > &v)
int operator>=(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & operator/(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & sinh(const ADvari< Double > &v)
ADvari< Double > & operator-(const ADvari< Double > &T)
ADvari< Double > & ADfn(Double f, int n, const IndepADvar< Double > *x, const Double *g)
ADvari< Double > & ADf1(Double f, Double g, const IndepADvar< Double > &x)
int operator<(const ADvari< Double > &L, const ADvari< Double > &R)
ADvari< Double > & acos(const ADvari< Double > &v)
ADvari< Double > & sin(const ADvari< Double > &v)
ADvari< Double > & cosh(const ADvari< Double > &v)
IndepADvar< Double > & ADvar_operatoreq(IndepADvar< Double > *, const ADvari< Double > &)
ADvari< Double > & cos(const ADvari< Double > &v)
ADvari< Double > & asin(const ADvari< Double > &v)
int operator!=(const ADvari< Double > &L, const ADvari< Double > &R)
Turn ADvar into a meta-function class usable with mpl::apply.
ADT_RAD IndepADvar< double > AI
ADT_RAD ADvari< double > Ai
Sacado::RadVec::ADvar< double > ADVar
void _uninit_f2c(void *x, int type, long len)