chkder.h
Go to the documentation of this file.
00001 #define chkder_log10e 0.43429448190325182765
00002 #define chkder_factor 100.
00003 
00004 namespace Eigen { 
00005 
00006 namespace internal {
00007 
00008 template<typename Scalar>
00009 void chkder(
00010         const Matrix< Scalar, Dynamic, 1 >  &x,
00011         const Matrix< Scalar, Dynamic, 1 >  &fvec,
00012         const Matrix< Scalar, Dynamic, Dynamic > &fjac,
00013         Matrix< Scalar, Dynamic, 1 >  &xp,
00014         const Matrix< Scalar, Dynamic, 1 >  &fvecp,
00015         int mode,
00016         Matrix< Scalar, Dynamic, 1 >  &err
00017         )
00018 {
00019     typedef DenseIndex Index;
00020 
00021     const Scalar eps = sqrt(NumTraits<Scalar>::epsilon());
00022     const Scalar epsf = chkder_factor * NumTraits<Scalar>::epsilon();
00023     const Scalar epslog = chkder_log10e * log(eps);
00024     Scalar temp;
00025 
00026     const Index m = fvec.size(), n = x.size();
00027 
00028     if (mode != 2) {
00029         /* mode = 1. */
00030         xp.resize(n);
00031         for (Index j = 0; j < n; ++j) {
00032             temp = eps * abs(x[j]);
00033             if (temp == 0.)
00034                 temp = eps;
00035             xp[j] = x[j] + temp;
00036         }
00037     }
00038     else {
00039         /* mode = 2. */
00040         err.setZero(m); 
00041         for (Index j = 0; j < n; ++j) {
00042             temp = abs(x[j]);
00043             if (temp == 0.)
00044                 temp = 1.;
00045             err += temp * fjac.col(j);
00046         }
00047         for (Index i = 0; i < m; ++i) {
00048             temp = 1.;
00049             if (fvec[i] != 0. && fvecp[i] != 0. && abs(fvecp[i] - fvec[i]) >= epsf * abs(fvec[i]))
00050                 temp = eps * abs((fvecp[i] - fvec[i]) / eps - err[i]) / (abs(fvec[i]) + abs(fvecp[i]));
00051             err[i] = 1.;
00052             if (temp > NumTraits<Scalar>::epsilon() && temp < eps)
00053                 err[i] = (chkder_log10e * log(temp) - epslog) / epslog;
00054             if (temp >= eps)
00055                 err[i] = 0.;
00056         }
00057     }
00058 }
00059 
00060 } // end namespace internal
00061 
00062 } // end namespace Eigen


win_eigen
Author(s): Daniel Stonier
autogenerated on Mon Oct 6 2014 12:24:12