我尝试在matlab中编写使用magma库,所以基本上我编写了一个mexfunction,它使用magma函数合并了c代码,然后将此mexfunction编译成mexa64文件,因此我可以在matlab中使用。
mexfunction或源c代码如下:(称为eig_magma)
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <cuda_runtime_api.h>
#include <cublas.h>
// includes, project
#include "flops.h"
#include "magma.h"
#include "magma_lapack.h"
#include "testings.h"
#include <math.h>
#include "mex.h"
#define A(i,j) A[i + j*lda]
extern "C"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
#define L_OUT plhs[0]
#define A_IN prhs[0]
#define S_OUT plhs[1]
magma_init();
real_Double_t gpu_perf, gpu_time;
double *h_A, *h_R, *h_work, *L,*S;
double *w;
magma_int_t *iwork;
magma_int_t N, n2, info, lwork, liwork, lda, aux_iwork[1];
double aux_work[1];
magma_vec_t jobz = MagmaVec;
magma_uplo_t uplo = MagmaLower;
magma_int_t i,j;
N = mxGetM(A_IN);
lda = N;
n2 = lda*N;
// query for workspace sizes
magma_dsyevd( jobz, uplo,
N, NULL, lda, NULL,
aux_work, -1,
aux_iwork, -1,
&info );
lwork = (magma_int_t) aux_work[0];
liwork = aux_iwork[0];
TESTING_MALLOC_CPU( h_A, double, n2);
h_A = (double *)mxGetData(A_IN);
L_OUT = mxCreateDoubleMatrix(N,N,mxREAL);
L = mxGetPr(L_OUT);
S_OUT = mxCreateDoubleMatrix(N,1,mxREAL);
S = mxGetPr(S_OUT);
//print and check
printf("print out h_A\n");
for(i=0;i<N;i++)
{
for(j=0;j<N;j++)
printf("%f\t",h_A[i+j*lda]);
printf("\n");
}
/* Allocate host memory for the matrix */
TESTING_MALLOC_CPU( w, double, N );
TESTING_MALLOC_CPU( iwork, magma_int_t, liwork );
TESTING_MALLOC_PIN( h_R, double, N*lda );
TESTING_MALLOC_PIN( h_work, double, lwork );
printf("allocate memory on cpu with pin!\n");
printf("copy h_A to h_R, and then print h_R\n");
memcpy( h_R, h_A, sizeof(double)*n2);
// lapackf77_dlacpy( MagmaUpperLowerStr, &N, &N, h_A, &lda,h_R, &lda );
// print and check
for(i=0;i<N;i++)
{
for(j=0;j<N;j++)
{
printf("%f\t",h_R[i+j*lda]);
}
printf("\n");
}
/* Performs operation using MAGA */
gpu_time = magma_wtime();
printf("start eig\n");
magma_dsyevd( jobz, uplo,
N, h_R, lda, w,
h_work, lwork,
iwork, liwork,
&info );
gpu_time = magma_wtime() - gpu_time;
printf("%d\n",info);
if (info != 0)
printf("magma_dsyevd returned error %d: %s.\n",
(int) info, magma_strerror( info ));
printf("convey the output\n");
memcpy(L,h_R,sizeof(double)*n2);
memcpy(S,w,sizeof(double)*N);
TESTING_FREE_CPU( w );
TESTING_FREE_CPU( iwork );
TESTING_FREE_PIN( h_R );
TESTING_FREE_PIN( h_work);
TESTING_FINALIZE();
}
因为使用 magma 需要 cuda 和 magma lib 我编写了 makefile 来将代码编译成 mex 文件:
# Paths where MAGMA, CUDA, and OpenBLAS are installed
MAGMADIR = /home/yuxin/magma-1.5.0
CUDADIR = /usr/local/cuda
MATLABHOME = /opt/share/MATLAB/R2012b.app/
CC = icpc
LD = icpc
CFLAGS = -Wall
LDFLAGS = -Wall -mkl=parallel
MEX_CFLAGS = -fPIC -ansi -pthread -DMX_COMPAT_32 \
-DMATLAB_MEX_FILE
MAGMA_CFLAGS := -DADD_ -DHAVE_CUBLAS -I$(MAGMADIR)/include -I$(CUDADIR)/include -I$(MAGMADIR)/testing
MAGMA_LIBS := -L$(MAGMADIR)/lib -L$(CUDADIR)/lib64 \
-lmagma -lcublas -lcudart
MEXLIBS := -L/opt/share/MATLAB/R2012b.app/bin/glnxa64 -lmx -lmex -lmat -lm -lstdc++
MEX_INCLU := -I$(MATLABHOME)/extern/include -I$(MATLABHOME)/simulink/include
MEXFLAGS := -shared
OBJECT = eig_magma.o
EXECUTABLE = eig_magma
$(EXECUTABLE): $(OBJECT)
$(CC) $(LDFLAGS) $(MEXFLAGS) $(OBJECT) -o $@ $(MAGMA_LIBS) $(MEXLIBS)
eig_magma.o: eig_magma.c
$(CC) $(CFLAGS) $(MAGMA_CFLAGS) $(MEX_INCLU) $(MEX_CFLAGS) -c $< -o $@
clean:
rm -rf eig_magma *.o eig_magma.mexa64
它可以编译成功,但是当我在 matlab 中运行 eig_magma 时,当 matlab 执行到
magma_dsyevd(jobz,uplo,N,h_R,lda,w,h_work,lwork,iwork,liwork,&info);
matlab 崩溃了......我也尝试过只用c 编写eig_magma,不使用mex 函数,它编译成功并且工作正常。
任何人都知道如何解决这个问题?
谢谢
玉心