ADMB Documentation: fvar_fn1.cpp (* $Id$) Source File

Go to the documentation of this file.
00001 /*
00002  * $Id$
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00011 #include "fvar.hpp"
00012 
00013 #include <stdio.h>
00014 #include <cmath>
00015 
00016 void gradfree(dlink *);
00017 
00018 //extern prevariable * FRETURN;
00019 //extern int RETURN_INDEX;
00020 //extern dlist * GRAD_LIST;          //js
00021 //extern grad_stack  * GRAD_STACK1;
00022 
00029 prevariable& exp(const prevariable& v1)
00030 {
00031   if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00032     gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00033 
00034   double tmp = ::exp(v1.v->x);
00035 
00036 #ifndef OPT_LIB
00037 
00040   #if !defined(__SUNPRO_CC) && !(defined(_MSC_VER) && (_MSC_VER <= 1700))
00041   if (!std::isfinite(tmp))
00042   {
00043     cerr << "Error: Result of \"exp(prevariable(" << value(v1) << ")) = "
00044          << tmp << "\" is not finite.\n";
00045     ad_exit(1);
00046   }
00047   #endif
00048 #endif
00049 
00050   gradient_structure::RETURN_PTR->v->x=tmp;
00051   gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00052     &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x),tmp);
00053 
00054   return *gradient_structure::RETURN_PTR;
00055 }
00056 
00061 prevariable& atan(const prevariable& v1)
00062     {
00063       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00064         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00065       gradient_structure::RETURN_PTR->v->x= ::atan(v1.v->x);
00066       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00067         &(gradient_structure::RETURN_PTR->v->x),
00068         &(v1.v->x),1./(1.+v1.v->x * v1.v->x) );
00069       return(*gradient_structure::RETURN_PTR);
00070     }
00071 
00076 prevariable& ldexp(const prevariable& v1, const int& exponent)
00077     {
00078       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00079         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00080       gradient_structure::RETURN_PTR->v->x=::ldexp(v1.v->x, exponent);
00081       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00082         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x),pow(2.0,exponent));
00083       return(*gradient_structure::RETURN_PTR);
00084     }
00085 
00090 prevariable& sqrt(const prevariable& v1)
00091     {
00092       double tmp=v1.v->x;
00093       if (tmp==0.0)
00094       {
00095         cerr << "Attempting to take the derivative of sqrt(prevariable x)"
00096          " at x=0\n";
00097         ad_exit(1);
00098       }
00099       tmp=::sqrt(tmp);
00100       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00101         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00102       gradient_structure::RETURN_PTR->v->x=tmp;
00103       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00104         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x),1./(2.*tmp));
00105       return(*gradient_structure::RETURN_PTR);
00106     }
00107 
00112 prevariable& sqr(const prevariable& v1)
00113     {
00114       double tmp=v1.v->x;
00115       if (tmp==0.0)
00116       {
00117         cerr << "Attempting to take the derivative of sqrt(prevariable x)"
00118          " at x=0\n";
00119         ad_exit(1);
00120       }
00121       tmp=::sqrt(tmp);
00122       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00123         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00124       gradient_structure::RETURN_PTR->v->x=tmp;
00125       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00126         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x),1./(2.*tmp));
00127       return(*gradient_structure::RETURN_PTR);
00128     }
00129 
00134 prevariable& tan(const prevariable& v1)
00135     {
00136       double t = ::tan(v1.v->x);
00137       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00138         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00139       gradient_structure::RETURN_PTR->v->x= t;
00140       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00141         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x), 1+t*t);
00142       return(*gradient_structure::RETURN_PTR);
00143     }
00144 
00149 prevariable& tanh(const prevariable& v1)
00150     {
00151       double t = ::tanh(v1.v->x);
00152       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00153         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00154       gradient_structure::RETURN_PTR->v->x= t;
00155       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00156         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x), 1-t*t);
00157       return(*gradient_structure::RETURN_PTR);
00158     }
00159 
00164 prevariable& acos(const prevariable& v1)
00165     {
00166       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00167         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00168       gradient_structure::RETURN_PTR->v->x=::acos(v1.v->x);
00169       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00170         &(gradient_structure::RETURN_PTR->v->x),
00171         &(v1.v->x),-1./::sqrt(1.- v1.v->x * v1.v->x));
00172       return(*gradient_structure::RETURN_PTR);
00173     }
00174 
00179 prevariable& asin(const prevariable& v1)
00180     {
00181       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00182         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00183       gradient_structure::RETURN_PTR->v->x=::asin(v1.v->x);
00184       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00185         &(gradient_structure::RETURN_PTR->v->x),
00186         &(v1.v->x),1./::sqrt(1.- v1.v->x * v1.v->x));
00187       return(*gradient_structure::RETURN_PTR);
00188     }
00189 
00194 prevariable& pow(const prevariable& v1, const prevariable& v2)
00195     {
00196       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00197         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00198       double x=::pow(v1.v->x,(v2.v->x)-1);
00199       double y=x* v1.v->x;
00200       gradient_structure::RETURN_PTR->v->x=y;
00201       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00202         &(gradient_structure::RETURN_PTR->v->x),
00203         &(v1.v->x), v2.v->x * x  ,&(v2.v->x),
00204         y * ::log(v1.v->x));
00205       return(*gradient_structure::RETURN_PTR);
00206     }
00207 
00212 prevariable& pow(const double u, const prevariable& v1)
00213     {
00214       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00215         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00216       double y=::pow(u,(v1.v->x));
00217 
00218       gradient_structure::RETURN_PTR->v->x=y;
00219       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00220         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x), y * ::log(u));
00221 
00222       return(*gradient_structure::RETURN_PTR);
00223     }
00224 
00229 prevariable& sinh(const prevariable& v1)
00230     {
00231       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00232         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00233       gradient_structure::RETURN_PTR->v->x=::sinh(v1.v->x);
00234       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00235         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x),::cosh(v1.v->x));
00236       return(*gradient_structure::RETURN_PTR);
00237     }
00238 
00243 prevariable& cosh(const prevariable& v1)
00244     {
00245       if (++gradient_structure::RETURN_PTR > gradient_structure::MAX_RETURN)
00246         gradient_structure::RETURN_PTR = gradient_structure::MIN_RETURN;
00247       gradient_structure::RETURN_PTR->v->x=::cosh(v1.v->x);
00248       gradient_structure::GRAD_STACK1->set_gradient_stack(default_evaluation,
00249         &(gradient_structure::RETURN_PTR->v->x), &(v1.v->x),::sinh(v1.v->x));
00250       return(*gradient_structure::RETURN_PTR);
00251     }
00252 
00257 prevariable& atan2(const prevariable& v1, const prevariable& v2)
00258 {
00259   if (value(v1) == 0 && value(v2) == 0)
00260   {
00261     cerr << "Error: The ADMB function \"atan2(y, x)\" is undefined "
00262     "for y and x equal zero.\n";
00263     ad_exit(1);
00264   }
00265   if (value(v1) == 0 && value(v2) > 0)
00266   {
00267     return atan(v1/v2);
00268   }
00269   dvariable x = (sqrt(v2 * v2 + v1 * v1) - v2)/v1;
00270   return atan(x) * 2.0;
00271 }
00272 
00277 prevariable& atan2(const prevariable& v1, const double v2)
00278 {
00279   if (value(v1) == 0 && v2 == 0)
00280   {
00281     cerr << "Error: The ADMB function \"atan2(y, x)\" is undefined "
00282     "for y and x equal zero.\n";
00283     ad_exit(1);
00284   }
00285   if (value(v1) == 0 && v2 > 0)
00286   {
00287     return atan(v1/v2);
00288   }
00289   dvariable x = (sqrt(v2 * v2 + v1 * v1) - v2)/v1;
00290   return atan(x) * 2.0;
00291 }
00292 
00297 prevariable& atan2(const double v1, const prevariable& v2)
00298 {
00299   if (v1 == 0 && value(v2) == 0)
00300   {
00301     cerr << "Error: The ADMB function \"atan2(y, x)\" is undefined "
00302     "for y and x equal zero.\n";
00303     ad_exit(1);
00304   }
00305   if (v1 == 0 && value(v2) > 0)
00306   {
00307     return atan(v1/v2);
00308   }
00309   dvariable x = (sqrt(v2 * v2 + v1 * v1) - v2)/v1;
00310   return atan(x) * 2.0;
00311 }