From 04bd16aeb391e67716344268cf0f43d1f31f180a Mon Sep 17 00:00:00 2001 From: cnlohr Date: Sat, 17 Mar 2018 14:21:20 -0400 Subject: Update dcl and test. --- redist/dclapack.h | 20 ++++++++++++++++---- redist/dclhelpers.c | 14 +++++++++++--- redist/dclhelpers.h | 8 +++++++- redist/test_dcl.c | 41 ++++++++++++++++++++++++++--------------- 4 files changed, 60 insertions(+), 23 deletions(-) (limited to 'redist') diff --git a/redist/dclapack.h b/redist/dclapack.h index 2823c5d..cae377b 100644 --- a/redist/dclapack.h +++ b/redist/dclapack.h @@ -37,9 +37,9 @@ /* * Returns the identity matrix (with size n x n, but width Ic) */ -#define IDENTITY(I, n) \ +#define IDENTITY(I, m, n) \ { \ - for (int _i = 0; _i < (n); _i++) { \ + for (int _i = 0; _i < (m); _i++) { \ for (int _j = 0; _j < _i; _j++) { \ _(I, _i, _j) = 0.0f; \ } \ @@ -50,6 +50,17 @@ } \ } + +/* + * Returns the identity matrix (with size n x n, but width Ic) + */ +#define ZERO(Z, m, n) \ + { \ + for (int _i = 0; _i < (m); _i++) \ + for (int _j = 0; _j < (n); _j++) \ + _(Z, _i, _j) = 0.0f; \ + } + /* * R = Transpose(A) * A is (m by n) @@ -79,7 +90,7 @@ _(U, _i, _j) = _(A, _i, _j); \ } \ } \ - IDENTITY(L, n); \ + IDENTITY(L, n, n); \ \ for (_i = 0; _i < (n)-1; _i++) { \ \ @@ -188,7 +199,7 @@ #define INV(Ainv,A,n,ORDER) { \ INV_SETUP(ORDER) \ int Piv[ORDER]; \ - IDENTITY(I,n); \ + IDENTITY(I,n,n); \ LU(L,U,A,Piv,n); \ PIVOT(Ipiv,I,Piv,n,n); \ L_SUB(C,L,Ipiv,n,n); \ @@ -258,6 +269,7 @@ PRINT(Ainv,n,n); \ for (int _k = 0; _k < n; _k++) { \ sum += _(A, _i, _k) * _(B, _k, _j); \ } \ + printf( "-> %d %d = %f %f %f %f[%d %d]\n", _i, _j, alpha,sum,beta, _(C,_i,_j), _i, _j ); \ _(R, _i, _j) = alpha * sum + beta * _(C, _i, _j); \ } \ } \ diff --git a/redist/dclhelpers.c b/redist/dclhelpers.c index 3e51fd2..5b956b0 100644 --- a/redist/dclhelpers.c +++ b/redist/dclhelpers.c @@ -5,15 +5,23 @@ #include #include #include +#include void dclPrint( const DCL_FLOAT * PMATRIX, int PMATRIXc, int n, int m ) { PRINT( PMATRIX, n, m ); } -void dclIdentity( DCL_FLOAT * I, int Ic, int n ) +void dclIdentity( DCL_FLOAT * I, int Ic, int m, int n ) { - IDENTITY( I, n ); + IDENTITY( I, m, n ); +} + + +/* Returns the zero matrix */ +void dclZero( DCL_FLOAT * Z, int Zc, int m, int n ) +{ + memset( Z, 0, m*n*sizeof(DCL_FLOAT) ); } void dclTransp( DCL_FLOAT * R, int Rc, const DCL_FLOAT * A, int Ac, int n, int m ) @@ -103,7 +111,7 @@ void dcldgemm( } else tb = B; - printf("%d %d %d\n", tac, tbc, Cc); + GMULADD(C, ta, tb, C, alpha, beta, m, n, k); } diff --git a/redist/dclhelpers.h b/redist/dclhelpers.h index 8779e1a..b4acec8 100644 --- a/redist/dclhelpers.h +++ b/redist/dclhelpers.h @@ -3,11 +3,17 @@ #define DCL_FLOAT FLT +//Use this macro to safely +#define DMS( m ) ((m)[0]), (sizeof((m)[0])/sizeof((m)[0][0])) + /* Prints matrix A of size[n][m] */ void dclPrint( const DCL_FLOAT * A, int Ac, int n, int m ); /* Returns the identity matrix */ -void dclIdentity( DCL_FLOAT * I, int Ic, int n ); +void dclIdentity( DCL_FLOAT * I, int Ic, int m, int n ); + +/* Returns the zero matrix */ +void dclZero( DCL_FLOAT * I, int Ic, int m, int n ); /* R = Transpose(A) A is (n by m) diff --git a/redist/test_dcl.c b/redist/test_dcl.c index 056acdf..6b0d870 100644 --- a/redist/test_dcl.c +++ b/redist/test_dcl.c @@ -12,7 +12,7 @@ int main() dclTransp( B[0], 2, A[0], 4, 2, 4 ); dclPrint( B[0], 2, 4, 2 ); - int i; + int i, j; for( i = 0; i < 8; i++ ) { printf( "%f\n", ((float*)(B[0]))[i] ); @@ -34,54 +34,65 @@ int main() { FLT A[3][4]; - dclIdentity(A[0], 4, 3); - dclPrint(A[0], 4, 3, 4); + dclIdentity( DMS(A), 3, 4); + dclPrint( DMS(A), 3, 4); FLT x[4][2] = { {7, -7}, {8, -8}, {9, -9}, {10, -10}, }; FLT R[4][2]; + dclZero( DMS(R), 4, 2 ); - printf("%p %p %p\n", A, x, R); // dclMul(R, 1, A[0], 4, x, 1, 4, 1, 3); dcldgemm(0, 0, 3, 4, 2, 1, A[0], 4, x[0], 2, 0, R[0], 2); - dclPrint(x[0], 2, 4, 2); - dclPrint(R[0], 2, 4, 2); + dclPrint(DMS(x), 4, 2); + dclPrint(DMS(R), 3, 2); for (int j = 0; j < 2; j++) { for (int i = 0; i < 3; i++) + { + printf( "[%d][%d]\n", i, j ); assert(R[i][j] == x[i][j]); + } assert(fabs(R[3][j]) < .0000001); } } - printf( "The following should be an identity matrix\n" ); - dclPrint( MM[0], 3, 3, 3 ); - +#if 0 //Currently failing test... { +// FLT em1[3][4]; +// FLT em2[4][2]; +// FLT emo[4][2]; + FLT em1[12][20]; FLT em2[20][12]; - FLT emo[12][12]; + FLT emo[20][12]; int x, y; for( y = 0; y < 12; y++ ) for( x = 0; x < 20; x++ ) - { em1[y][x] = (rand()%1000)/1000.0; - em2[x][y] = (rand()%1000)/1000.0; - } + + for( y = 0; y < 20; y++ ) + for( x = 0; x < 12; x++ ) + em2[y][x] = (rand()%1000)/1000.0; int m = 12; int n = 20; int k = 12; - dcldgemm( 0, 0, m, n, k, 0.1, em1[0], 20, em2[0], 12, .1, emo[0], 12 ); - dclPrint( emo[0], 12, 12, 12 ); + + dclPrint( em1[0], 20, 12, 20 ); + dclPrint( em2[0], 12, 20, 12 ); + + dcldgemm( 0, 0, m, n, k, 1.0, em1[0], 20, em2[0], 12, .1, emo[0], 12 ); + dclPrint( emo[0], 12, 20, 12 ); } +#endif } -- cgit v1.2.3