Reputation: 1
does C++ allow using external functions in template specification implementations? What I wanna do is make a wrapper to BLAS/LAPACK (the CBLAS doesn't cut it for me, I need a couple of functions that are not in there and are in the external BLAS), somewhat like this
extern "C" {
void saxpy_(int* n, const float* const sa, const float* const sx, int* incx, float* sy, int* incy);
void daxpy_(int* n, const double* const sa, const double* const sx, int* incx, double* sy, int* incy);
}
namespace matrix {
template<class T>
void ax_plus_y(const T a, const Matrix<T, Dynamic, Dynamic>& x, Matrix<T, Dynamic, Dynamic>& y) {
y += a * x;
}
template<>
void ax_plus_y<float>(const float a, const Matrix<float, Dynamic, 1>& x, Matrix<float, Dynamic, 1>& y) {
int n = x.size();
int incx = x.innerStride();
int incy = y.innerStride();
saxpy_(&n, &a, x.data, &incx, y.data(), &incy);
}
template<>
void ax_plus_y<double>(const double a, const Matrix<double, Dynamic, 1>& x, Matrix<double, Dynamic, 1>& y) {
int n = x.size();
int incx = x.innerStride();
int incy = y.innerStride();
daxpy_(&n, &a, x.data, &incx, y.data(), &incy);
}
}
That's just a snippet of what I need, but is this possible at all? I can't really use only one function, 'cause I need the different external calls for each data type.
P.S. The Matrix
is from Eigen, shouldn't be relevant to concept.
EDIT: To clarify, this is what g++ throws up with:
error: template-id ‘ax_plus_y<float>’ for ‘void matrix::ax_plus_y(float, const Eigen::Matrix<float, -0x00000000000000001, 1>&, Eigen::Matrix<float, -0x00000000000000001, 1>&)’ does not match any template declaration
Upvotes: 1
Views: 507
Reputation: 42805
The problem I see is that you declare three things:
template<class T>
void ax_plus_y(const T a, const Matrix<T, Dynamic, Dynamic>& x,
Matrix<T, Dynamic, Dynamic>& y);
template<>
void ax_plus_y<float>(const float a, const Matrix<float, Dynamic, 1>& x,
Matrix<float, Dynamic, 1>& y);
template<>
void ax_plus_y<double>(const double a, const Matrix<double, Dynamic, 1>& x,
Matrix<double, Dynamic, 1>& y);
The function signatures of the specializations, which use Matrix<double, Dynamic, 1>
, do not match the general function signature, which uses Matrix<double, Dynamic, Dynamic>
, so the compiler will not associate the specializations with the general template. I don't think that the compiler is going to jump through implicit-typecast hoops to make this work for you.
EDIT: As to your general "is it possible" question: Yes, but it's going to take a lot of effort, and you might want to consider writing some kind of code-generating program to keep you from going insane.
I've done this years ago, but not with Matrix
or Vector
classes, just with const T*, size_t
pairs. My goal was legibility (BLAS/LAPACK function names are more about brevity than clarity), support for other libraries such as FFTW and SSE, and support for types outside those packages. The idea was that you'd use the general routines and it would generate the best code for your platform. Unfortunately, I can't open-source it and it's a decade out of support.
In your case, you're going to have to write something like:
template<class T, long long D1=Dynamic, long long D2=Dynamic>
void ax_plus_y(const T a, const Matrix<T, D1, D2>& x, Matrix<T, D1, D2>& y) {
y += a * x;
}
template<long long D1=Dynamic>
void ax_plus_y<float, D1, 1>(const float a, const Matrix<float, D1, 1>& x, Matrix<float, D1, 1>& y) {
int n = x.size();
int incx = x.innerStride();
int incy = y.innerStride();
saxpy_(&n, &a, x.data, &incx, y.data(), &incy);
}
template<long long D2=Dynamic>
void ax_plus_y<float, 1, D2>(const float a, const Matrix<float, 1, D2>& x, Matrix<float, 1, D2>& y) {
int n = x.size();
int incx = x.outerStride();
int incy = y.outerStride();
saxpy_(&n, &a, x.data, &incx, y.data(), &incy);
}
and repeat for double
.
Upvotes: 0
Reputation: 1318
There is no problem with extern functions in template functions. See e.g. the following snippet that works fine:
#include <iostream>
#include <math.h>
template<typename T>
void doit(T number)
{
std::cout << cos(number) << std::endl;
}
template<>
void doit<float>(float number)
{
std::cout << "a float: " << cos(number) << std::endl;
}
int main()
{
doit((double)0.0);
doit((float)1.0);
}
Maybe it's a problem that the first template defintion is in namespace "matrix" while the others aren't?
Upvotes: 0