lots of fixes for complex numbers. unittests now test complex numbers.

tvmet-1.7.1/include/tvmet/Traits.h

@@ -56,6 +56,7 @@ struct Traits : public TraitsBase<T>
   using Base::imag;
   using Base::conj;
   using Base::sqrt;
+  using Base::isLessThan_nonfuzzy;
   
   static value_type random()
   {
@@ -74,7 +75,7 @@ struct Traits : public TraitsBase<T>
     if(isFloat())
       return(abs(a) <= abs(b) * epsilon());
     else
-      return(a==0);
+      return(a==static_cast<value_type>(0));
   }
   
   /**
@@ -105,13 +106,12 @@ struct Traits : public TraitsBase<T>
     *
     * Full story: returns a <= b || isApprox(a, b);
     */
-  static bool isLessThan( argument_type a, argument_type b )
+  static bool isLessThan(argument_type a, argument_type b)
   {
-    assert(!isComplex());
     if(isFloat())
-      return(a <= b || isApprox(a, b));
+      return(isLessThan_nonfuzzy(a,b) || isApprox(a, b));
     else
-      return(a<=b);
+      return(isLessThan_nonfuzzy(a,b));
   }
 
 };

tvmet-1.7.1/testsuite/main.h

@@ -23,6 +23,8 @@
 #include <stdlib.h>
 #include <time.h>
 
+#define EIGEN_USE_COMPLEX
+
 #ifdef EIGEN_USE_COMPLEX
 #include <complex>
 #endif

tvmet-1.7.1/testsuite/testtraits.cpp

@@ -25,12 +25,26 @@
 
 template<typename T> struct TestTraits
 {
-  void real()
+  void real_imag_conj_abs()
   {
     T x = Traits<T>::random();
-    typedef typename Traits<T>::real_type real_type;
+    typedef typename Traits<T>::real_type real;
-    real_type r = Traits<T>::real(x);
+    real r = Traits<T>::real(x);
-    TEST_APPROX(r, x);
+    real i = Traits<T>::imag(x);
+    T c = Traits<T>::conj(x);
+    real a = Traits<T>::abs(x);
+    
+    // a must be real
+    TEST_APPROX(a, Traits<real>::real(a));
+    TEST_APPROX(a, Traits<T>::real(a));
+    TEST_ZERO(Traits<real>::imag(a));
+    TEST_ZERO(Traits<T>::imag(a));
+    
+    // check Pythagora's formula
+    if(Traits<T>::isFloat() || !Traits<T>::isComplex()) TEST_APPROX(r*r + i*i, a*a);
+    
+    // check complex conjugation
+    TEST_APPROX(-i, Traits<T>::imag(c));
   }
 
   void imag()
@@ -60,20 +74,18 @@ template<typename T> struct TestTraits
   
   void sqrt()
   {
+    // only test compilation here
     T x = Traits<T>::random();
-    T a = Traits<T>::abs(x);
+    Traits<T>::sqrt(x);
-    T b = Traits<T>::sqrt(a);
-    // T could be an integer type, so b*b=a is not necessarily true
-    TEST_LESSTHAN(b*b, a);
-    TEST_LESSTHAN(a, (b+1)*(b+1));
   }
   
   void isApprox()
   {
     T x = Traits<T>::random();
+    T e = T(Traits<T>::epsilon()) / T(10);
     TEST(Traits<T>::isApprox(x,x));
-    TEST(Traits<T>::isApprox(x,x+Traits<T>::epsilon()/10));
+    TEST(Traits<T>::isApprox(x,x+e));
-    TEST(!Traits<T>::isApprox(x,x+1));
+    TEST(!Traits<T>::isApprox(x,x+T(1)));
   }
   
   void isNegligible()
@@ -88,6 +100,7 @@ template<typename T> struct TestTraits
     TEST(!Traits<T>::isNegligible(one, x));
   }
   
+  
   void isZero()
   {
     T zero(0), one(1), x = Traits<T>::random(), y = Traits<T>::random();
@@ -98,18 +111,22 @@ template<typename T> struct TestTraits
   
   void isLessThan()
   {
+    if(Traits<T>::isComplex()) {
+      T x = Traits<T>::random(), y = Traits<T>::random();
+      TEST(!Traits<T>::isLessThan(x,y));
+    }
+    else {
       T one(1), x = Traits<T>::random();
+      T e = T(Traits<T>::epsilon()) / T(10);
       TEST(Traits<T>::isLessThan(x, x+one));
       TEST(!Traits<T>::isLessThan(x+one, x));
-    TEST(Traits<T>::isLessThan(x, x+Traits<T>::epsilon()/10));
+      TEST(Traits<T>::isLessThan(x, x+e));
+    }
   }
   
   TestTraits()
   {
-    real();
+    real_imag_conj_abs();
-    imag();
-    conj();
-    abs();
     sqrt();
     isApprox();
     isNegligible();
@@ -123,4 +140,7 @@ void TvmetTestSuite::testTraits()
   TestTraits<int>();
   TestTraits<float>();
   TestTraits<double>();
+  TestTraits<std::complex<int> >();
+  TestTraits<std::complex<float> >();
+  TestTraits<std::complex<double> >();
 }