#pragma once #ifndef BERTINI_PYTHON_EIGENPY_INTERACTION_HPP #define BERTINI_PYTHON_EIGENPY_INTERACTION_HPP #include "python_common.hpp" #include #include #include // this code derived from // https://github.com/stack-of-tasks/eigenpy/issues/365 // where I asked about using custom types, and @jcarpent responded with a discussion // of an application of this in Pinnochio, a library for rigid body dynamics. namespace eigenpy { namespace internal { // template specialization for real numbers template <> struct getitem { using NumT = bertini::mpfr_float; static PyObject* run(void* data, void* /* arr */) { NumT & mpfr_scalar = *static_cast(data); auto & backend = mpfr_scalar.backend(); if(backend.data()[0]._mpfr_d == 0) // If the mpfr_scalar is not initialized, we have to init it. { mpfr_scalar = NumT(0); } boost::python::object m(boost::ref(mpfr_scalar)); Py_INCREF(m.ptr()); return m.ptr(); } }; // a template specialization for complex numbers template <> struct getitem { using NumT = bertini::mpfr_complex; static PyObject* run(void* data, void* /* arr */) { NumT & mpfr_scalar = *static_cast(data); auto & backend = mpfr_scalar.backend(); if(backend.data()[0].re->_mpfr_d == 0) // If the mpfr_scalar is not initialized, we have to init it. { mpfr_scalar = NumT(0); } boost::python::object m(boost::ref(mpfr_scalar)); Py_INCREF(m.ptr()); return m.ptr(); } }; } // namespace internal // i lifted this from EigenPy and adapted it, basically removing the calls for the comparitors. template void registerUfunct_without_comparitors(){ const int type_code = Register::getTypeCode(); PyObject *numpy_str; #if PY_MAJOR_VERSION >= 3 numpy_str = PyUnicode_FromString("numpy"); #else numpy_str = PyString_FromString("numpy"); #endif PyObject *numpy; numpy = PyImport_Import(numpy_str); Py_DECREF(numpy_str); import_ufunc(); // Matrix multiply { int types[3] = {type_code, type_code, type_code}; std::stringstream ss; ss << "return result of multiplying two matrices of "; ss << bp::type_info(typeid(Scalar)).name(); PyUFuncObject *ufunc = (PyUFuncObject *)PyObject_GetAttrString(numpy, "matmul"); if (!ufunc) { std::stringstream ss; ss << "Impossible to define matrix_multiply for given type " << bp::type_info(typeid(Scalar)).name() << std::endl; eigenpy::Exception(ss.str()); } if (PyUFunc_RegisterLoopForType((PyUFuncObject *)ufunc, type_code, &internal::gufunc_matrix_multiply, types, 0) < 0) { std::stringstream ss; ss << "Impossible to register matrix_multiply for given type " << bp::type_info(typeid(Scalar)).name() << std::endl; eigenpy::Exception(ss.str()); } Py_DECREF(ufunc); } // Binary operators EIGENPY_REGISTER_BINARY_UFUNC(add, type_code, Scalar, Scalar, Scalar); EIGENPY_REGISTER_BINARY_UFUNC(subtract, type_code, Scalar, Scalar, Scalar); EIGENPY_REGISTER_BINARY_UFUNC(multiply, type_code, Scalar, Scalar, Scalar); EIGENPY_REGISTER_BINARY_UFUNC(divide, type_code, Scalar, Scalar, Scalar); // Comparison operators EIGENPY_REGISTER_BINARY_UFUNC(equal, type_code, Scalar, Scalar, bool); EIGENPY_REGISTER_BINARY_UFUNC(not_equal, type_code, Scalar, Scalar, bool); //these are commented out because the comparisons are NOT defined for complex types!! // EIGENPY_REGISTER_BINARY_UFUNC(greater, type_code, Scalar, Scalar, bool); // EIGENPY_REGISTER_BINARY_UFUNC(less, type_code, Scalar, Scalar, bool); // EIGENPY_REGISTER_BINARY_UFUNC(greater_equal, type_code, Scalar, Scalar, bool); // EIGENPY_REGISTER_BINARY_UFUNC(less_equal, type_code, Scalar, Scalar, bool); // Unary operators EIGENPY_REGISTER_UNARY_UFUNC(negative, type_code, Scalar, Scalar); Py_DECREF(numpy); } } // namespace eigenpy namespace bertini{ namespace python{ void EnableEigenPy(); }} // namespaces #endif // include guard