Skip to content

Commit

Permalink
Still debugging multiple processes
Browse files Browse the repository at this point in the history
  • Loading branch information
@Heatherms27
Heatherms27 committed Jan 2, 2024
1 parent 775a1d8 commit 11a2166
Show file tree
Hide file tree
Showing 6 changed files with 488 additions and 332 deletions.
9 changes: 5 additions & 4 deletions examples/ex_eigs_dpar.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -98,8 +98,9 @@ int main (int argc, char *argv[]) {
primme.aNorm = 1.0;
primme.printLevel = 4;

primme.maxBlockSize = 2;
primme.expansionParams.expansion = primme_expansion_fullLanczos;
primme.maxBlockSize = 1;
//primme.expansionParams.expansion = primme_expansion_fullLanczos;
primme.expansionParams.expansion = primme_expansion_davidson;
/* Set method to solve the problem */
//primme_set_method(PRIMME_DYNAMIC, &primme);
primme_set_method(PRIMME_DEFAULT_MIN_MATVECS, &primme);
Expand Down Expand Up @@ -147,8 +148,8 @@ int main (int argc, char *argv[]) {
evals[i], rnorms[i]);
}
fprintf(primme.outputFile, " %d eigenpairs converged\n", primme.initSize);
fprintf(primme.outputFile, "Tolerance : %-22.15E\n",
primme.aNorm*primme.eps);
fprintf(primme.outputFile, "Tolerance : %-22.15E\t (%-22.4E tolerance x %-22.4E norm)\n",
primme.aNorm*primme.eps, primme.eps, primme.aNorm);
fprintf(primme.outputFile, "Iterations: %-" PRIMME_INT_P "\n",
primme.stats.numOuterIterations);
fprintf(primme.outputFile, "Restarts : %-" PRIMME_INT_P "\n", primme.stats.numRestarts);
Expand Down
20 changes: 12 additions & 8 deletions examples/ex_eigs_dparTests.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -117,13 +117,6 @@ int main (int argc, char *argv[]) {
* n % numProcs */
primme.globalSumReal = par_GlobalSum;

/* Set problem matrix */
ierr = CSR_to_PETSc_Matrix(argv[1], &A, &primme, &ierr); CHKERRQ(ierr);
primme.matrix = &A;
primme.matrixMatvec = PETScMatvec;
/* Function that implements the matrix-vector product
A*x for solving the problem A*x = l*x */

/* Default Parameters */
primme.maxBasisSize = 100;
primme.numEvals = 10; /* Number of wanted eigenpairs */
Expand All @@ -134,8 +127,10 @@ int main (int argc, char *argv[]) {
primme.eps = 1e-1; /* ||r|| <= eps * ||matrix|| */
primme.aNorm = 1.0;
primme.target = primme_largest;
//primme.maxMatvecs = 1000;
primme.maxMatvecs = 1000;
primme.maxOuterIterations = 1000;

primme.dynamicMethodSwitch = 0;
primme.expansionParams.expansion = primme_expansion_davidson;
primme.projectionParams.projection = primme_proj_default;

Expand All @@ -147,6 +142,8 @@ int main (int argc, char *argv[]) {
return 0;
}
if(strcmp(argv[i], "-basisSize") == 0) primme.maxBasisSize = atoi(argv[i+1]);
if(strcmp(argv[i], "-maxIters") == 0) primme.maxOuterIterations = atoi(argv[i+1]);
if(strcmp(argv[i], "-maxMatvecs") == 0) primme.maxMatvecs = atoi(argv[i+1]);
if(strcmp(argv[i], "-numEvals") == 0) primme.numEvals = atoi(argv[i+1]);
if(strcmp(argv[i], "-blockSize") == 0) primme.maxBlockSize = atoi(argv[i+1]);
if(strcmp(argv[i], "-printLevel") == 0) primme.printLevel = atoi(argv[i+1]);
Expand Down Expand Up @@ -193,6 +190,13 @@ int main (int argc, char *argv[]) {
} /* End residual */
} /* End command line argument options */

/* Set problem matrix */
ierr = CSR_to_PETSc_Matrix(argv[1], &A, &primme, &ierr); CHKERRQ(ierr);
primme.matrix = &A;
primme.matrixMatvec = PETScMatvec;
/* Function that implements the matrix-vector product
A*x for solving the problem A*x = l*x */

/* Set method to solve the problem */
primme_set_method(PRIMME_DEFAULT_MIN_MATVECS, &primme);
//primme_set_method(PRIMME_DYNAMIC, &primme);
Expand Down
193 changes: 30 additions & 163 deletions examples/ex_eigs_dseq.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,21 +38,26 @@
#include <complex.h>
#include "primme.h" /* header file is required to run primme */

void LaplacianMatrixMatvec(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *ierr);
//void LaplacianMatrixMatvec(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *ierr);
void DiagonalMatrixMatvec(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *ierr);
void LaplacianApplyPreconditioner(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *ierr);
//void LaplacianApplyPreconditioner(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *ierr);
void DiagonalApplyPreconditioner(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *ierr);

#ifndef min
# define min(a, b) ((a) < (b) ? (a) : (b))
#endif

int main (int argc, char *argv[]) {

/* Solver arrays and parameters */
double *evals; /* Array with the computed eigenvalues */
double *rnorms; /* Array with the computed eigenpairs residual norms */
complex double *evecs; /* Array with the computed eigenvectors;
double *evecs; /* Array with the computed eigenvectors;
first vector starts in evecs[0],
second vector starts in evecs[primme.n],
third vector starts in evecs[primme.n*2]... */
primme_params primme;

/* PRIMME configuration struct */
double targetShifts[1];

Expand All @@ -71,17 +76,17 @@ int main (int argc, char *argv[]) {
/* Set problem parameters */
primme.n = 1000; /* set problem dimension */
primme.numEvals = 10; /* Number of wanted eigenpairs */
primme.eps = 1e-5; /* ||r|| <= eps * ||matrix|| */
primme.eps = .1; /* ||r|| <= eps * ||matrix|| */
primme.target = primme_largest;
/* Wanted the smallest eigenvalues */

/* Set preconditioner (optional) */
primme.applyPreconditioner = DiagonalApplyPreconditioner;
primme.correctionParams.precondition = 1;
//primme.applyPreconditioner = DiagonalApplyPreconditioner;
//primme.correctionParams.precondition = 1;
primme.projectionParams.projection = primme_proj_sketched;

/* Set advanced parameters if you know what are you doing (optional) */
//primme.minRestartSize = 100;
primme.minRestartSize = 100;
primme.maxBasisSize = 750;
primme.initSize = 0;
primme.locking = 0;
Expand All @@ -90,7 +95,7 @@ int main (int argc, char *argv[]) {
primme.printLevel = 4;

primme.maxBlockSize = 1;
primme.expansionParams.expansion = primme_expansion_fullLanczos;
primme.expansionParams.expansion = primme_expansion_davidson;
/* Set method to solve the problem */
//primme_set_method(PRIMME_DYNAMIC, &primme);
primme_set_method(PRIMME_DEFAULT_MIN_MATVECS, &primme);
Expand All @@ -103,11 +108,11 @@ int main (int argc, char *argv[]) {

/* Allocate space for converged Ritz values and residual norms */
evals = (double*)malloc(primme.numEvals*sizeof(double));
evecs = (complex double*)malloc(primme.n*primme.numEvals*sizeof(complex double));
evecs = (double*)malloc(primme.n*primme.numEvals*sizeof(double));
rnorms = (double*)malloc(primme.numEvals*sizeof(double));

/* Call primme */
ret = zprimme(evals, evecs, rnorms, &primme);
ret = dprimme(evals, evecs, rnorms, &primme);

if (ret != 0) {
fprintf(primme.outputFile,
Expand Down Expand Up @@ -146,149 +151,6 @@ int main (int argc, char *argv[]) {
}


/* Note that d/zprimme can be called more than once before call primme_free. */
/* Find the 5 eigenpairs closest to .5 */
// primme.numTargetShifts = 1;
// targetShifts[0] = .5;
// primme.targetShifts = targetShifts;
// primme.target = primme_closest_abs;
// primme.numEvals = 5;
// primme.initSize = 0; /* primme.initSize may be not zero after a d/zprimme;
// so set it to zero to avoid the already converged eigenvectors
// being used as initial vectors. */
//
// /* Call primme */
// ret = zprimme(evals, evecs, rnorms, &primme);
//
// if (ret != 0) {
// fprintf(primme.outputFile,
// "Error: primme returned with nonzero exit status: %d \n",ret);
// return -1;
// }
//
// /* Reporting (optional) */
// for (i=0; i < primme.initSize; i++) {
// fprintf(primme.outputFile, "Eval[%d]: %-22.15E rnorm: %-22.15E\n", i+1,
// evals[i], rnorms[i]);
// }
// fprintf(primme.outputFile, " %d eigenpairs converged\n", primme.initSize);
// fprintf(primme.outputFile, "Tolerance : %-22.15E\n",
// primme.aNorm*primme.eps);
// fprintf(primme.outputFile, "Iterations: %-" PRIMME_INT_P "\n",
// primme.stats.numOuterIterations);
// fprintf(primme.outputFile, "Restarts : %-" PRIMME_INT_P "\n", primme.stats.numRestarts);
// fprintf(primme.outputFile, "Matvecs : %-" PRIMME_INT_P "\n", primme.stats.numMatvecs);
// fprintf(primme.outputFile, "Preconds : %-" PRIMME_INT_P "\n", primme.stats.numPreconds);
// if (primme.stats.lockingIssue) {
// fprintf(primme.outputFile, "\nA locking problem has occurred.\n");
// fprintf(primme.outputFile,
// "Some eigenpairs do not have a residual norm less than the tolerance.\n");
// fprintf(primme.outputFile,
// "However, the subspace of evecs is accurate to the required tolerance.\n");
// }
//
// switch (primme.dynamicMethodSwitch) {
// case -1: fprintf(primme.outputFile,
// "Recommended method for next run: DEFAULT_MIN_MATVECS\n"); break;
// case -2: fprintf(primme.outputFile,
// "Recommended method for next run: DEFAULT_MIN_TIME\n"); break;
// case -3: fprintf(primme.outputFile,
// "Recommended method for next run: DYNAMIC (close call)\n"); break;
// }
//
//
// /* Perturb the 5 approximate eigenvectors in evecs and used them as initial solution.
// This time the solver should converge faster than the last one. */
// for (i=0; i<primme.n*5; i++)
// evecs[i] += rand()/(double)RAND_MAX*1e-4;
// primme.initSize = 5;
// primme.numEvals = 5;
//
// /* Call primme */
// ret = zprimme(evals, evecs, rnorms, &primme);
//
// if (ret != 0) {
// fprintf(primme.outputFile,
// "Error: primme returned with nonzero exit status: %d \n",ret);
// return -1;
// }
//
// /* Reporting (optional) */
// for (i=0; i < primme.initSize; i++) {
// fprintf(primme.outputFile, "Eval[%d]: %-22.15E rnorm: %-22.15E\n", i+1,
// evals[i], rnorms[i]);
// }
// fprintf(primme.outputFile, " %d eigenpairs converged\n", primme.initSize);
// fprintf(primme.outputFile, "Tolerance : %-22.15E\n",
// primme.aNorm*primme.eps);
// fprintf(primme.outputFile, "Iterations: %-" PRIMME_INT_P "\n",
// primme.stats.numOuterIterations);
// fprintf(primme.outputFile, "Restarts : %-" PRIMME_INT_P "\n", primme.stats.numRestarts);
// fprintf(primme.outputFile, "Matvecs : %-" PRIMME_INT_P "\n", primme.stats.numMatvecs);
// fprintf(primme.outputFile, "Preconds : %-" PRIMME_INT_P "\n", primme.stats.numPreconds);
// if (primme.stats.lockingIssue) {
// fprintf(primme.outputFile, "\nA locking problem has occurred.\n");
// fprintf(primme.outputFile,
// "Some eigenpairs do not have a residual norm less than the tolerance.\n");
// fprintf(primme.outputFile,
// "However, the subspace of evecs is accurate to the required tolerance.\n");
// }
//
// switch (primme.dynamicMethodSwitch) {
// case -1: fprintf(primme.outputFile,
// "Recommended method for next run: DEFAULT_MIN_MATVECS\n"); break;
// case -2: fprintf(primme.outputFile,
// "Recommended method for next run: DEFAULT_MIN_TIME\n"); break;
// case -3: fprintf(primme.outputFile,
// "Recommended method for next run: DYNAMIC (close call)\n"); break;
// }
//
//
// /* Find the next 5 eigenpairs closest to .5 */
// primme.initSize = 0;
// primme.numEvals = 5;
// primme.numOrthoConst = 5; /* solver will find solutions orthogonal to the already
// 5 approximate eigenvectors in evecs */
//
// /* Call primme */
// ret = zprimme(evals, evecs, rnorms, &primme);
//
// if (ret != 0) {
// fprintf(primme.outputFile,
// "Error: primme returned with nonzero exit status: %d \n",ret);
// return -1;
// }
//
// /* Reporting (optional) */
// for (i=0; i < primme.initSize; i++) {
// fprintf(primme.outputFile, "Eval[%d]: %-22.15E rnorm: %-22.15E\n", i+1,
// evals[i], rnorms[i]);
// }
// fprintf(primme.outputFile, " %d eigenpairs converged\n", primme.initSize);
// fprintf(primme.outputFile, "Tolerance : %-22.15E\n",
// primme.aNorm*primme.eps);
// fprintf(primme.outputFile, "Iterations: %-" PRIMME_INT_P "\n",
// primme.stats.numOuterIterations);
// fprintf(primme.outputFile, "Restarts : %-" PRIMME_INT_P "\n", primme.stats.numRestarts);
// fprintf(primme.outputFile, "Matvecs : %-" PRIMME_INT_P "\n", primme.stats.numMatvecs);
// fprintf(primme.outputFile, "Preconds : %-" PRIMME_INT_P "\n", primme.stats.numPreconds);
// if (primme.stats.lockingIssue) {
// fprintf(primme.outputFile, "\nA locking problem has occurred.\n");
// fprintf(primme.outputFile,
// "Some eigenpairs do not have a residual norm less than the tolerance.\n");
// fprintf(primme.outputFile,
// "However, the subspace of evecs is accurate to the required tolerance.\n");
// }
//
// switch (primme.dynamicMethodSwitch) {
// case -1: fprintf(primme.outputFile,
// "Recommended method for next run: DEFAULT_MIN_MATVECS\n"); break;
// case -2: fprintf(primme.outputFile,
// "Recommended method for next run: DEFAULT_MIN_TIME\n"); break;
// case -3: fprintf(primme.outputFile,
// "Recommended method for next run: DYNAMIC (close call)\n"); break;
// }
//
primme_free(&primme);
free(evals);
free(evecs);
Expand Down Expand Up @@ -332,21 +194,25 @@ void LaplacianMatrixMatvec(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, i
void DiagonalMatrixMatvec(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *ldy, int *blockSize, primme_params *primme, int *err) {

int i; /* vector index, from 0 to *blockSize-1*/
int row; /* Laplacian matrix row index, from 0 to matrix dimension */
complex double *xvec; /* pointer to i-th input vector x */
complex double *yvec; /* pointer to i-th output vector y */
int row; /* local matrix row index, from 0 to nLocal */
/* In this example, row0 is the global index of the first local row */
int row0 = primme->n / primme->numProcs * primme->procID +
min(primme->n % primme->numProcs, primme->procID);
double *xvec; /* pointer to i-th input vector x */
double *yvec; /* pointer to i-th output vector y */

for (i = 0; i < *blockSize; i++) {
xvec = (complex double *)x + *ldx*i;
yvec = (complex double *)y + *ldy*i;
for (row = 0; row < primme->n; row++) {
yvec[row] = 0.0;
yvec[row] += (complex double)((row+1)*(row+1))*xvec[row];
}
for (i=0; i<*blockSize; i++) {
xvec = (double *)x + *ldx*i;
yvec = (double *)y + *ldy*i;
for (row = 0; row < primme->nLocal; row++) {
double this_row = (row + row0 + 1);
yvec[row] = this_row * this_row * xvec[row];
}
}
*err = 0;
}


/* This performs Y = M^{-1} * X, where
- X, input dense matrix of size primme.n x blockSize;
Expand Down Expand Up @@ -387,3 +253,4 @@ void DiagonalApplyPreconditioner(void *x, PRIMME_INT *ldx, void *y, PRIMME_INT *
}
*ierr = 0;
}

0 comments on commit 11a2166

Please sign in to comment.