Support std::complex<T>.

14 files changed, 90 insertions, 12 deletions

diff --git a/src/vZ/Traits.hpp b/src/vZ/Traits.hpp
index 2c1496e..1c8f082 100644
--- a/src/vZ/Traits.hpp
+++ b/src/vZ/Traits.hpp
@@ -21,6 +21,8 @@
 #ifndef VZ_TRAITS_HPP
 #define VZ_TRAITS_HPP
 
+#include <complex>
+
 namespace vZ
 {
   // Traits class
@@ -36,6 +38,17 @@ namespace vZ
   private:
     Traits();
   };
+
+  // Specialization for std::complex<T>
+  template <typename T>
+  class Traits<std::complex<T> >
+  {
+  public:
+    typedef typename Traits<T>::Scalar Scalar;
+
+  private:
+    Traits();
+  };
 }
 
 #endif // VZ_TRAITS_HPP
diff --git a/tests/BS23.cpp b/tests/BS23.cpp
index cdd2b9a..c2eb345 100644
--- a/tests/BS23.cpp
+++ b/tests/BS23.cpp
@@ -51,7 +51,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 4.8e-6) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/CK45.cpp b/tests/CK45.cpp
index 62f4870..c040436 100644
--- a/tests/CK45.cpp
+++ b/tests/CK45.cpp
@@ -51,7 +51,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 2.7e-6) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/Complex.cpp b/tests/Complex.cpp
new file mode 100644
index 0000000..df517b0
--- /dev/null
+++ b/tests/Complex.cpp
@@ -0,0 +1,63 @@
+/*************************************************************************
+ * Copyright (C) 2010 Tavian Barnes <tavianator@gmail.com>               *
+ *                                                                       *
+ * This file is part of The vZ Test Suite.                               *
+ *                                                                       *
+ * The vZ Test Suite is free software; you can redistribute it and/or    *
+ * modify it under the terms of the GNU General Public License as        *
+ * published by the Free Software Foundation; either version 3 of the    *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The vZ Test Suite is distributed in the hope that it will be useful,  *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of        *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU     *
+ * General Public License for more details.                              *
+ *                                                                       *
+ * You should have received a copy of the GNU General Public License     *
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>. *
+ *************************************************************************/
+
+#include "vZ.hpp"
+#include <complex>
+#include <cmath>
+#include <cstdlib>
+#include <iostream>
+#include <iomanip>
+
+// y' = -y (y == C*exp(t))
+std::complex<double>
+f(double t, std::complex<double> y)
+{
+  return y;
+}
+
+int
+main()
+{
+  vZ::GenericDP45Integrator<std::complex<double> > integrator(f);
+  integrator.tol(1e-6)
+            .y(std::complex<double>(1.0, 0.0))
+            .x(0.0)
+            .h(0.06);
+
+  integrator.integrate(2.0);
+
+  std::complex<double> actual = integrator.y();
+  std::complex<double> expected(std::exp(2.0), 0.0);
+
+  std::cout << std::setprecision(10)
+            << "Numerical:  " << actual << std::endl
+            << "Expected:   " << expected  << std::endl
+            << "h:          " << integrator.h() << std::endl
+            << "Iterations: " << integrator.iterations() << std::endl
+            << "Rejections: " << integrator.rejections() << std::endl;
+
+  double error = std::abs(expected - actual)/std::abs(expected);
+  if (error > 6.0e-7) {
+    std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
+    return EXIT_FAILURE;
+  } else {
+    std::cout << "Error:      " << 100.0*error << "%" << std::endl;
+    return EXIT_SUCCESS;
+  }
+}
diff --git a/tests/DP45.cpp b/tests/DP45.cpp
index 5d6cb43..12c3202 100644
--- a/tests/DP45.cpp
+++ b/tests/DP45.cpp
@@ -51,7 +51,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 6.0e-7) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/EquationSystem.cpp b/tests/EquationSystem.cpp
index b1fed00..0a929aa 100644
--- a/tests/EquationSystem.cpp
+++ b/tests/EquationSystem.cpp
@@ -59,7 +59,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 6.0e-7) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/Euler.cpp b/tests/Euler.cpp
index d5de734..fa2d97c 100644
--- a/tests/Euler.cpp
+++ b/tests/Euler.cpp
@@ -48,7 +48,7 @@ main()
             << "Expected:   " << expected  << std::endl
             << "Iterations: " << integrator.iterations() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 0.01) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/HE12.cpp b/tests/HE12.cpp
index 20b37f2..fd2fb56 100644
--- a/tests/HE12.cpp
+++ b/tests/HE12.cpp
@@ -51,7 +51,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 8.7e-7) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/Heun.cpp b/tests/Heun.cpp
index 739a42f..d9e4cbe 100644
--- a/tests/Heun.cpp
+++ b/tests/Heun.cpp
@@ -48,7 +48,7 @@ main()
             << "Expected:   " << expected  << std::endl
             << "Iterations: " << integrator.iterations() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 1.4e-4) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/Makefile.am b/tests/Makefile.am
index 2091b15..52e3064 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -30,7 +30,8 @@ check_PROGRAMS = Euler-test                                                    \
                  DP45-test                                                     \
                  Vector-test                                                   \
                  EquationSystem-test                                           \
-                 EquationSystem-Vector-test
+                 EquationSystem-Vector-test                                    \
+                 Complex-test
 TESTS          = $(check_PROGRAMS)
 
 Euler_test_SOURCES                 = Euler.cpp
@@ -45,3 +46,4 @@ DP45_test_SOURCES                  = DP45.cpp
 Vector_test_SOURCES                = Vector.cpp
 EquationSystem_test_SOURCES        = EquationSystem.cpp
 EquationSystem_Vector_test_SOURCES = EquationSystem-Vector.cpp
+Complex_test_SOURCES               = Complex.cpp
diff --git a/tests/Midpoint.cpp b/tests/Midpoint.cpp
index 385adc4..f667ca5 100644
--- a/tests/Midpoint.cpp
+++ b/tests/Midpoint.cpp
@@ -48,7 +48,7 @@ main()
             << "Expected:   " << expected  << std::endl
             << "Iterations: " << integrator.iterations() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 1.4e-4) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/RK4.cpp b/tests/RK4.cpp
index 8e2f43e..75ec46a 100644
--- a/tests/RK4.cpp
+++ b/tests/RK4.cpp
@@ -48,7 +48,7 @@ main()
             << "Expected:   " << expected  << std::endl
             << "Iterations: " << integrator.iterations() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 4.2e-8) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/RKF45.cpp b/tests/RKF45.cpp
index 61f7949..59bad82 100644
--- a/tests/RKF45.cpp
+++ b/tests/RKF45.cpp
@@ -51,7 +51,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 1.7e-6) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;
diff --git a/tests/Vector.cpp b/tests/Vector.cpp
index e1ef4dd..e3b660c 100644
--- a/tests/Vector.cpp
+++ b/tests/Vector.cpp
@@ -51,7 +51,7 @@ main()
             << "Iterations: " << integrator.iterations() << std::endl
             << "Rejections: " << integrator.rejections() << std::endl;
 
-  double error = std::fabs(expected - actual)/expected;
+  double error = std::abs(expected - actual)/expected;
   if (error > 1.5e-6) {
     std::cerr << "Error:      " << 100.0*error << "%" << std::endl;
     return EXIT_FAILURE;