CppAD: A C++ Algorithmic Differentiation Package 20110419
|
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