/* testGMMM.c */
#include "../InpMtx.h"
#include "../../Drand.h"
/*--------------------------------------------------------------------*/
int
main ( int argc, char *argv[] )
/*
------------------------------------------------------
generate a random matrix and test the InpMtx_*_gmmm*()
matrix-matrix multiply methods.
the output is a matlab file to test correctness.
created -- 98nov14, cca
------------------------------------------------------
*/
{
DenseMtx *X, *Y ;
double alpha[2], beta[2] ;
double alphaImag, alphaReal, betaImag, betaReal ;
Drand *drand ;
int dataType, msglvl, ncolA, nitem, nrhs, nrowA, nrowX,
nrowY, seed, coordType, rc, symflag, transposeflag ;
InpMtx *A ;
FILE *msgFile ;
if ( argc != 16 ) {
fprintf(stdout,
"\n\n %% usage : %s msglvl msgFile symflag coordType transpose"
"\n %% nrow ncol nent nrhs seed "
"\n %% alphaReal alphaImag betaReal betaImag"
"\n %% msglvl -- message level"
"\n %% msgFile -- message file"
"\n %% dataType -- type of matrix entries"
"\n %% 1 -- real"
"\n %% 2 -- complex"
"\n %% symflag -- symmetry flag"
"\n %% 0 -- symmetric"
"\n %% 1 -- hermitian"
"\n %% 2 -- nonsymmetric"
"\n %% coordType -- storage mode"
"\n %% 1 -- by rows"
"\n %% 2 -- by columns"
"\n %% 3 -- by chevrons, (requires nrow = ncol)"
"\n %% transpose -- transpose flag"
"\n %% 0 -- Y := beta * Y + alpha * A * X"
"\n %% 1 -- Y := beta * Y + alpha * A^H * X, nonsymmetric only"
"\n %% 2 -- Y := beta * Y + alpha * A^T * X, nonsymmetric only"
"\n %% nrowA -- number of rows in A"
"\n %% ncolA -- number of columns in A"
"\n %% nitem -- number of items"
"\n %% nrhs -- number of right hand sides"
"\n %% seed -- random number seed"
"\n %% alphaReal -- y := beta*y + alpha*A*x"
"\n %% alphaImag -- y := beta*y + alpha*A*x"
"\n %% betaReal -- y := beta*y + alpha*A*x"
"\n %% betaImag -- y := beta*y + alpha*A*x"
"\n", argv[0]) ;
return(0) ;
}
msglvl = atoi(argv[1]) ;
if ( strcmp(argv[2], "stdout") == 0 ) {
msgFile = stdout ;
} else if ( (msgFile = fopen(argv[2], "a")) == NULL ) {
fprintf(stderr, "\n fatal error in %s"
"\n unable to open file %s\n",
argv[0], argv[2]) ;
return(-1) ;
}
dataType = atoi(argv[3]) ;
symflag = atoi(argv[4]) ;
coordType = atoi(argv[5]) ;
transposeflag = atoi(argv[6]) ;
nrowA = atoi(argv[7]) ;
ncolA = atoi(argv[8]) ;
nitem = atoi(argv[9]) ;
nrhs = atoi(argv[10]) ;
seed = atoi(argv[11]) ;
alphaReal = atof(argv[12]) ;
alphaImag = atof(argv[13]) ;
betaReal = atof(argv[14]) ;
betaImag = atof(argv[15]) ;
fprintf(msgFile,
"\n %% %s "
"\n %% msglvl -- %d"
"\n %% msgFile -- %s"
"\n %% dataType -- %d"
"\n %% symflag -- %d"
"\n %% coordType -- %d"
"\n %% transposeflag -- %d"
"\n %% nrowA -- %d"
"\n %% ncolA -- %d"
"\n %% nitem -- %d"
"\n %% nrhs -- %d"
"\n %% seed -- %d"
"\n %% alphaReal -- %e"
"\n %% alphaImag -- %e"
"\n %% betaReal -- %e"
"\n %% betaImag -- %e"
"\n",
argv[0], msglvl, argv[2], dataType, symflag, coordType,
transposeflag, nrowA, ncolA, nitem, nrhs, seed,
alphaReal, alphaImag, betaReal, betaImag) ;
fflush(msgFile) ;
if ( dataType != 1 && dataType != 2 ) {
fprintf(stderr, "\n invalid value %d for dataType\n", dataType) ;
exit(-1) ;
}
if ( symflag != 0 && symflag != 1 && symflag != 2 ) {
fprintf(stderr, "\n invalid value %d for symflag\n", symflag) ;
exit(-1) ;
}
if ( coordType != 1 && coordType != 2 && coordType != 3 ) {
fprintf(stderr,
"\n invalid value %d for coordType\n", coordType) ;
exit(-1) ;
}
if ( transposeflag < 0
|| transposeflag > 2 ) {
fprintf(stderr, "\n error, transposeflag = %d, must be 0, 1 or 2",
transposeflag) ;
exit(-1) ;
}
if ( (transposeflag == 1 && symflag != 2)
|| (transposeflag == 2 && symflag != 2) ) {
fprintf(stderr, "\n error, transposeflag = %d, symflag = %d",
transposeflag, symflag) ;
exit(-1) ;
}
if ( transposeflag == 1 && dataType != 2 ) {
fprintf(stderr, "\n error, transposeflag = %d, dataType = %d",
transposeflag, dataType) ;
exit(-1) ;
}
if ( symflag == 1 && dataType != 2 ) {
fprintf(stderr,
"\n symflag = 1 (hermitian), dataType != 2 (complex)") ;
exit(-1) ;
}
if ( nrowA <= 0 || ncolA <= 0 || nitem <= 0 ) {
fprintf(stderr,
"\n invalid value: nrow = %d, ncol = %d, nitem = %d",
nrowA, ncolA, nitem) ;
exit(-1) ;
}
if ( symflag < 2 && nrowA != ncolA ) {
fprintf(stderr,
"\n invalid data: symflag = %d, nrow = %d, ncol = %d",
symflag, nrowA, ncolA) ;
exit(-1) ;
}
alpha[0] = alphaReal ;
alpha[1] = alphaImag ;
beta[0] = betaReal ;
beta[1] = betaImag ;
drand = Drand_new() ;
Drand_setSeed(drand, seed) ;
Drand_setUniform(drand, -1.0, 1.0) ;
/*
----------------------------
initialize the matrix object
and fill with random entries
----------------------------
*/
A = InpMtx_new() ;
InpMtx_init(A, coordType, dataType, 0, 0) ;
rc = InpMtx_randomMatrix(A, dataType, coordType, INPMTX_BY_VECTORS,
nrowA, ncolA, symflag, 1, nitem, seed) ;
if ( rc != 1 ) {
fprintf(stderr, "\n error return %d from InpMtx_randomMatrix()", rc);
exit(-1) ;
}
/*
-------------------------------------------
write the assembled matrix to a matlab file
-------------------------------------------
*/
InpMtx_writeForMatlab(A, "A", msgFile) ;
if ( symflag == 0 ) {
fprintf(msgFile,
"\n for k = 1:%d"
"\n for j = k+1:%d"
"\n A(j,k) = A(k,j) ;"
"\n end"
"\n end", nrowA, ncolA) ;
} else if ( symflag == 1 ) {
fprintf(msgFile,
"\n for k = 1:%d"
"\n for j = k+1:%d"
"\n A(j,k) = ctranspose(A(k,j)) ;"
"\n end"
"\n end", nrowA, ncolA) ;
}
/*
-------------------------------
generate dense matrices X and Y
-------------------------------
*/
if ( transposeflag == 0 ) {
nrowX = ncolA ;
nrowY = nrowA ;
} else {
nrowX = nrowA ;
nrowY = ncolA ;
}
X = DenseMtx_new() ;
Y = DenseMtx_new() ;
DenseMtx_init(X, dataType, 0, 0, nrowX, nrhs, 1, nrowX) ;
DenseMtx_fillRandomEntries(X, drand) ;
DenseMtx_init(Y, dataType, 0, 0, nrowY, nrhs, 1, nrowY) ;
DenseMtx_fillRandomEntries(Y, drand) ;
fprintf(msgFile, "\n X = zeros(%d,%d) ;", nrowX, nrhs) ;
DenseMtx_writeForMatlab(X, "X", msgFile) ;
fprintf(msgFile, "\n Y = zeros(%d,%d) ;", nrowY, nrhs) ;
DenseMtx_writeForMatlab(Y, "Y", msgFile) ;
/*
----------------------------------
perform the matrix-matrix multiply
----------------------------------
*/
fprintf(msgFile, "\n beta = %20.12e + %20.2e*i;", beta[0], beta[1]);
fprintf(msgFile, "\n alpha = %20.12e + %20.2e*i;", alpha[0], alpha[1]);
fprintf(msgFile, "\n Z = zeros(%d,1) ;", nrowY) ;
if ( transposeflag == 0 ) {
if ( symflag == 0 ) {
InpMtx_sym_gmmm(A, beta, Y, alpha, X) ;
} else if ( symflag == 1 ) {
InpMtx_herm_gmmm(A, beta, Y, alpha, X) ;
} else if ( symflag == 2 ) {
InpMtx_nonsym_gmmm(A, beta, Y, alpha, X) ;
}
DenseMtx_writeForMatlab(Y, "Z", msgFile) ;
fprintf(msgFile, "\n maxerr = max(Z - beta*Y - alpha*A*X) ") ;
fprintf(msgFile, "\n") ;
} else if ( transposeflag == 1 ) {
InpMtx_nonsym_gmmm_H(A, beta, Y, alpha, X) ;
DenseMtx_writeForMatlab(Y, "Z", msgFile) ;
fprintf(msgFile,
"\n maxerr = max(Z - beta*Y - alpha*ctranspose(A)*X) ") ;
fprintf(msgFile, "\n") ;
} else if ( transposeflag == 2 ) {
InpMtx_nonsym_gmmm_T(A, beta, Y, alpha, X) ;
DenseMtx_writeForMatlab(Y, "Z", msgFile) ;
fprintf(msgFile,
"\n maxerr = max(Z - beta*Y - alpha*transpose(A)*X) ") ;
fprintf(msgFile, "\n") ;
}
/*
------------------------
free the working storage
------------------------
*/
InpMtx_free(A) ;
DenseMtx_free(X) ;
DenseMtx_free(Y) ;
Drand_free(drand) ;
fclose(msgFile) ;
return(1) ; }
/*--------------------------------------------------------------------*/
syntax highlighted by Code2HTML, v. 0.9.1