CppAD: A C++ Algorithmic Differentiation Package
20130102
|
00001 /* $Id$ */ 00002 # ifndef CPPAD_BASE_COMPLEX_INCLUDED 00003 # define CPPAD_BASE_COMPLEX_INCLUDED 00004 /* -------------------------------------------------------------------------- 00005 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-12 Bradley M. Bell 00006 00007 CppAD is distributed under multiple licenses. This distribution is under 00008 the terms of the 00009 Eclipse Public License Version 1.0. 00010 00011 A copy of this license is included in the COPYING file of this distribution. 00012 Please visit http://www.coin-or.org/CppAD/ for information on other licenses. 00013 -------------------------------------------------------------------------- */ 00014 # include <limits> 00015 # include <complex> 00016 00017 // needed before one can use CPPAD_ASSERT_FIRST_CALL_NOT_PARALLEL 00018 # include <cppad/thread_alloc.hpp> 00019 00020 /* 00021 $begin base_complex.hpp$$ 00022 $spell 00023 eps 00024 abs_geq 00025 Rel 00026 Lt Le Eq Ge Gt 00027 imag 00028 gcc 00029 isnan 00030 cppad.hpp 00031 sqrt 00032 exp 00033 cos 00034 std 00035 const 00036 CppAD 00037 Op 00038 inline 00039 enum 00040 undef 00041 acos 00042 asin 00043 atan 00044 erf 00045 Cond 00046 namespace 00047 bool 00048 $$ 00049 00050 $index complex, double Base$$ 00051 $index Base, double complex$$ 00052 $index double, complex Base$$ 00053 00054 $section Enable use of AD<Base> where Base is std::complex<double>$$ 00055 00056 $children% 00057 example/complex_poly.cpp% 00058 example/not_complex_ad.cpp 00059 %$$ 00060 00061 $head Example$$ 00062 The file $cref complex_poly.cpp$$ contains an example use of 00063 $code std::complex<double>$$ type for a CppAD $icode Base$$ type. 00064 It returns true if it succeeds and false otherwise. 00065 00066 $head See Also$$ 00067 The file $cref not_complex_ad.cpp$$ contains an example using 00068 complex arithmetic where the function is not complex differentiable. 00069 00070 $head Include Order$$ 00071 This file is included before $code <cppad/cppad.hpp>$$ 00072 so it is necessary to define the error handler 00073 in addition to including 00074 $cref/base_require.hpp/base_require/Include Order/$$ 00075 $codep */ 00076 # include <limits> 00077 # include <complex> 00078 # include <cppad/base_require.hpp> 00079 # include <cppad/local/cppad_assert.hpp> 00080 00081 /* $$ 00082 00083 $head CondExpOp$$ 00084 The type $code std::complex<double>$$ does not supports the 00085 $code <$$, $code <=$$, $code ==$$, $code >=$$, and $code >$$ operators; see 00086 $cref/not ordered/base_cond_exp/CondExpTemplate/Not Ordered/$$. 00087 Hence its $code CondExpOp$$ function is defined by 00088 $codep */ 00089 namespace CppAD { 00090 inline std::complex<double> CondExpOp( 00091 enum CppAD::CompareOp cop , 00092 const std::complex<double> &left , 00093 const std::complex<double> &right , 00094 const std::complex<double> &trueCase , 00095 const std::complex<double> &falseCase ) 00096 { CppAD::ErrorHandler::Call( 00097 true , __LINE__ , __FILE__ , 00098 "std::complex<float> CondExpOp(...)", 00099 "Error: cannot use CondExp with a complex type" 00100 ); 00101 return std::complex<double>(0); 00102 } 00103 } 00104 /* $$ 00105 00106 $head CondExpRel$$ 00107 The $cref/CPPAD_COND_EXP_REL/base_cond_exp/CondExpRel/$$ macro invocation 00108 $codep */ 00109 namespace CppAD { 00110 CPPAD_COND_EXP_REL( std::complex<double> ) 00111 } 00112 /* $$ 00113 used $code CondExpOp$$ above to 00114 define $codei%CondExp%Rel%$$ for $code std::complex<double>$$ arguments 00115 and $icode%Rel%$$ equal to 00116 $code Lt$$, $code Le$$, $code Eq$$, $code Ge$$, and $code Gt$$. 00117 00118 $head EqualOpSeq$$ 00119 Complex numbers do not carry operation sequence information. 00120 Thus they are equal in this sense if and only if there values are equal. 00121 $codep */ 00122 namespace CppAD { 00123 inline bool EqualOpSeq( 00124 const std::complex<double> &x , 00125 const std::complex<double> &y ) 00126 { return x == y; 00127 } 00128 } 00129 /* $$ 00130 00131 $head Identical$$ 00132 Complex numbers do not carry operation sequence information. 00133 Thus they are all parameters so the identical functions just check values. 00134 $codep */ 00135 namespace CppAD { 00136 inline bool IdenticalPar(const std::complex<double> &x) 00137 { return true; } 00138 inline bool IdenticalZero(const std::complex<double> &x) 00139 { return (x == std::complex<double>(0., 0.) ); } 00140 inline bool IdenticalOne(const std::complex<double> &x) 00141 { return (x == std::complex<double>(1., 0.) ); } 00142 inline bool IdenticalEqualPar( 00143 const std::complex<double> &x, const std::complex<double> &y) 00144 { return (x == y); } 00145 } 00146 /* $$ 00147 00148 $head Ordered$$ 00149 Complex types do not support comparison operators, 00150 $codep */ 00151 # undef CPPAD_USER_MACRO 00152 # define CPPAD_USER_MACRO(Fun) \ 00153 inline bool Fun(const std::complex<double>& x) \ 00154 { CppAD::ErrorHandler::Call( \ 00155 true , __LINE__ , __FILE__ , \ 00156 #Fun"(x)", \ 00157 "Error: cannot use " #Fun " with x complex<double> " \ 00158 ); \ 00159 return false; \ 00160 } 00161 namespace CppAD { 00162 CPPAD_USER_MACRO(LessThanZero) 00163 CPPAD_USER_MACRO(LessThanOrZero) 00164 CPPAD_USER_MACRO(GreaterThanOrZero) 00165 CPPAD_USER_MACRO(GreaterThanZero) 00166 inline bool abs_geq( 00167 const std::complex<double>& x , 00168 const std::complex<double>& y ) 00169 { return std::abs(x) >= std::abs(y); } 00170 } 00171 /* $$ 00172 00173 $head Integer$$ 00174 The implementation of this function must agree 00175 with the CppAD user specifications for complex arguments to the 00176 $cref/Integer/Integer/x/Complex Types/$$ function: 00177 $codep */ 00178 namespace CppAD { 00179 inline int Integer(const std::complex<double> &x) 00180 { return static_cast<int>( x.real() ); } 00181 } 00182 /* $$ 00183 00184 $head isnan$$ 00185 The gcc 4.1.1 complier defines the function 00186 $codei% 00187 int std::complex<double>::isnan( std::complex<double> %z% ) 00188 %$$ 00189 (which is not specified in the C++ 1998 standard ISO/IEC 14882). 00190 This causes an ambiguity between the function above and the CppAD 00191 $cref/isnan/nan/$$ template function. 00192 We avoid this ambiguity by defining a non-template version of 00193 this function in the CppAD namespace. 00194 $codep */ 00195 namespace CppAD { 00196 inline bool isnan(const std::complex<double>& z) 00197 { CPPAD_ASSERT_FIRST_CALL_NOT_PARALLEL; 00198 static const double nan = std::numeric_limits<double>::quiet_NaN(); 00199 return (z != z) | (z.real() == nan) | (z.imag() == nan); 00200 } 00201 } 00202 /* $$ 00203 00204 $head Valid Unary Math$$ 00205 The following macro invocations define the standard unary 00206 math functions that are valid with complex arguments and are 00207 required to use $code AD< std::complex<double> >$$. 00208 $codep */ 00209 namespace CppAD { 00210 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, cos) 00211 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, cosh) 00212 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, exp) 00213 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, log) 00214 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, sin) 00215 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, sinh) 00216 CPPAD_STANDARD_MATH_UNARY(std::complex<double>, sqrt) 00217 } 00218 /* $$ 00219 00220 $head Invalid Unary Math$$ 00221 The following macro definition and invocations define the standard unary 00222 math functions that are invalid with complex arguments and are 00223 required to use $code AD< std::complex<double> >$$. 00224 $codep */ 00225 # undef CPPAD_USER_MACRO 00226 # define CPPAD_USER_MACRO(Fun) \ 00227 inline std::complex<double> Fun(const std::complex<double>& x) \ 00228 { CppAD::ErrorHandler::Call( \ 00229 true , __LINE__ , __FILE__ , \ 00230 #Fun"(x)", \ 00231 "Error: cannot use " #Fun " with x complex<double> " \ 00232 ); \ 00233 return std::complex<double>(0); \ 00234 } 00235 namespace CppAD { 00236 CPPAD_USER_MACRO(abs) 00237 CPPAD_USER_MACRO(acos) 00238 CPPAD_USER_MACRO(asin) 00239 CPPAD_USER_MACRO(atan) 00240 CPPAD_USER_MACRO(sign) 00241 } 00242 /* $$ 00243 00244 $head pow $$ 00245 The following defines a $code CppAD::pow$$ function that 00246 is required to use $code AD< std::complex<double> >$$: 00247 $codep */ 00248 namespace CppAD { 00249 inline std::complex<double> pow( 00250 const std::complex<double> &x , 00251 const std::complex<double> &y ) 00252 { return std::pow(x, y); } 00253 } 00254 /*$$ 00255 00256 $head limits$$ 00257 The following defines the numeric limits functions 00258 $code epsilon$$, $code min$$, and $code max$$ for the type 00259 $code std::complex<double>$$. 00260 It also defines the deprecated $code epsilon$$ function: 00261 $codep */ 00262 namespace CppAD { 00263 template <> 00264 class numeric_limits< std::complex<double> > { 00265 public: 00266 // machine epsilon 00267 static std::complex<double> epsilon(void) 00268 { double eps = std::numeric_limits<double>::epsilon(); 00269 return std::complex<double>(eps, 0.0); 00270 } 00271 // minimum positive normalized value 00272 static std::complex<double> min(void) 00273 { double min = std::numeric_limits<double>::min(); 00274 return std::complex<double>(min, 0.0); 00275 } 00276 // maximum finite value 00277 static std::complex<double> max(void) 00278 { double max = std::numeric_limits<double>::max(); 00279 return std::complex<double>(max, 0.0); 00280 } 00281 }; 00282 // deprecated machine epsilon 00283 template <> 00284 inline std::complex<double> epsilon< std::complex<double> > (void) 00285 { return numeric_limits< std::complex<double> >::epsilon(); } 00286 } 00287 /* $$ 00288 $end 00289 */ 00290 # undef CPPAD_USER_MACRO_ONE 00291 # define CPPAD_USER_MACRO_ONE(Fun) \ 00292 inline bool Fun(const std::complex<float>& x) \ 00293 { CppAD::ErrorHandler::Call( \ 00294 true , __LINE__ , __FILE__ , \ 00295 #Fun"(x)", \ 00296 "Error: cannot use " #Fun " with x complex<float> " \ 00297 ); \ 00298 return false; \ 00299 } 00300 # undef CPPAD_USER_MACRO_TWO 00301 # define CPPAD_USER_MACRO_TWO(Fun) \ 00302 inline std::complex<float> Fun(const std::complex<float>& x) \ 00303 { CppAD::ErrorHandler::Call( \ 00304 true , __LINE__ , __FILE__ , \ 00305 #Fun"(x)", \ 00306 "Error: cannot use " #Fun " with x complex<float> " \ 00307 ); \ 00308 return std::complex<float>(0); \ 00309 } 00310 namespace CppAD { 00311 // CondExpOp ------------------------------------------------------ 00312 inline std::complex<float> CondExpOp( 00313 enum CppAD::CompareOp cop , 00314 const std::complex<float> &left , 00315 const std::complex<float> &right , 00316 const std::complex<float> &trueCase , 00317 const std::complex<float> &falseCase ) 00318 { CppAD::ErrorHandler::Call( 00319 true , __LINE__ , __FILE__ , 00320 "std::complex<float> CondExpOp(...)", 00321 "Error: cannot use CondExp with a complex type" 00322 ); 00323 return std::complex<float>(0); 00324 } 00325 // CondExpRel -------------------------------------------------------- 00326 CPPAD_COND_EXP_REL( std::complex<float> ) 00327 // EqualOpSeq ----------------------------------------------------- 00328 inline bool EqualOpSeq( 00329 const std::complex<float> &x , 00330 const std::complex<float> &y ) 00331 { return x == y; 00332 } 00333 // Identical ------------------------------------------------------ 00334 inline bool IdenticalPar(const std::complex<float> &x) 00335 { return true; } 00336 inline bool IdenticalZero(const std::complex<float> &x) 00337 { return (x == std::complex<float>(0., 0.) ); } 00338 inline bool IdenticalOne(const std::complex<float> &x) 00339 { return (x == std::complex<float>(1., 0.) ); } 00340 inline bool IdenticalEqualPar( 00341 const std::complex<float> &x, const std::complex<float> &y) 00342 { return (x == y); } 00343 // Ordered -------------------------------------------------------- 00344 CPPAD_USER_MACRO_ONE(LessThanZero) 00345 CPPAD_USER_MACRO_ONE(LessThanOrZero) 00346 CPPAD_USER_MACRO_ONE(GreaterThanOrZero) 00347 CPPAD_USER_MACRO_ONE(GreaterThanZero) 00348 inline bool abs_geq( 00349 const std::complex<float>& x , 00350 const std::complex<float>& y ) 00351 { return std::abs(x) >= std::abs(y); } 00352 // Integer ------------------------------------------------------ 00353 inline int Integer(const std::complex<float> &x) 00354 { return static_cast<int>( x.real() ); } 00355 // isnan ------------------------------------------------------------- 00356 inline bool isnan(const std::complex<float>& z) 00357 { CPPAD_ASSERT_FIRST_CALL_NOT_PARALLEL; 00358 static const float nan = std::numeric_limits<float>::quiet_NaN(); 00359 return (z != z) | (z.real() == nan) | (z.imag() == nan); 00360 } 00361 // Valid standard math functions -------------------------------- 00362 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, cos) 00363 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, cosh) 00364 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, exp) 00365 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, log) 00366 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, sin) 00367 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, sinh) 00368 CPPAD_STANDARD_MATH_UNARY(std::complex<float>, sqrt) 00369 // Invalid standrd math functions ------------------------------- 00370 CPPAD_USER_MACRO_TWO(abs) 00371 CPPAD_USER_MACRO_TWO(acos) 00372 CPPAD_USER_MACRO_TWO(asin) 00373 CPPAD_USER_MACRO_TWO(atan) 00374 CPPAD_USER_MACRO_TWO(sign) 00375 // The pow function 00376 inline std::complex<float> pow( 00377 const std::complex<float> &x , 00378 const std::complex<float> &y ) 00379 { return std::pow(x, y); } 00380 // numeric_limits ------------------------------------------------- 00381 template <> 00382 class numeric_limits< std::complex<float> > { 00383 public: 00384 /// machine epsilon 00385 static std::complex<float> epsilon(void) 00386 { float eps = std::numeric_limits<float>::epsilon(); 00387 return std::complex<float>(eps, 0.0); 00388 } 00389 /// minimum positive normalized value 00390 static std::complex<float> min(void) 00391 { float min = std::numeric_limits<float>::min(); 00392 return std::complex<float>(min, 0.0); 00393 } 00394 /// maximum finite value 00395 static std::complex<float> max(void) 00396 { float max = std::numeric_limits<float>::max(); 00397 return std::complex<float>(max, 0.0); 00398 } 00399 }; 00400 template <> 00401 inline std::complex<float> epsilon< std::complex<float> >(void) 00402 { return numeric_limits< std::complex<float> >::epsilon(); } 00403 } 00404 00405 // undefine macros only used by this file 00406 # undef CPPAD_USER_MACRO 00407 # undef CPPAD_USER_MACRO_ONE 00408 # undef CPPAD_USER_MACRO_TWO 00409 00410 # endif