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


libicr
Author(s): Robert Krug
autogenerated on Mon Jan 6 2014 11:32:31