CppAD: A C++ Algorithmic Differentiation Package 20110419
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-07 Bradley M. Bell
00007 
00008 CppAD is distributed under multiple licenses. This distribution is under
00009 the terms of the 
00010                     Common 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 $syntax%%b% = NearEqual(%x%, %y%, %r%, %a%)%$$
00046 
00047 
00048 $head Purpose$$
00049 Returns true,
00050 if $italic x$$ and $italic y$$ are nearly equal,
00051 and false otherwise.
00052 
00053 $head x$$
00054 The argument $italic x$$ 
00055 has one of the following possible prototypes
00056 $syntax%
00057         const %Type%               &%x%,
00058         const std::complex<%Type%> &%x%, 
00059 %$$
00060 
00061 $head y$$
00062 The argument $italic y$$ 
00063 has one of the following possible prototypes
00064 $syntax%
00065         const %Type%               &%y%,
00066         const std::complex<%Type%> &%y%, 
00067 %$$
00068 
00069 $head r$$
00070 The relative error criteria $italic r$$ has prototype
00071 $syntax%
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 $italic a$$ has prototype
00082 $syntax%
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 $italic b$$ has prototype
00093 $syntax%
00094         bool %b%
00095 %$$
00096 If either $italic x$$ or $italic 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 $italic Type$$ must be a
00104 $xref/NumericType/$$.
00105 The routine $xref/CheckNumericType/$$ will generate
00106 an error message if this is not the case.
00107 In addition, the following operations must be defined objects
00108 $italic a$$ and $italic b$$ of type $italic Type$$:
00109 $table
00110 $bold Operation$$     $cnext 
00111         $bold Description$$ $rnext
00112 $syntax%%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 $xref/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