add a proof of concept autodiff jacobian helper class based on adolc

with unit test and FindAdolc cmake module
This commit is contained in:
Gael Guennebaud 2009-02-27 16:19:13 +00:00
parent 170128770a
commit 40774c625e
8 changed files with 341 additions and 91 deletions

View File

@ -92,13 +92,13 @@ set(INCLUDE_INSTALL_DIR
add_subdirectory(Eigen)
add_subdirectory(doc)
if(EIGEN_BUILD_TESTS)
include(CTest)
add_subdirectory(test)
endif(EIGEN_BUILD_TESTS)
add_subdirectory(doc)
add_subdirectory(unsupported)
if(EIGEN_BUILD_DEMOS)
@ -108,3 +108,7 @@ endif(EIGEN_BUILD_DEMOS)
if(EIGEN_BUILD_BTL)
add_subdirectory(bench/btl)
endif(EIGEN_BUILD_BTL)
if(EIGEN_BUILD_TESTS)
ei_testing_print_summary()
endif(EIGEN_BUILD_TESTS)

View File

@ -10,6 +10,11 @@ macro(ei_add_target_property target prop value)
endmacro(ei_add_target_property)
macro(ei_add_property prop value)
get_property(previous GLOBAL PROPERTY ${prop})
set_property(GLOBAL PROPERTY ${prop} "${previous}${value}")
endmacro(ei_add_property)
# Macro to add a test
#
# the unique parameter testname must correspond to a file
@ -73,7 +78,90 @@ macro(ei_add_test testname)
if(WIN32)
add_test(${testname} "${targetname}")
else(WIN32)
add_test(${testname} "${CMAKE_CURRENT_SOURCE_DIR}/runtest.sh" "${testname}")
add_test(${testname} "${Eigen_SOURCE_DIR}/test/runtest.sh" "${testname}")
endif(WIN32)
endmacro(ei_add_test)
# print a summary of the different options
macro(ei_testing_print_summary)
message("************************************************************")
message("*** Eigen's unit tests configuration summary ***")
message("************************************************************")
get_property(EIGEN_TESTING_SUMMARY GLOBAL PROPERTY EIGEN_TESTING_SUMMARY)
get_property(EIGEN_TESTED_BACKENDS GLOBAL PROPERTY EIGEN_TESTED_BACKENDS)
get_property(EIGEN_MISSING_BACKENDS GLOBAL PROPERTY EIGEN_MISSING_BACKENDS)
message("Enabled backends: ${EIGEN_TESTED_BACKENDS}")
message("Disabled backends: ${EIGEN_MISSING_BACKENDS}")
if(EIGEN_TEST_SSE2)
message("SSE2: ON")
else(EIGEN_TEST_SSE2)
message("SSE2: AUTO")
endif(EIGEN_TEST_SSE2)
if(EIGEN_TEST_SSE3)
message("SSE3: ON")
else(EIGEN_TEST_SSE3)
message("SSE3: AUTO")
endif(EIGEN_TEST_SSE3)
if(EIGEN_TEST_SSSE3)
message("SSSE3: ON")
else(EIGEN_TEST_SSSE3)
message("SSSE3: AUTO")
endif(EIGEN_TEST_SSSE3)
if(EIGEN_TEST_ALTIVEC)
message("Altivec: ON")
else(EIGEN_TEST_ALTIVEC)
message("Altivec: AUTO")
endif(EIGEN_TEST_ALTIVEC)
if(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("Explicit vec: OFF")
else(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("Explicit vec: AUTO")
endif(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("\n${EIGEN_TESTING_SUMMARY}")
# message("CXX: ${CMAKE_CXX_COMPILER}")
# if(CMAKE_COMPILER_IS_GNUCXX)
# execute_process(COMMAND ${CMAKE_CXX_COMPILER} --version COMMAND head -n 1 OUTPUT_VARIABLE EIGEN_CXX_VERSION_STRING OUTPUT_STRIP_TRAILING_WHITESPACE)
# message("CXX_VERSION: ${EIGEN_CXX_VERSION_STRING}")
# endif(CMAKE_COMPILER_IS_GNUCXX)
# message("CXX_FLAGS: ${CMAKE_CXX_FLAGS}")
# message("Sparse lib flags: ${SPARSE_LIBS}")
message("************************************************************")
endmacro(ei_testing_print_summary)
macro(ei_init_testing)
define_property(GLOBAL PROPERTY EIGEN_TESTED_BACKENDS BRIEF_DOCS " " FULL_DOCS " ")
define_property(GLOBAL PROPERTY EIGEN_MISSING_BACKENDS BRIEF_DOCS " " FULL_DOCS " ")
define_property(GLOBAL PROPERTY EIGEN_TESTING_SUMMARY BRIEF_DOCS " " FULL_DOCS " ")
set_property(GLOBAL PROPERTY EIGEN_TESTED_BACKENDS "")
set_property(GLOBAL PROPERTY EIGEN_MISSING_BACKENDS "")
set_property(GLOBAL PROPERTY EIGEN_TESTING_SUMMARY "")
endmacro(ei_init_testing)
if(CMAKE_COMPILER_IS_GNUCXX)
if(CMAKE_SYSTEM_NAME MATCHES Linux)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -g2")
set(CMAKE_CXX_FLAGS_RELWITHDEBINFO "${CMAKE_CXX_FLAGS_RELWITHDEBINFO} -O2 -g2")
set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} -fno-inline-functions")
set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -O0 -g2")
endif(CMAKE_SYSTEM_NAME MATCHES Linux)
set(EI_OFLAG "-O2")
# MSVC fails with:
# cl : Command line warning D9025 : overriding '/Od' with '/O2'
# cl : Command line error D8016 : '/RTC1' and '/O2' command-line options are incompatible
# elseif(MSVC)
# set(EI_OFLAG "/O2")
else(CMAKE_COMPILER_IS_GNUCXX)
set(EI_OFLAG "")
endif(CMAKE_COMPILER_IS_GNUCXX)

20
cmake/FindAdolc.cmake Normal file
View File

@ -0,0 +1,20 @@
if (ADOLC_INCLUDES AND ADOLC_LIBRARIES)
set(ADOLC_FIND_QUIETLY TRUE)
endif (ADOLC_INCLUDES AND ADOLC_LIBRARIES)
find_path(ADOLC_INCLUDES
NAMES
adolc/adouble.h
PATHS
$ENV{ADOLCDIR}
${INCLUDE_INSTALL_DIR}
)
find_library(ADOLC_LIBRARIES adolc PATHS $ENV{ADOLCDIR} ${LIB_INSTALL_DIR})
include(FindPackageHandleStandardArgs)
find_package_handle_standard_args(ADOLC DEFAULT_MSG
ADOLC_INCLUDES ADOLC_LIBRARIES)
mark_as_advanced(ADOLC_INCLUDES ADOLC_LIBRARIES)

View File

@ -1,4 +1,8 @@
include(EigenTesting)
enable_testing()
ei_init_testing()
find_package(GSL)
if(GSL_FOUND AND GSL_VERSION_MINOR LESS 9)
set(GSL_FOUND "")
@ -6,22 +10,22 @@ endif(GSL_FOUND AND GSL_VERSION_MINOR LESS 9)
if(GSL_FOUND)
add_definitions("-DHAS_GSL" ${GSL_DEFINITIONS})
include_directories(${GSL_INCLUDE_DIR})
ei_add_property(EIGEN_TESTED_BACKENDS "GSL, ")
else(GSL_FOUND)
ei_add_property(EIGEN_MISSING_BACKENDS "GSL, ")
set(GSL_LIBRARIES " ")
endif(GSL_FOUND)
set(SPARSE_LIBS "")
set(EIGEN_TESTED_BACKENDS "")
set(EIGEN_MISSING_BACKENDS "")
find_package(Taucs)
if(TAUCS_FOUND)
add_definitions("-DEIGEN_TAUCS_SUPPORT")
include_directories(${TAUCS_INCLUDES})
set(SPARSE_LIBS ${SPARSE_LIBS} ${TAUCS_LIBRARIES})
set(EIGEN_TESTED_BACKENDS ${EIGEN_TESTED_BACKENDS} Taucs)
ei_add_property(EIGEN_TESTED_BACKENDS "Taucs, ")
else(TAUCS_FOUND)
set(EIGEN_MISSING_BACKENDS ${EIGEN_MISSING_BACKENDS} Taucs)
ei_add_property(EIGEN_MISSING_BACKENDS "Taucs, ")
endif(TAUCS_FOUND)
find_package(Cholmod)
@ -29,9 +33,9 @@ if(CHOLMOD_FOUND)
add_definitions("-DEIGEN_CHOLMOD_SUPPORT")
include_directories(${CHOLMOD_INCLUDES})
set(SPARSE_LIBS ${SPARSE_LIBS} ${CHOLMOD_LIBRARIES})
set(EIGEN_TESTED_BACKENDS ${EIGEN_TESTED_BACKENDS} Cholmod)
ei_add_property(EIGEN_TESTED_BACKENDS "Cholmod, ")
else(CHOLMOD_FOUND)
set(EIGEN_MISSING_BACKENDS ${EIGEN_MISSING_BACKENDS} Cholmod)
ei_add_property(EIGEN_MISSING_BACKENDS "Cholmod, ")
endif(CHOLMOD_FOUND)
find_package(Umfpack)
@ -39,9 +43,9 @@ if(UMFPACK_FOUND)
add_definitions("-DEIGEN_UMFPACK_SUPPORT")
include_directories(${UMFPACK_INCLUDES})
set(SPARSE_LIBS ${SPARSE_LIBS} ${UMFPACK_LIBRARIES})
set(EIGEN_TESTED_BACKENDS ${EIGEN_TESTED_BACKENDS} UmfPack)
ei_add_property(EIGEN_TESTED_BACKENDS "UmfPack, ")
else(UMFPACK_FOUND)
set(EIGEN_MISSING_BACKENDS ${EIGEN_MISSING_BACKENDS} UmfPack)
ei_add_property(EIGEN_MISSING_BACKENDS "UmfPack, ")
endif(UMFPACK_FOUND)
find_package(SuperLU)
@ -49,18 +53,18 @@ if(SUPERLU_FOUND)
add_definitions("-DEIGEN_SUPERLU_SUPPORT")
include_directories(${SUPERLU_INCLUDES})
set(SPARSE_LIBS ${SPARSE_LIBS} ${SUPERLU_LIBRARIES})
set(EIGEN_TESTED_BACKENDS ${EIGEN_TESTED_BACKENDS} SuperLU)
ei_add_property(EIGEN_TESTED_BACKENDS "SuperLU, ")
else(SUPERLU_FOUND)
set(EIGEN_MISSING_BACKENDS ${EIGEN_MISSING_BACKENDS} SuperLU)
ei_add_property(EIGEN_MISSING_BACKENDS "SuperLU, ")
endif(SUPERLU_FOUND)
find_package(GoogleHash)
if(GOOGLEHASH_FOUND)
add_definitions("-DEIGEN_GOOGLEHASH_SUPPORT")
include_directories(${GOOGLEHASH_INCLUDES})
set(EIGEN_TESTED_BACKENDS ${EIGEN_TESTED_BACKENDS} GoogleHash)
ei_add_property(EIGEN_TESTED_BACKENDS "GoogleHash, ")
else(GOOGLEHASH_FOUND)
set(EIGEN_MISSING_BACKENDS ${EIGEN_MISSING_BACKENDS} GoogleHash)
ei_add_property(EIGEN_MISSING_BACKENDS "GoogleHash, ")
endif(GOOGLEHASH_FOUND)
option(EIGEN_TEST_NOQT "Disable Qt support in unit tests" OFF)
@ -69,33 +73,12 @@ if(NOT EIGEN_TEST_NOQT)
endif(NOT EIGEN_TEST_NOQT)
if(QT4_FOUND)
include(${QT_USE_FILE})
set(EIGEN_TESTED_BACKENDS ${EIGEN_TESTED_BACKENDS} "Qt4 support")
ei_add_property(EIGEN_TESTED_BACKENDS "Qt4 support, ")
else(QT4_FOUND)
set(EIGEN_MISSING_BACKENDS ${EIGEN_MISSING_BACKENDS} "Qt4 support")
ei_add_property(EIGEN_MISSING_BACKENDS "Qt4 support, ")
endif(QT4_FOUND)
if(CMAKE_COMPILER_IS_GNUCXX)
if(CMAKE_SYSTEM_NAME MATCHES Linux)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -g2")
set(CMAKE_CXX_FLAGS_RELWITHDEBINFO "${CMAKE_CXX_FLAGS_RELWITHDEBINFO} -O2 -g2")
set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} -fno-inline-functions")
set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -O0 -g2")
endif(CMAKE_SYSTEM_NAME MATCHES Linux)
set(EI_OFLAG "-O2")
# MSVC fails with:
# cl : Command line warning D9025 : overriding '/Od' with '/O2'
# cl : Command line error D8016 : '/RTC1' and '/O2' command-line options are incompatible
# elseif(MSVC)
# set(EI_OFLAG "/O2")
else(CMAKE_COMPILER_IS_GNUCXX)
set(EI_OFLAG "")
endif(CMAKE_COMPILER_IS_GNUCXX)
include(EigenTesting)
enable_testing()
if(TEST_LIB)
add_definitions("-DEIGEN_EXTERN_INSTANTIATIONS=1")
endif(TEST_LIB)
@ -154,55 +137,13 @@ ei_add_test(sparse_product)
ei_add_test(sparse_solvers " " "${SPARSE_LIBS}")
ei_add_test(reverse)
# print a summary of the different options
message("************************************************************")
message("*** Eigen's unit tests configuration summary ***")
message("************************************************************")
if(GSL_FOUND)
message("Comparison with GSL: ON")
else(GSL_FOUND)
message("Comparison with GSL: OFF")
endif(GSL_FOUND)
message("Enabled backends: ${EIGEN_TESTED_BACKENDS}")
message("Disabled backends: ${EIGEN_MISSING_BACKENDS}")
if(EIGEN_TEST_SSE2)
message("SSE2: ON")
else(EIGEN_TEST_SSE2)
message("SSE2: AUTO")
endif(EIGEN_TEST_SSE2)
if(EIGEN_TEST_SSE3)
message("SSE3: ON")
else(EIGEN_TEST_SSE3)
message("SSE3: AUTO")
endif(EIGEN_TEST_SSE3)
if(EIGEN_TEST_SSSE3)
message("SSSE3: ON")
else(EIGEN_TEST_SSSE3)
message("SSSE3: AUTO")
endif(EIGEN_TEST_SSSE3)
if(EIGEN_TEST_ALTIVEC)
message("Altivec: ON")
else(EIGEN_TEST_ALTIVEC)
message("Altivec: AUTO")
endif(EIGEN_TEST_ALTIVEC)
if(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("Explicit vec: OFF")
else(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("Explicit vec: AUTO")
endif(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("CXX: ${CMAKE_CXX_COMPILER}")
ei_add_property(EIGEN_TESTING_SUMMARY "CXX: ${CMAKE_CXX_COMPILER}\n")
if(CMAKE_COMPILER_IS_GNUCXX)
execute_process(COMMAND ${CMAKE_CXX_COMPILER} --version COMMAND head -n 1 OUTPUT_VARIABLE EIGEN_CXX_VERSION_STRING OUTPUT_STRIP_TRAILING_WHITESPACE)
message("CXX_VERSION: ${EIGEN_CXX_VERSION_STRING}")
ei_add_property(EIGEN_TESTING_SUMMARY "CXX_VERSION: ${EIGEN_CXX_VERSION_STRING}\n")
endif(CMAKE_COMPILER_IS_GNUCXX)
message("CXX_FLAGS: ${CMAKE_CXX_FLAGS}")
message("Sparse lib flags: ${SPARSE_LIBS}")
ei_add_property(EIGEN_TESTING_SUMMARY "CXX_FLAGS: ${CMAKE_CXX_FLAGS}\n")
ei_add_property(EIGEN_TESTING_SUMMARY "Sparse lib flags: ${SPARSE_LIBS}\n")
message("************************************************************")

View File

@ -2,7 +2,6 @@
add_subdirectory(Eigen)
if(EIGEN_BUILD_TESTS)
include(CTest)
add_subdirectory(test)
endif(EIGEN_BUILD_TESTS)

View File

@ -98,9 +98,9 @@ namespace adtl {
}
namespace Eigen { namespace unsupported { /*@}*/ } }
namespace Eigen {
template<> struct EigenNumTraits<adtl::adouble>
template<> struct NumTraits<adtl::adouble>
{
typedef adtl::adouble Real;
typedef adtl::adouble FloatingPoint;
@ -113,4 +113,60 @@ template<> struct EigenNumTraits<adtl::adouble>
};
};
template<typename Functor> class AdolcForwardJacobian : public Functor
{
typedef adtl::adouble ActiveScalar;
public:
AdolcForwardJacobian() : Functor() {}
AdolcForwardJacobian(const Functor& f) : Functor(f) {}
// forward constructors
template<typename T0>
AdolcForwardJacobian(const T0& a0) : Functor(a0) {}
template<typename T0, typename T1>
AdolcForwardJacobian(const T0& a0, const T1& a1) : Functor(a0, a1) {}
template<typename T0, typename T1, typename T2>
AdolcForwardJacobian(const T0& a0, const T1& a1, const T1& a2) : Functor(a0, a1, a2) {}
typedef typename Functor::InputType InputType;
typedef typename Functor::ValueType ValueType;
typedef typename Functor::JacobianType JacobianType;
typedef Matrix<ActiveScalar, InputType::SizeAtCompileTime, 1> ActiveInput;
typedef Matrix<ActiveScalar, ValueType::SizeAtCompileTime, 1> ActiveValue;
void operator() (const InputType& x, ValueType* v, JacobianType* _jac) const
{
ei_assert(v!=0);
if (!_jac)
{
Functor::operator()(x, v);
return;
}
JacobianType& jac = *_jac;
ActiveInput ax = x.template cast<ActiveScalar>();
ActiveValue av(jac.rows());
for (int j=0; j<jac.cols(); j++)
for (int i=0; i<jac.cols(); i++)
ax[i].setADValue(j, i==j ? 1 : 0);
Functor::operator()(ax, &av);
for (int i=0; i<jac.rows(); i++)
{
(*v)[i] = av[i].getValue();
for (int j=0; j<jac.cols(); j++)
jac.coeffRef(i,j) = av[i].getADValue(j);
}
}
protected:
};
}
#endif // EIGEN_ADLOC_FORWARD

View File

@ -1,6 +1,16 @@
enable_testing()
include(EigenTesting)
# ei_add_test(foo "CXXFLAGS" "libs")
enable_testing()
find_package(Adolc)
include_directories(../../test)
if(ADOLC_FOUND)
include_directories(${ADOLC_INCLUDES})
ei_add_property(EIGEN_TESTED_BACKENDS "Adolc")
ei_add_test(forward_adolc " " ${ADOLC_LIBRARIES})
else(ADOLC_FOUND)
ei_add_property(EIGEN_MISSING_BACKENDS "Adolc")
endif(ADOLC_FOUND)

View File

@ -0,0 +1,132 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra. Eigen itself is part of the KDE project.
//
// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
//
// Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 3 of the License, or (at your option) any later version.
//
// Alternatively, 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 2 of
// the License, or (at your option) any later version.
//
// Eigen 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 Lesser General Public License or the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License and a copy of the GNU General Public License along with
// Eigen. If not, see <http://www.gnu.org/licenses/>.
#include "main.h"
#define NUMBER_DIRECTIONS 16
#include <unsupported/Eigen/AdolcForward>
int adtl::ADOLC_numDir;
template<typename _Scalar, int NX, int NY>
struct TestFunc1
{
typedef _Scalar Scalar;
enum {
InputsAtCompileTime = NX,
ValuesAtCompileTime = NY
};
typedef Matrix<Scalar,InputsAtCompileTime,1> InputType;
typedef Matrix<Scalar,ValuesAtCompileTime,1> ValueType;
typedef Matrix<Scalar,ValuesAtCompileTime,InputsAtCompileTime> JacobianType;
template<typename T>
void operator() (const Matrix<T,InputsAtCompileTime,1>& x, Matrix<T,ValuesAtCompileTime,1>* _v) const
{
Matrix<T,ValuesAtCompileTime,1>& v = *_v;
v[0] = 2 * x[0] * x[0] + x[0] * x[1];
v[1] = 3 * x[1] * x[0] + 0.5 * x[1] * x[1];
if(NX>2)
{
v[0] += 0.5 * x[2];
v[1] += x[2];
}
if(NY>2)
{
v[2] = 3 * x[1] * x[0] * x[0];
}
if (NX>2 && NY>2)
v[2] *= x[2];
}
void operator() (const InputType& x, ValueType* v, JacobianType* _j) const
{
(*this)(x, v);
if(_j)
{
JacobianType& j = *_j;
j(0,0) = 4 * x[0] + x[1];
j(1,0) = 3 * x[1];
j(0,1) = x[0];
j(1,1) = 3 * x[0] + 2 * 0.5 * x[1];
if (NX>2)
{
j(0,2) = 0.5;
j(1,2) = 1;
}
if(NY>2)
{
j(2,0) = 3 * x[1] * 2 * x[0];
j(2,1) = 3 * x[0] * x[0];
}
if (NX>2 && NY>2)
{
j(2,0) *= x[2];
j(2,1) *= x[2];
j(2,2) = 3 * x[1] * x[0] * x[0];
j(2,2) = 3 * x[1] * x[0] * x[0];
}
}
}
};
template<typename Func> void adolc_forward_jacobian(const Func& f)
{
typename Func::InputType x = Func::InputType::Random();
typename Func::ValueType y, yref;
typename Func::JacobianType j, jref;
jref.setZero();
yref.setZero();
f(x,&yref,&jref);
// std::cerr << y.transpose() << "\n\n";;
// std::cerr << j << "\n\n";;
j.setZero();
y.setZero();
AdolcForwardJacobian<Func> autoj(f);
autoj(x, &y, &j);
// std::cerr << y.transpose() << "\n\n";;
// std::cerr << j << "\n\n";;
VERIFY_IS_APPROX(y, yref);
VERIFY_IS_APPROX(j, jref);
}
void test_forward_adolc()
{
adtl::ADOLC_numDir = NUMBER_DIRECTIONS;
for(int i = 0; i < g_repeat; i++) {
CALL_SUBTEST(( adolc_forward_jacobian(TestFunc1<double,2,2>()) ));
CALL_SUBTEST(( adolc_forward_jacobian(TestFunc1<double,2,3>()) ));
CALL_SUBTEST(( adolc_forward_jacobian(TestFunc1<double,3,2>()) ));
CALL_SUBTEST(( adolc_forward_jacobian(TestFunc1<double,3,3>()) ));
}
}