CppAD: A C++ Algorithmic Differentiation Package  20130102
near_equal.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_NEAR_EQUAL_INCLUDED
00003 # define CPPAD_NEAR_EQUAL_INCLUDED
00004 
00005 /* --------------------------------------------------------------------------
00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-12 Bradley M. Bell
00007 
00008 CppAD is distributed under multiple licenses. This distribution is under
00009 the terms of the 
00010                     Eclipse Public License Version 1.0.
00011 
00012 A copy of this license is included in the COPYING file of this distribution.
00013 Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
00014 -------------------------------------------------------------------------- */
00015 
00016 /*
00017 $begin NearEqual$$
00018 $spell
00019      cppad.hpp
00020      sqrt
00021      cout
00022      endl
00023      Microsoft
00024      std
00025      Cpp
00026      namespace
00027      const
00028      bool
00029 $$
00030 
00031 $section Determine if Two Values Are Nearly Equal$$
00032 
00033 $index NearEqual$$
00034 $index equal, near$$
00035 $index absolute, difference$$
00036 $index relative, difference$$
00037 $index difference, absolute$$
00038 $index difference, relative$$
00039 
00040 $head Syntax$$
00041 
00042 $code # include <cppad/near_equal.hpp>$$
00043 $pre
00044 $$
00045 $icode%b% = NearEqual(%x%, %y%, %r%, %a%)%$$
00046 
00047 
00048 $head Purpose$$
00049 Returns true,
00050 if $icode x$$ and $icode y$$ are nearly equal,
00051 and false otherwise.
00052 
00053 $head x$$
00054 The argument $icode x$$ 
00055 has one of the following possible prototypes
00056 $codei%
00057      const %Type%               &%x%,
00058      const std::complex<%Type%> &%x%, 
00059 %$$
00060 
00061 $head y$$
00062 The argument $icode y$$ 
00063 has one of the following possible prototypes
00064 $codei%
00065      const %Type%               &%y%,
00066      const std::complex<%Type%> &%y%, 
00067 %$$
00068 
00069 $head r$$
00070 The relative error criteria $icode r$$ has prototype
00071 $codei%
00072      const %Type% &%r%
00073 %$$
00074 It must be greater than or equal to zero.
00075 The relative error condition is defined as:
00076 $latex \[
00077      | x - y | \leq r ( |x| + |y| ) 
00078 \] $$
00079 
00080 $head a$$
00081 The absolute error criteria $icode a$$ has prototype
00082 $codei%
00083      const %Type% &%a%
00084 %$$
00085 It must be greater than or equal to zero.
00086 The absolute error condition is defined as:
00087 $latex \[
00088      | x - y | \leq a
00089 \] $$
00090 
00091 $head b$$
00092 The return value $icode b$$ has prototype
00093 $codei%
00094      bool %b%
00095 %$$
00096 If either $icode x$$ or $icode y$$ is infinite or not a number, 
00097 the return value is false.
00098 Otherwise, if either the relative or absolute error 
00099 condition (defined above) is satisfied, the return value is true.
00100 Otherwise, the return value is false.
00101 
00102 $head Type$$
00103 The type $icode Type$$ must be a
00104 $cref NumericType$$.
00105 The routine $cref CheckNumericType$$ will generate
00106 an error message if this is not the case.
00107 In addition, the following operations must be defined objects
00108 $icode a$$ and $icode b$$ of type $icode Type$$:
00109 $table
00110 $bold Operation$$     $cnext 
00111      $bold Description$$ $rnext
00112 $icode%a% <= %b%$$  $cnext 
00113      less that or equal operator (returns a $code bool$$ object)
00114 $tend
00115 
00116 $head Include Files$$
00117 The file $code cppad/near_equal.hpp$$ is included by $code cppad/cppad.hpp$$
00118 but it can also be included separately with out the rest of 
00119 the $code CppAD$$ routines.
00120 
00121 $head Example$$
00122 $children%
00123      example/near_equal.cpp
00124 %$$
00125 The file $cref near_equal.cpp$$ contains an example
00126 and test of $code NearEqual$$.
00127 It return true if it succeeds and false otherwise.
00128 
00129 $head Exercise$$
00130 $index exercise, NearEqual$$
00131 Create and run a program that contains the following code:
00132 $codep
00133      using std::complex;
00134      using std::cout;
00135      using std::endl;
00136 
00137      complex<double> one(1., 0), i(0., 1);
00138      complex<double> x = one / i;
00139      complex<double> y = - i;
00140      double          r = 1e-12;
00141      double          a = 0;
00142      bool           ok = CppAD::NearEqual(x, y, r, a);
00143      if( ok )
00144           cout << "Ok"    << endl;
00145      else cout << "Error" << endl;
00146 $$
00147 
00148 $end
00149 
00150 */
00151 
00152 # include <complex>
00153 # include <cppad/local/cppad_assert.hpp>
00154 # include <cppad/check_numeric_type.hpp>
00155 
00156 namespace CppAD { // Begin CppAD namespace
00157 
00158 // determine if both x and y are finite values (z1 and z2 are zero).
00159 template <class Type>
00160 bool near_equal_isfinite(
00161      const Type &z1, const Type &z2, const Type &x , const Type &y)
00162 {    Type infinity = Type(1) / z1;
00163      Type nan      = z1      / z2;
00164 
00165      // handle bug where some compilers return true for nan == nan
00166      bool xNan = ( x != x || x == nan );
00167      bool yNan = ( y != y || y == nan );
00168 
00169      // infinite cases
00170      bool xInf = (x == infinity   || x == - infinity);
00171      bool yInf = (x == infinity   || x == - infinity);
00172      
00173      return ! (xNan | yNan | xInf | yInf);
00174 }
00175 
00176 template <class Type>
00177 bool NearEqual(const Type &x, const Type &y, const Type &r, const Type &a)
00178 {
00179      CheckNumericType<Type>();
00180      Type zero(0);
00181 
00182      CPPAD_ASSERT_KNOWN(
00183           zero <= r,
00184           "Error in NearEqual: relative error is less than zero"
00185      );
00186      CPPAD_ASSERT_KNOWN(
00187           zero <= a,
00188           "Error in NearEqual: absolute error is less than zero"
00189      );
00190 
00191      // check for special cases
00192      if( ! CppAD::near_equal_isfinite(zero, zero, x, y) )
00193           return false;
00194 
00195      Type ax = x;
00196      if( ax <= zero )
00197           ax = - ax;
00198 
00199      Type ay = y;
00200      if( ay <= zero )
00201           ay = - ay;
00202 
00203      Type ad = x - y;
00204      if( ad <= zero )
00205           ad = - ad;
00206 
00207      if( ad <= a )
00208           return true;
00209 
00210      if( ad <= r * (ax + ay) )
00211           return true;
00212 
00213      return false;
00214 }
00215 
00216 template <class Type>
00217 bool NearEqual(
00218      const std::complex<Type> &x , 
00219      const std::complex<Type> &y , 
00220      const              Type  &r , 
00221      const              Type  & a )
00222 {
00223      CheckNumericType<Type>();
00224      Type zero(0);
00225 
00226      CPPAD_ASSERT_KNOWN(
00227           zero <= r,
00228           "Error in NearEqual: relative error is less than zero"
00229      );
00230      CPPAD_ASSERT_KNOWN(
00231           zero <= a,
00232           "Error in NearEqual: absolute error is less than zero"
00233      );
00234 
00235      // check for special cases
00236      if( ! CppAD::near_equal_isfinite(zero, zero, x.real(), x.imag()) )
00237           return false;
00238      if( ! CppAD::near_equal_isfinite(zero, zero, y.real(), y.imag()) )
00239           return false;
00240 
00241      std::complex<Type> d = x - y;
00242 
00243      Type ad = std::abs(d);
00244      if( ad <= a )
00245           return true;
00246 
00247      Type ax = std::abs(x);
00248      Type ay = std::abs(y);
00249      if( ad <= r * (ax + ay) )
00250           return true;
00251 
00252      return false;
00253 }
00254 
00255 template <class Type>
00256 bool NearEqual(
00257      const std::complex<Type> &x , 
00258      const              Type  &y , 
00259      const              Type  &r , 
00260      const              Type  & a )
00261 {
00262      return NearEqual(x, std::complex<Type>(y, Type(0)), r, a);
00263 }
00264 
00265 template <class Type>
00266 bool NearEqual(
00267      const              Type  &x , 
00268      const std::complex<Type> &y , 
00269      const              Type  &r , 
00270      const              Type  & a )
00271 {
00272      return NearEqual(std::complex<Type>(x, Type(0)), y, r, a);
00273 }
00274 
00275 } // END CppAD namespace
00276 
00277 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines