|
|
|
@ -4,70 +4,70 @@ |
|
|
|
#include <algorithm> |
|
|
|
|
|
|
|
extern "C"{ |
|
|
|
DLLEXPORT float s_matrix_norm(char norm, int m, int n, float a[], float work[]) |
|
|
|
DLLEXPORT float s_matrix_norm(char norm, MKL_INT m, MKL_INT n, float a[], float work[]) |
|
|
|
{ |
|
|
|
return slange_(&norm, &m, &n, a, &m, work); |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT double d_matrix_norm(char norm, int m, int n, double a[], double work[]) |
|
|
|
DLLEXPORT double d_matrix_norm(char norm, MKL_INT m, MKL_INT n, double a[], double work[]) |
|
|
|
{ |
|
|
|
return dlange_(&norm, &m, &n, a, &m, work); |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT float c_matrix_norm(char norm, int m, int n, MKL_Complex8 a[], float work[]) |
|
|
|
DLLEXPORT float c_matrix_norm(char norm, MKL_INT m, MKL_INT n, MKL_Complex8 a[], float work[]) |
|
|
|
{ |
|
|
|
return clange_(&norm, &m, &n, a, &m, work); |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT double z_matrix_norm(char norm, int m, int n, MKL_Complex16 a[], double work[]) |
|
|
|
DLLEXPORT double z_matrix_norm(char norm, MKL_INT m, MKL_INT n, MKL_Complex16 a[], double work[]) |
|
|
|
{ |
|
|
|
return zlange_(&norm, &m, &n, a, &m, work); |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_lu_factor(int m, float a[], int ipiv[]) |
|
|
|
DLLEXPORT MKL_INT s_lu_factor(MKL_INT m, float a[], MKL_INT ipiv[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
sgetrf_(&m,&m,a,&m,ipiv,&info); |
|
|
|
for(int i = 0; i < m; ++i ){ |
|
|
|
for(MKL_INT i = 0; i < m; ++i ){ |
|
|
|
ipiv[i] -= 1; |
|
|
|
} |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_lu_factor(int m, double a[], int ipiv[]) |
|
|
|
DLLEXPORT MKL_INT d_lu_factor(MKL_INT m, double a[], MKL_INT ipiv[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
dgetrf_(&m,&m,a,&m,ipiv,&info); |
|
|
|
for(int i = 0; i < m; ++i ){ |
|
|
|
for(MKL_INT i = 0; i < m; ++i ){ |
|
|
|
ipiv[i] -= 1; |
|
|
|
} |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_lu_factor(int m, MKL_Complex8 a[], int ipiv[]) |
|
|
|
DLLEXPORT MKL_INT c_lu_factor(MKL_INT m, MKL_Complex8 a[], MKL_INT ipiv[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
cgetrf_(&m,&m,a,&m,ipiv,&info); |
|
|
|
for(int i = 0; i < m; ++i ){ |
|
|
|
for(MKL_INT i = 0; i < m; ++i ){ |
|
|
|
ipiv[i] -= 1; |
|
|
|
} |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_lu_factor(int m, MKL_Complex16 a[], int ipiv[]) |
|
|
|
DLLEXPORT MKL_INT z_lu_factor(MKL_INT m, MKL_Complex16 a[], MKL_INT ipiv[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
zgetrf_(&m,&m,a,&m,ipiv,&info); |
|
|
|
for(int i = 0; i < m; ++i ){ |
|
|
|
for(MKL_INT i = 0; i < m; ++i ){ |
|
|
|
ipiv[i] -= 1; |
|
|
|
} |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_lu_inverse(int n, float a[], float work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT s_lu_inverse(MKL_INT n, float a[], float work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
sgetrf_(&n,&n,a,&n,ipiv,&info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -80,10 +80,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_lu_inverse(int n, double a[], double work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT d_lu_inverse(MKL_INT n, double a[], double work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
dgetrf_(&n,&n,a,&n,ipiv,&info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -96,10 +96,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_lu_inverse(int n, MKL_Complex8 a[], MKL_Complex8 work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT c_lu_inverse(MKL_INT n, MKL_Complex8 a[], MKL_Complex8 work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
cgetrf_(&n,&n,a,&n,ipiv,&info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -112,10 +112,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_lu_inverse(int n, MKL_Complex16 a[], MKL_Complex16 work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT z_lu_inverse(MKL_INT n, MKL_Complex16 a[], MKL_Complex16 work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
zgetrf_(&n,&n,a,&n,ipiv,&info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -128,13 +128,13 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_lu_inverse_factored(int n, float a[], int ipiv[], float work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT s_lu_inverse_factored(MKL_INT n, float a[], MKL_INT ipiv[], float work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int i; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
sgetri_(&n,a,&n,ipiv,work,&lwork,&info); |
|
|
|
|
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
@ -143,14 +143,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_lu_inverse_factored(int n, double a[], int ipiv[], double work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT d_lu_inverse_factored(MKL_INT n, double a[], MKL_INT ipiv[], double work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int i; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
|
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
dgetri_(&n,a,&n,ipiv,work,&lwork,&info); |
|
|
|
|
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
@ -159,14 +159,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_lu_inverse_factored(int n, MKL_Complex8 a[], int ipiv[], MKL_Complex8 work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT c_lu_inverse_factored(MKL_INT n, MKL_Complex8 a[], MKL_INT ipiv[], MKL_Complex8 work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int i; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
|
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
cgetri_(&n,a,&n,ipiv,work,&lwork,&info); |
|
|
|
|
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
@ -175,14 +175,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_lu_inverse_factored(int n, MKL_Complex16 a[], int ipiv[], MKL_Complex16 work[], int lwork) |
|
|
|
DLLEXPORT MKL_INT z_lu_inverse_factored(MKL_INT n, MKL_Complex16 a[], MKL_INT ipiv[], MKL_Complex16 work[], MKL_INT lwork) |
|
|
|
{ |
|
|
|
int i; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
|
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
zgetri_(&n,a,&n,ipiv,work,&lwork,&info); |
|
|
|
|
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
@ -191,10 +191,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_lu_solve_factored(int n, int nrhs, float a[], int ipiv[], float b[]) |
|
|
|
DLLEXPORT MKL_INT s_lu_solve_factored(MKL_INT n, MKL_INT nrhs, float a[], MKL_INT ipiv[], float b[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
int i; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
@ -207,10 +207,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_lu_solve_factored(int n, int nrhs, double a[], int ipiv[], double b[]) |
|
|
|
DLLEXPORT MKL_INT d_lu_solve_factored(MKL_INT n, MKL_INT nrhs, double a[], MKL_INT ipiv[], double b[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
int i; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
@ -223,10 +223,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_lu_solve_factored(int n, int nrhs, MKL_Complex8 a[], int ipiv[], MKL_Complex8 b[]) |
|
|
|
DLLEXPORT MKL_INT c_lu_solve_factored(MKL_INT n, MKL_INT nrhs, MKL_Complex8 a[], MKL_INT ipiv[], MKL_Complex8 b[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
int i; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
@ -239,10 +239,10 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_lu_solve_factored(int n, int nrhs, MKL_Complex16 a[], int ipiv[], MKL_Complex16 b[]) |
|
|
|
DLLEXPORT MKL_INT z_lu_solve_factored(MKL_INT n, MKL_INT nrhs, MKL_Complex16 a[], MKL_INT ipiv[], MKL_Complex16 b[]) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
int i; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_INT i; |
|
|
|
for(i = 0; i < n; ++i ){ |
|
|
|
ipiv[i] += 1; |
|
|
|
} |
|
|
|
@ -255,13 +255,13 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_lu_solve(int n, int nrhs, float a[], float b[]) |
|
|
|
DLLEXPORT MKL_INT s_lu_solve(MKL_INT n, MKL_INT nrhs, float a[], float b[]) |
|
|
|
{ |
|
|
|
float* clone = new float[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(float)); |
|
|
|
|
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
sgetrf_(&n, &n, clone, &n, ipiv, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -277,13 +277,13 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_lu_solve(int n, int nrhs, double a[], double b[]) |
|
|
|
DLLEXPORT MKL_INT d_lu_solve(MKL_INT n, MKL_INT nrhs, double a[], double b[]) |
|
|
|
{ |
|
|
|
double* clone = new double[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(double)); |
|
|
|
|
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
dgetrf_(&n, &n, clone, &n, ipiv, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -299,13 +299,13 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_lu_solve(int n, int nrhs, MKL_Complex8 a[], MKL_Complex8 b[]) |
|
|
|
DLLEXPORT MKL_INT c_lu_solve(MKL_INT n, MKL_INT nrhs, MKL_Complex8 a[], MKL_Complex8 b[]) |
|
|
|
{ |
|
|
|
MKL_Complex8* clone = new MKL_Complex8[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(MKL_Complex8)); |
|
|
|
|
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
cgetrf_(&n, &n, clone, &n, ipiv, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -321,13 +321,13 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_lu_solve(int n, int nrhs, MKL_Complex16 a[], MKL_Complex16 b[]) |
|
|
|
DLLEXPORT MKL_INT z_lu_solve(MKL_INT n, MKL_INT nrhs, MKL_Complex16 a[], MKL_Complex16 b[]) |
|
|
|
{ |
|
|
|
MKL_Complex16* clone = new MKL_Complex16[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(MKL_Complex16)); |
|
|
|
|
|
|
|
int* ipiv = new int[n]; |
|
|
|
int info = 0; |
|
|
|
MKL_INT* ipiv = new MKL_INT[n]; |
|
|
|
MKL_INT info = 0; |
|
|
|
zgetrf_(&n, &n, clone, &n, ipiv, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -343,14 +343,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_cholesky_factor(int n, float a[]){ |
|
|
|
DLLEXPORT MKL_INT s_cholesky_factor(MKL_INT n, float a[]){ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
spotrf_(&uplo, &n, a, &n, &info); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
int index = i * n; |
|
|
|
for (int j = 0; j < n && i > j; ++j) |
|
|
|
MKL_INT index = i * n; |
|
|
|
for (MKL_INT j = 0; j < n && i > j; ++j) |
|
|
|
{ |
|
|
|
a[index + j] = 0; |
|
|
|
} |
|
|
|
@ -358,14 +358,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_cholesky_factor(int n, double* a){ |
|
|
|
DLLEXPORT MKL_INT d_cholesky_factor(MKL_INT n, double* a){ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
dpotrf_(&uplo, &n, a, &n, &info); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
int index = i * n; |
|
|
|
for (int j = 0; j < n && i > j; ++j) |
|
|
|
MKL_INT index = i * n; |
|
|
|
for (MKL_INT j = 0; j < n && i > j; ++j) |
|
|
|
{ |
|
|
|
a[index + j] = 0; |
|
|
|
} |
|
|
|
@ -373,15 +373,15 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_cholesky_factor(int n, MKL_Complex8 a[]){ |
|
|
|
DLLEXPORT MKL_INT c_cholesky_factor(MKL_INT n, MKL_Complex8 a[]){ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_Complex8 zero = {0.0f, 0.0f}; |
|
|
|
cpotrf_(&uplo, &n, a, &n, &info); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
int index = i * n; |
|
|
|
for (int j = 0; j < n && i > j; ++j) |
|
|
|
MKL_INT index = i * n; |
|
|
|
for (MKL_INT j = 0; j < n && i > j; ++j) |
|
|
|
{ |
|
|
|
a[index + j] = zero; |
|
|
|
} |
|
|
|
@ -389,15 +389,15 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_cholesky_factor(int n, MKL_Complex16 a[]){ |
|
|
|
DLLEXPORT MKL_INT z_cholesky_factor(MKL_INT n, MKL_Complex16 a[]){ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_Complex16 zero = {0.0, 0.0}; |
|
|
|
zpotrf_(&uplo, &n, a, &n, &info); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
int index = i * n; |
|
|
|
for (int j = 0; j < n && i > j; ++j) |
|
|
|
MKL_INT index = i * n; |
|
|
|
for (MKL_INT j = 0; j < n && i > j; ++j) |
|
|
|
{ |
|
|
|
a[index + j] = zero; |
|
|
|
} |
|
|
|
@ -405,12 +405,12 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_cholesky_solve(int n, int nrhs, float a[], float b[]) |
|
|
|
DLLEXPORT MKL_INT s_cholesky_solve(MKL_INT n, MKL_INT nrhs, float a[], float b[]) |
|
|
|
{ |
|
|
|
float* clone = new float[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(float)); |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
spotrf_(&uplo, &n, clone, &n, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -422,12 +422,12 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_cholesky_solve(int n, int nrhs, double a[], double b[]) |
|
|
|
DLLEXPORT MKL_INT d_cholesky_solve(MKL_INT n, MKL_INT nrhs, double a[], double b[]) |
|
|
|
{ |
|
|
|
double* clone = new double[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(double)); |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
dpotrf_(&uplo, &n, clone, &n, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -439,12 +439,12 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_cholesky_solve(int n, int nrhs, MKL_Complex8 a[], MKL_Complex8 b[]) |
|
|
|
DLLEXPORT MKL_INT c_cholesky_solve(MKL_INT n, MKL_INT nrhs, MKL_Complex8 a[], MKL_Complex8 b[]) |
|
|
|
{ |
|
|
|
MKL_Complex8* clone = new MKL_Complex8[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(MKL_Complex8)); |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
cpotrf_(&uplo, &n, clone, &n, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -456,12 +456,12 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_cholesky_solve(int n, int nrhs, MKL_Complex16 a[], MKL_Complex16 b[]) |
|
|
|
DLLEXPORT MKL_INT z_cholesky_solve(MKL_INT n, MKL_INT nrhs, MKL_Complex16 a[], MKL_Complex16 b[]) |
|
|
|
{ |
|
|
|
MKL_Complex16* clone = new MKL_Complex16[n*n]; |
|
|
|
std::memcpy(clone, a, n*n*sizeof(MKL_Complex16)); |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
zpotrf_(&uplo, &n, clone, &n, &info); |
|
|
|
|
|
|
|
if (info != 0){ |
|
|
|
@ -473,46 +473,46 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_cholesky_solve_factored(int n, int nrhs, float a[], float b[]) |
|
|
|
DLLEXPORT MKL_INT s_cholesky_solve_factored(MKL_INT n, MKL_INT nrhs, float a[], float b[]) |
|
|
|
{ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
spotrs_(&uplo, &n, &nrhs, a, &n, b, &n, &info); |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_cholesky_solve_factored(int n, int nrhs, double a[], double b[]) |
|
|
|
DLLEXPORT MKL_INT d_cholesky_solve_factored(MKL_INT n, MKL_INT nrhs, double a[], double b[]) |
|
|
|
{ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
dpotrs_(&uplo, &n, &nrhs, a, &n, b, &n, &info); |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_cholesky_solve_factored(int n, int nrhs, MKL_Complex8 a[], MKL_Complex8 b[]) |
|
|
|
DLLEXPORT MKL_INT c_cholesky_solve_factored(MKL_INT n, MKL_INT nrhs, MKL_Complex8 a[], MKL_Complex8 b[]) |
|
|
|
{ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
cpotrs_(&uplo, &n, &nrhs, a, &n, b, &n, &info); |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_cholesky_solve_factored(int n, int nrhs, MKL_Complex16 a[], MKL_Complex16 b[]) |
|
|
|
DLLEXPORT MKL_INT z_cholesky_solve_factored(MKL_INT n, MKL_INT nrhs, MKL_Complex16 a[], MKL_Complex16 b[]) |
|
|
|
{ |
|
|
|
char uplo = 'L'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
zpotrs_(&uplo, &n, &nrhs, a, &n, b, &n, &info); |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_qr_factor(int m, int n, float r[], float tau[], float q[], float work[], int len) |
|
|
|
DLLEXPORT MKL_INT s_qr_factor(MKL_INT m, MKL_INT n, float r[], float tau[], float q[], float work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
sgeqrf_(&m, &n, r, &m, tau, work, &len, &info); |
|
|
|
|
|
|
|
for (int i = 0; i < m; ++i) |
|
|
|
for (MKL_INT i = 0; i < m; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < m && j < n; ++j) |
|
|
|
for (MKL_INT j = 0; j < m && j < n; ++j) |
|
|
|
{ |
|
|
|
if (i > j) |
|
|
|
{ |
|
|
|
@ -534,14 +534,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_qr_factor(int m, int n, double r[], double tau[], double q[], double work[], int len) |
|
|
|
DLLEXPORT MKL_INT d_qr_factor(MKL_INT m, MKL_INT n, double r[], double tau[], double q[], double work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
dgeqrf_(&m, &n, r, &m, tau, work, &len, &info); |
|
|
|
|
|
|
|
for (int i = 0; i < m; ++i) |
|
|
|
for (MKL_INT i = 0; i < m; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < m && j < n; ++j) |
|
|
|
for (MKL_INT j = 0; j < m && j < n; ++j) |
|
|
|
{ |
|
|
|
if (i > j) |
|
|
|
{ |
|
|
|
@ -563,14 +563,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_qr_factor(int m, int n, MKL_Complex8 r[], MKL_Complex8 tau[], MKL_Complex8 q[], MKL_Complex8 work[], int len) |
|
|
|
DLLEXPORT MKL_INT c_qr_factor(MKL_INT m, MKL_INT n, MKL_Complex8 r[], MKL_Complex8 tau[], MKL_Complex8 q[], MKL_Complex8 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
cgeqrf_(&m, &n, r, &m, tau, work, &len, &info); |
|
|
|
|
|
|
|
for (int i = 0; i < m; ++i) |
|
|
|
for (MKL_INT i = 0; i < m; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < m && j < n; ++j) |
|
|
|
for (MKL_INT j = 0; j < m && j < n; ++j) |
|
|
|
{ |
|
|
|
if (i > j) |
|
|
|
{ |
|
|
|
@ -592,14 +592,14 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_qr_factor(int m, int n, MKL_Complex16 r[], MKL_Complex16 tau[], MKL_Complex16 q[], MKL_Complex16 work[], int len) |
|
|
|
DLLEXPORT MKL_INT z_qr_factor(MKL_INT m, MKL_INT n, MKL_Complex16 r[], MKL_Complex16 tau[], MKL_Complex16 q[], MKL_Complex16 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
zgeqrf_(&m, &n, r, &m, tau, work, &len, &info); |
|
|
|
|
|
|
|
for (int i = 0; i < m; ++i) |
|
|
|
for (MKL_INT i = 0; i < m; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < m && j < n; ++j) |
|
|
|
for (MKL_INT j = 0; j < m && j < n; ++j) |
|
|
|
{ |
|
|
|
if (i > j) |
|
|
|
{ |
|
|
|
@ -621,9 +621,9 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_qr_solve(int m, int n, int bn, float r[], float b[], float x[], float work[], int len) |
|
|
|
DLLEXPORT MKL_INT s_qr_solve(MKL_INT m, MKL_INT n, MKL_INT bn, float r[], float b[], float x[], float work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
float* clone_r = new float[m*n]; |
|
|
|
std::memcpy(clone_r, r, m*n*sizeof(float)); |
|
|
|
|
|
|
|
@ -644,9 +644,9 @@ extern "C"{ |
|
|
|
char tran = 'T'; |
|
|
|
sormqr_(&side, &tran, &m, &bn, &n, clone_r, &m, tau, clone_b, &m, work, &len, &info); |
|
|
|
cblas_strsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, 1.0, clone_r, m, clone_b, m); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -658,9 +658,9 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_qr_solve(int m, int n, int bn, double r[], double b[], double x[], double work[], int len) |
|
|
|
DLLEXPORT MKL_INT d_qr_solve(MKL_INT m, MKL_INT n, MKL_INT bn, double r[], double b[], double x[], double work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
double* clone_r = new double[m*n]; |
|
|
|
std::memcpy(clone_r, r, m*n*sizeof(double)); |
|
|
|
|
|
|
|
@ -682,9 +682,9 @@ extern "C"{ |
|
|
|
|
|
|
|
dormqr_(&side, &tran, &m, &bn, &n, clone_r, &m, tau, clone_b, &m, work, &len, &info); |
|
|
|
cblas_dtrsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, 1.0, clone_r, m, clone_b, m); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -696,9 +696,9 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_qr_solve(int m, int n, int bn, MKL_Complex8 r[], MKL_Complex8 b[], MKL_Complex8 x[], MKL_Complex8 work[], int len) |
|
|
|
DLLEXPORT MKL_INT c_qr_solve(MKL_INT m, MKL_INT n, MKL_INT bn, MKL_Complex8 r[], MKL_Complex8 b[], MKL_Complex8 x[], MKL_Complex8 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_Complex8* clone_r = new MKL_Complex8[m*n]; |
|
|
|
std::memcpy(clone_r, r, m*n*sizeof(MKL_Complex8)); |
|
|
|
|
|
|
|
@ -722,9 +722,9 @@ extern "C"{ |
|
|
|
MKL_Complex8 one = {1.0, 0.0}; |
|
|
|
cblas_ctrsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, &one, clone_r, m, clone_b, m); |
|
|
|
|
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -736,9 +736,9 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_qr_solve(int m, int n, int bn, MKL_Complex16 r[], MKL_Complex16 b[], MKL_Complex16 x[], MKL_Complex16 work[], int len) |
|
|
|
DLLEXPORT MKL_INT z_qr_solve(MKL_INT m, MKL_INT n, MKL_INT bn, MKL_Complex16 r[], MKL_Complex16 b[], MKL_Complex16 x[], MKL_Complex16 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_Complex16* clone_r = new MKL_Complex16[m*n]; |
|
|
|
std::memcpy(clone_r, r, m*n*sizeof(MKL_Complex16)); |
|
|
|
|
|
|
|
@ -762,9 +762,9 @@ extern "C"{ |
|
|
|
MKL_Complex16 one = {1.0, 0.0}; |
|
|
|
cblas_ztrsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, &one, clone_r, m, clone_b, m); |
|
|
|
|
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -776,20 +776,20 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_qr_solve_factored(int m, int n, int bn, float r[], float b[], float tau[], float x[], float work[], int len) |
|
|
|
DLLEXPORT MKL_INT s_qr_solve_factored(MKL_INT m, MKL_INT n, MKL_INT bn, float r[], float b[], float tau[], float x[], float work[], MKL_INT len) |
|
|
|
{ |
|
|
|
char side ='L'; |
|
|
|
char tran = 'T'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
|
|
|
|
float* clone_b = new float[m*bn]; |
|
|
|
std::memcpy(clone_b, b, m*bn*sizeof(float)); |
|
|
|
|
|
|
|
sormqr_(&side, &tran, &m, &bn, &n, r, &m, tau, clone_b, &m, work, &len, &info); |
|
|
|
cblas_strsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, 1.0, r, m, clone_b, m); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -799,20 +799,20 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_qr_solve_factored(int m, int n, int bn, double r[], double b[], double tau[], double x[], double work[], int len) |
|
|
|
DLLEXPORT MKL_INT d_qr_solve_factored(MKL_INT m, MKL_INT n, MKL_INT bn, double r[], double b[], double tau[], double x[], double work[], MKL_INT len) |
|
|
|
{ |
|
|
|
char side ='L'; |
|
|
|
char tran = 'T'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
|
|
|
|
double* clone_b = new double[m*bn]; |
|
|
|
std::memcpy(clone_b, b, m*bn*sizeof(double)); |
|
|
|
|
|
|
|
dormqr_(&side, &tran, &m, &bn, &n, r, &m, tau, clone_b, &m, work, &len, &info); |
|
|
|
cblas_dtrsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, 1.0, r, m, clone_b, m); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -822,11 +822,11 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_qr_solve_factored(int m, int n, int bn, MKL_Complex8 r[], MKL_Complex8 b[], MKL_Complex8 tau[], MKL_Complex8 x[], MKL_Complex8 work[], int len) |
|
|
|
DLLEXPORT MKL_INT c_qr_solve_factored(MKL_INT m, MKL_INT n, MKL_INT bn, MKL_Complex8 r[], MKL_Complex8 b[], MKL_Complex8 tau[], MKL_Complex8 x[], MKL_Complex8 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
char side ='L'; |
|
|
|
char tran = 'C'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
|
|
|
|
MKL_Complex8* clone_b = new MKL_Complex8[m*bn]; |
|
|
|
std::memcpy(clone_b, b, m*bn*sizeof(MKL_Complex8)); |
|
|
|
@ -834,9 +834,9 @@ extern "C"{ |
|
|
|
cunmqr_(&side, &tran, &m, &bn, &n, r, &m, tau, clone_b, &m, work, &len, &info); |
|
|
|
MKL_Complex8 one = {1.0f, 0.0f}; |
|
|
|
cblas_ctrsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, &one, r, m, clone_b, m); |
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -846,11 +846,11 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_qr_solve_factored(int m, int n, int bn, MKL_Complex16 r[], MKL_Complex16 b[], MKL_Complex16 tau[], MKL_Complex16 x[], MKL_Complex16 work[], int len) |
|
|
|
DLLEXPORT MKL_INT z_qr_solve_factored(MKL_INT m, MKL_INT n, MKL_INT bn, MKL_Complex16 r[], MKL_Complex16 b[], MKL_Complex16 tau[], MKL_Complex16 x[], MKL_Complex16 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
char side ='L'; |
|
|
|
char tran = 'C'; |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
|
|
|
|
MKL_Complex16* clone_b = new MKL_Complex16[m*bn]; |
|
|
|
std::memcpy(clone_b, b, m*bn*sizeof(MKL_Complex16)); |
|
|
|
@ -859,9 +859,9 @@ extern "C"{ |
|
|
|
MKL_Complex16 one = {1.0, 0.0}; |
|
|
|
cblas_ztrsm(CblasColMajor, CblasLeft, CblasUpper, CblasNoTrans, CblasNonUnit, n, bn, &one, r, m, clone_b, m); |
|
|
|
|
|
|
|
for (int i = 0; i < n; ++i) |
|
|
|
for (MKL_INT i = 0; i < n; ++i) |
|
|
|
{ |
|
|
|
for (int j = 0; j < bn; ++j) |
|
|
|
for (MKL_INT j = 0; j < bn; ++j) |
|
|
|
{ |
|
|
|
x[j * n + i] = clone_b[j * m + i]; |
|
|
|
} |
|
|
|
@ -871,32 +871,32 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int s_svd_factor(bool compute_vectors, int m, int n, float a[], float s[], float u[], float v[], float work[], int len) |
|
|
|
DLLEXPORT MKL_INT s_svd_factor(bool compute_vectors, MKL_INT m, MKL_INT n, float a[], float s[], float u[], float v[], float work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
char job = compute_vectors ? 'A' : 'N'; |
|
|
|
sgesvd_(&job, &job, &m, &n, a, &m, s, u, &m, v, &n, work, &len, &info); |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int d_svd_factor(bool compute_vectors, int m, int n, double a[], double s[], double u[], double v[], double work[], int len) |
|
|
|
DLLEXPORT MKL_INT d_svd_factor(bool compute_vectors, MKL_INT m, MKL_INT n, double a[], double s[], double u[], double v[], double work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
MKL_INT info = 0; |
|
|
|
char job = compute_vectors ? 'A' : 'N'; |
|
|
|
dgesvd_(&job, &job, &m, &n, a, &m, s, u, &m, v, &n, work, &len, &info); |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int c_svd_factor(bool compute_vectors, int m, int n, MKL_Complex8 a[], MKL_Complex8 s[], MKL_Complex8 u[], MKL_Complex8 v[], MKL_Complex8 work[], int len) |
|
|
|
DLLEXPORT MKL_INT c_svd_factor(bool compute_vectors, MKL_INT m, MKL_INT n, MKL_Complex8 a[], MKL_Complex8 s[], MKL_Complex8 u[], MKL_Complex8 v[], MKL_Complex8 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
int dim_s = std::min(m,n); |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_INT dim_s = std::min(m,n); |
|
|
|
float* rwork = new float[5 * dim_s]; |
|
|
|
float* s_local = new float[dim_s]; |
|
|
|
char job = compute_vectors ? 'A' : 'N'; |
|
|
|
cgesvd_(&job, &job, &m, &n, a, &m, s_local, u, &m, v, &n, work, &len, rwork, &info); |
|
|
|
|
|
|
|
for(int index = 0; index < dim_s; ++index){ |
|
|
|
for(MKL_INT index = 0; index < dim_s; ++index){ |
|
|
|
MKL_Complex8 value = {s_local[index], 0.0f}; |
|
|
|
s[index] = value; |
|
|
|
} |
|
|
|
@ -906,16 +906,16 @@ extern "C"{ |
|
|
|
return info; |
|
|
|
} |
|
|
|
|
|
|
|
DLLEXPORT int z_svd_factor(bool compute_vectors, int m, int n, MKL_Complex16 a[], MKL_Complex16 s[], MKL_Complex16 u[], MKL_Complex16 v[], MKL_Complex16 work[], int len) |
|
|
|
DLLEXPORT MKL_INT z_svd_factor(bool compute_vectors, MKL_INT m, MKL_INT n, MKL_Complex16 a[], MKL_Complex16 s[], MKL_Complex16 u[], MKL_Complex16 v[], MKL_Complex16 work[], MKL_INT len) |
|
|
|
{ |
|
|
|
int info = 0; |
|
|
|
int dim_s = std::min(m,n); |
|
|
|
MKL_INT info = 0; |
|
|
|
MKL_INT dim_s = std::min(m,n); |
|
|
|
double* rwork = new double[5 * std::min(m, n)]; |
|
|
|
double* s_local = new double[dim_s]; |
|
|
|
char job = compute_vectors ? 'A' : 'N'; |
|
|
|
zgesvd_(&job, &job, &m, &n, a, &m, s_local, u, &m, v, &n, work, &len, rwork, &info); |
|
|
|
|
|
|
|
for(int index = 0; index < dim_s; ++index){ |
|
|
|
for(MKL_INT index = 0; index < dim_s; ++index){ |
|
|
|
MKL_Complex16 value = {s_local[index], 0.0f}; |
|
|
|
s[index] = value; |
|
|
|
} |
|
|
|
|