Resolve two old TODOs in jet.h
1. Use using directive to pull functions from the std:: namespace
instead of defining a forwarding function.
2. Remove Eigen 2.0 related definitions.
Change-Id: If5f24ef740c17bc4db300fa04f4a2e8f809970d2
diff --git a/include/ceres/jet.h b/include/ceres/jet.h
index 2106357..3c53088 100644
--- a/include/ceres/jet.h
+++ b/include/ceres/jet.h
@@ -435,35 +435,33 @@
// This is necessary because we want to use the same name (e.g. 'sqrt') for
// double-valued and Jet-valued functions, but we are not allowed to put
// Jet-valued functions inside namespace std.
-//
-// TODO(keir): Switch to "using".
-inline double abs (double x) { return std::abs(x); }
-inline double log (double x) { return std::log(x); }
-inline double exp (double x) { return std::exp(x); }
-inline double sqrt (double x) { return std::sqrt(x); }
-inline double cos (double x) { return std::cos(x); }
-inline double acos (double x) { return std::acos(x); }
-inline double sin (double x) { return std::sin(x); }
-inline double asin (double x) { return std::asin(x); }
-inline double tan (double x) { return std::tan(x); }
-inline double atan (double x) { return std::atan(x); }
-inline double sinh (double x) { return std::sinh(x); }
-inline double cosh (double x) { return std::cosh(x); }
-inline double tanh (double x) { return std::tanh(x); }
-inline double floor (double x) { return std::floor(x); }
-inline double ceil (double x) { return std::ceil(x); }
-inline double pow (double x, double y) { return std::pow(x, y); }
-inline double atan2(double y, double x) { return std::atan2(y, x); }
-inline double cbrt (double x) { return std::cbrt(x); }
-inline double exp2 (double x) { return std::exp2(x); }
-inline double log2 (double x) { return std::log2(x); }
-inline double hypot(double x, double y) { return std::hypot(x, y); }
-inline double fmax(double x, double y) { return std::fmax(x, y); }
-inline double fmin(double x, double y) { return std::fmin(x, y); }
-inline double isfinite(double x) { return std::isfinite(x); }
-inline double isinf(double x) { return std::isinf(x); }
-inline double isnan(double x) { return std::isnan(x); }
-inline double isnormal(double x) { return std::isnormal(x); }
+using std::abs;
+using std::acos;
+using std::asin;
+using std::atan;
+using std::atan2;
+using std::cbrt;
+using std::ceil;
+using std::cos;
+using std::cosh;
+using std::exp;
+using std::exp2;
+using std::floor;
+using std::fmax;
+using std::fmin;
+using std::hypot;
+using std::isfinite;
+using std::isinf;
+using std::isnan;
+using std::isnormal;
+using std::log;
+using std::log2;
+using std::pow;
+using std::sin;
+using std::sinh;
+using std::sqrt;
+using std::tan;
+using std::tanh;
// Legacy names from pre-C++11 days.
inline bool IsFinite (double x) { return std::isfinite(x); }
@@ -899,38 +897,6 @@
return Jet<T, N>(tmp1, tmp2 * f.v + tmp3 * g.v);
}
-// Define the helper functions Eigen needs to embed Jet types.
-//
-// NOTE(keir): machine_epsilon() and precision() are missing, because they don't
-// work with nested template types (e.g. where the scalar is itself templated).
-// Among other things, this means that decompositions of Jet's does not work,
-// for example
-//
-// Matrix<Jet<T, N> ... > A, x, b;
-// ...
-// A.solve(b, &x)
-//
-// does not work and will fail with a strange compiler error.
-//
-// TODO(keir): This is an Eigen 2.0 limitation that is lifted in 3.0. When we
-// switch to 3.0, also add the rest of the specialization functionality.
-template<typename T, int N> inline const Jet<T, N>& ei_conj(const Jet<T, N>& x) { return x; } // NOLINT
-template<typename T, int N> inline const Jet<T, N>& ei_real(const Jet<T, N>& x) { return x; } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_imag(const Jet<T, N>& ) { return Jet<T, N>(0.0); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_abs (const Jet<T, N>& x) { return fabs(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_abs2(const Jet<T, N>& x) { return x * x; } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_sqrt(const Jet<T, N>& x) { return sqrt(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_exp (const Jet<T, N>& x) { return exp(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_log (const Jet<T, N>& x) { return log(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_sin (const Jet<T, N>& x) { return sin(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_cos (const Jet<T, N>& x) { return cos(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_tan (const Jet<T, N>& x) { return tan(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_atan(const Jet<T, N>& x) { return atan(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_sinh(const Jet<T, N>& x) { return sinh(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_cosh(const Jet<T, N>& x) { return cosh(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_tanh(const Jet<T, N>& x) { return tanh(x); } // NOLINT
-template<typename T, int N> inline Jet<T, N> ei_pow (const Jet<T, N>& x, Jet<T, N> y) { return pow(x, y); } // NOLINT
-
// Note: This has to be in the ceres namespace for argument dependent lookup to
// function correctly. Otherwise statements like CHECK_LE(x, 2.0) fail with
// strange compile errors.
diff --git a/internal/ceres/jet_test.cc b/internal/ceres/jet_test.cc
index e15fd67..3c6e23c 100644
--- a/internal/ceres/jet_test.cc
+++ b/internal/ceres/jet_test.cc
@@ -30,13 +30,14 @@
#include "ceres/jet.h"
+#include <Eigen/Dense>
#include <algorithm>
#include <cmath>
-#include "glog/logging.h"
-#include "gtest/gtest.h"
#include "ceres/stringprintf.h"
#include "ceres/test_util.h"
+#include "glog/logging.h"
+#include "gtest/gtest.h"
#define VL VLOG(1)
@@ -760,6 +761,61 @@
EXPECT_FALSE(IsNaN(a));
}
+// The following test ensures that Jets have all the appropriate Eigen
+// related traits so that they can be used as part of matrix
+// decompositions.
+TEST(Jet, FullRankEigenLLTSolve) {
+ Eigen::Matrix<J, 3, 3> A;
+ Eigen::Matrix<J, 3, 1> b, x;
+ for (int i = 0; i < 3; ++i) {
+ for (int j = 0; j < 3; ++j) {
+ A(i,j) = MakeJet(0.0, i, j * j);
+ }
+ b(i) = MakeJet(i, i, i);
+ x(i) = MakeJet(0.0, 0.0, 0.0);
+ A(i,i) = MakeJet(1.0, i, i * i);
+ }
+ x = A.llt().solve(b);
+ for (int i = 0; i < 3; ++i) {
+ EXPECT_EQ(x(i).a, b(i).a);
+ }
+}
+
+TEST(Jet, FullRankEigenLDLTSolve) {
+ Eigen::Matrix<J, 3, 3> A;
+ Eigen::Matrix<J, 3, 1> b, x;
+ for (int i = 0; i < 3; ++i) {
+ for (int j = 0; j < 3; ++j) {
+ A(i,j) = MakeJet(0.0, i, j * j);
+ }
+ b(i) = MakeJet(i, i, i);
+ x(i) = MakeJet(0.0, 0.0, 0.0);
+ A(i,i) = MakeJet(1.0, i, i * i);
+ }
+ x = A.ldlt().solve(b);
+ for (int i = 0; i < 3; ++i) {
+ EXPECT_EQ(x(i).a, b(i).a);
+ }
+}
+
+TEST(Jet, FullRankEigenLUSolve) {
+ Eigen::Matrix<J, 3, 3> A;
+ Eigen::Matrix<J, 3, 1> b, x;
+ for (int i = 0; i < 3; ++i) {
+ for (int j = 0; j < 3; ++j) {
+ A(i,j) = MakeJet(0.0, i, j * j);
+ }
+ b(i) = MakeJet(i, i, i);
+ x(i) = MakeJet(0.0, 0.0, 0.0);
+ A(i,i) = MakeJet(1.0, i, i * i);
+ }
+
+ x = A.lu().solve(b);
+ for (int i = 0; i < 3; ++i) {
+ EXPECT_EQ(x(i).a, b(i).a);
+ }
+}
+
// ScalarBinaryOpTraits is only supported on Eigen versions >= 3.3
#if EIGEN_VERSION_AT_LEAST(3, 3, 0)
TEST(JetTraitsTest, MatrixScalarUnaryOps) {