Commit dff0a8bb authored by Pavel Kus's avatar Pavel Kus

single/double unification in check_correctness

of _GEMM and P_GEMM calls
parent fe6f28a3
......@@ -35,6 +35,7 @@
#undef cublas_PRECISION_TRMM
#undef cublas_PRECISION_GEMV
#undef cublas_PRECISION_SYMV
#undef scal_PRECISION_GEMM
#undef PRECISION_SUFFIX
#undef CONST_0_0
......@@ -87,6 +88,7 @@
#define cublas_PRECISION_TRMM cublas_DTRMM
#define cublas_PRECISION_GEMV cublas_DGEMV
#define cublas_PRECISION_SYMV cublas_DSYMV
#define scal_PRECISION_GEMM PDGEMM
#define CONST_0_0 0.0_rk8
#define CONST_0_5 0.5_rk8
#define CONST_1_0 1.0_rk8
......@@ -136,6 +138,7 @@
#define cublas_PRECISION_TRMM cublas_STRMM
#define cublas_PRECISION_GEMV cublas_SGEMV
#define cublas_PRECISION_SYMV cublas_SSYMV
#define scal_PRECISION_GEMM PSGEMM
#define CONST_0_0 0.0_rk4
#define CONST_0_5 0.5_rk4
#define CONST_1_0 1.0_rk4
......@@ -189,6 +192,7 @@
#undef cublas_PRECISION_TRMM
#undef cublas_PRECISION_GEMV
#undef cublas_PRECISION_SYMV
#undef scal_PRECISION_GEMM
#undef PRECISION_SUFFIX
#undef MPI_COMPLEX_PRECISION
#undef MPI_MATH_DATATYPE_PRECISION
......@@ -251,6 +255,7 @@
#define cublas_PRECISION_TRMM cublas_ZTRMM
#define cublas_PRECISION_GEMV cublas_ZGEMV
#define cublas_PRECISION_SYMV cublas_ZSYMV
#define scal_PRECISION_GEMM PZGEMM
#define MPI_COMPLEX_PRECISION MPI_DOUBLE_COMPLEX
#define MPI_MATH_DATATYPE_PRECISION MPI_DOUBLE_COMPLEX
#define MPI_COMPLEX_EXPLICIT_PRECISION MPI_COMPLEX16
......@@ -308,6 +313,7 @@
#define cublas_PRECISION_TRMM cublas_CTRMM
#define cublas_PRECISION_GEMV cublas_CGEMV
#define cublas_PRECISION_SYMV cublas_CSYMV
#define scal_PRECISION_GEMM PCGEMM
#define MPI_COMPLEX_PRECISION MPI_COMPLEX
#define MPI_MATH_DATATYPE_PRECISION MPI_COMPLEX
#define MPI_COMPLEX_EXPLICIT_PRECISION MPI_COMPLEX8
......
......@@ -104,45 +104,19 @@
#if REALCASE == 1
#ifdef WITH_MPI
#ifdef DOUBLE_PRECISION_REAL
call pdgemm('N', 'N', na, nev, na, 1.0_rk8, as, 1, 1, sc_desc, &
z, 1, 1, sc_desc, 0.0_rk8, tmp1, 1, 1, sc_desc)
#else
call psgemm('N', 'N', na, nev, na, 1.0_rk4, as, 1, 1, sc_desc, &
z, 1, 1, sc_desc, 0.0_rk4, tmp1, 1, 1, sc_desc)
#endif
call scal_PRECISION_GEMM('N', 'N', na, nev, na, CONST_1_0, as, 1, 1, sc_desc, &
z, 1, 1, sc_desc, CONST_0_0, tmp1, 1, 1, sc_desc)
#else /* WITH_MPI */
#ifdef DOUBLE_PRECISION_REAL
call dgemm('N','N',na,nev,na,1.0_rk8,as,na,z,na,0.0_rk8,tmp1,na)
#else
call sgemm('N','N',na,nev,na,1.0_rk4,as,na,z,na,0.0_rk4,tmp1,na)
#endif
call PRECISION_GEMM('N','N',na,nev,na,CONST_1_0,as,na,z,na,CONST_0_0,tmp1,na)
#endif /* WITH_MPI */
#endif /* REALCASE */
#if COMPLEXCASE == 1
#ifdef WITH_MPI
#ifdef DOUBLE_PRECISION_COMPLEX
call pzgemm('N', 'N', na, nev, na, CONE, as, 1, 1, sc_desc, &
z, 1, 1, sc_desc, CZERO, tmp1, 1, 1, sc_desc)
#else
call pcgemm('N', 'N', na, nev, na, CONE, as, 1, 1, sc_desc, &
call scal_PRECISION_GEMM('N', 'N', na, nev, na, CONE, as, 1, 1, sc_desc, &
z, 1, 1, sc_desc, CZERO, tmp1, 1, 1, sc_desc)
#endif
#else /* WITH_MPI */
#ifdef DOUBLE_PRECISION_COMPLEX
call zgemm('N','N',na,nev,na,CONE,as,na,z,na,CZERO,tmp1,na)
#else
call cgemm('N','N',na,nev,na,CONE,as,na,z,na,CZERO,tmp1,na)
#endif
call PRECISION_GEMM('N','N',na,nev,na,CONE,as,na,z,na,CZERO,tmp1,na)
#endif /* WITH_MPI */
#endif /* COMPLEXCASE */
......@@ -318,48 +292,21 @@
#if REALCASE == 1
#ifdef WITH_MPI
#ifdef DOUBLE_PRECISION_REAL
call pdgemm('T', 'N', nev, nev, na, 1.0_rk8, z, 1, 1, sc_desc, &
z, 1, 1, sc_desc, 0.0_rk8, tmp1, 1, 1, sc_desc)
#else
call psgemm('T', 'N', nev, nev, na, 1.0_rk4, z, 1, 1, sc_desc, &
z, 1, 1, sc_desc, 0.0_rk4, tmp1, 1, 1, sc_desc)
#endif
call scal_PRECISION_GEMM('T', 'N', nev, nev, na, CONST_1_0, z, 1, 1, sc_desc, &
z, 1, 1, sc_desc, CONST_0_0, tmp1, 1, 1, sc_desc)
#else /* WITH_MPI */
#ifdef DOUBLE_PRECISION_REAL
call dgemm('T','N',nev,nev,na,1.0_rk8,z,na, &
z,na,0.0_rk8,tmp1,na)
#else
call sgemm('T','N',nev,nev,na,1.0_rk4,z,na, &
z,na,0.0_rk4,tmp1,na)
#endif
call PRECISION_GEMM('T','N',nev,nev,na,CONST_1_0,z,na, &
z,na,CONST_0_0,tmp1,na)
#endif /* WITH_MPI */
#endif /* REALCASE */
#if COMPLEXCASE == 1
#ifdef WITH_MPI
#ifdef DOUBLE_PRECISION_COMPLEX
call pzgemm('C', 'N', nev, nev, na, CONE, z, 1, 1, sc_desc, &
z, 1, 1, sc_desc, CZERO, tmp1, 1, 1, sc_desc)
#else
call pcgemm('C', 'N', nev, nev, na, CONE, z, 1, 1, sc_desc, &
call scal_PRECISION_GEMM('C', 'N', nev, nev, na, CONE, z, 1, 1, sc_desc, &
z, 1, 1, sc_desc, CZERO, tmp1, 1, 1, sc_desc)
#endif
#else /* WITH_MPI */
#ifdef DOUBLE_PRECISION_COMPLEX
call zgemm('C','N',nev,nev,na,CONE,z,na,z,na,CZERO,tmp1,na)
#else
call cgemm('C','N',nev,nev,na,CONE,z,na,z,na,CZERO,tmp1,na)
#endif
call PRECISION_GEMM('C','N',nev,nev,na,CONE,z,na,z,na,CZERO,tmp1,na)
#endif /* WITH_MPI */
#endif /* COMPLEXCASE */
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment