Solver parameters for eigenvalue routines
Solver parameters for eigenvalue routines
I'm writing C code to solve a sparse, hermitian eigenvalue problem. I'm working from the provided 'example_sparse.c' which shows how to use routines for solving Ax=b. At the part where a solver and preconditioner are chosen the author comments, alluding to some mysterious 'documentation' which I haven't found yet. Where can I find this documentation or can I get a quick rundown of the solver parameters that the zlobpcg routine needs?
Re: Solver parameters for eigenvalue routines
Documentation is under magma/docs/html/index.html
It is also available from the Documentation link on MAGMA's website,
http://icl.cs.utk.edu/magma/
http://icl.cs.utk.edu/projectsfiles/magma/doxygen/
For sparse, there is an Overview page (the second one  I'll fix that duplicate name).
For lobpcg, look under Routines > Sparse > Sparse eigenvalues > Symmetric/Hermitian eigenvalues > * precision > magma_*lobpcg.
for whatever precision you need.
Feel free to ask again if you can't find what you need.
mark
It is also available from the Documentation link on MAGMA's website,
http://icl.cs.utk.edu/magma/
http://icl.cs.utk.edu/projectsfiles/magma/doxygen/
For sparse, there is an Overview page (the second one  I'll fix that duplicate name).
For lobpcg, look under Routines > Sparse > Sparse eigenvalues > Symmetric/Hermitian eigenvalues > * precision > magma_*lobpcg.
for whatever precision you need.
Feel free to ask again if you can't find what you need.
mark
Re: Solver parameters for eigenvalue routines
I'm looking for the specifics of the solver parameters, the corresponding section in example_sparse.c is
But since this is for solving Ax=b I'm not sure these parameters are applicable to zlobpcg. Most of these parameters are pretty self explanatory but some are less clear and their allowed values I can only guess at.
The zlobpcg page (http://icl.cs.utk.edu/projectsfiles/mag ... a8c3ba9f51) only mentions that the parameters are required but not what they do.
Code: Select all
// Choose a solver, preconditioner, etc.  see documentation for options.
opts.solver_par.solver = Magma_PIDRMERGE;
opts.solver_par.restart = 8;
opts.solver_par.maxiter = 1000;
opts.solver_par.rtol = 1e10;
opts.solver_par.maxiter = 1000;
opts.precond_par.solver = Magma_ILU;
opts.precond_par.levels = 0;
opts.precond_par.trisolver = Magma_CUSOLVE;
The zlobpcg page (http://icl.cs.utk.edu/projectsfiles/mag ... a8c3ba9f51) only mentions that the parameters are required but not what they do.

 Posts: 90
 Joined: Tue Sep 02, 2014 5:44 pm
Re: Solver parameters for eigenvalue routines
Dear Olafur,
LOBPCG is an iterative eigenvalue solver that can be used to compute a certain number of the largest/smallest eigenvalues. It is a Krylovbased method, and you can add a preconditioner to it. Also the stopping criterion and the total number of iterations are valid. In other words, a setup for LOBPCG could look like:
opts.solver_par.solver = Magma_LOBPCG; // choose an LOBPCG solver
opts.solver_par.num_eigenvalues = 4; // number of eigenvalues you want to compute
opts.solver_par.maxiter = 1000; // max number of iterations
opts.solver_par.rtol = 1e10; // stopping criterion  relative accuracy of first eigenvalue
opts.precond_par.solver = Magma_ILU; // preconditioner
opts.precond_par.levels = 0; // ILU(0)  no fillin
opts.precond_par.trisolver = Magma_CUSOLVE; //exact triangular solves
Maybe this helps?
Hartwig
LOBPCG is an iterative eigenvalue solver that can be used to compute a certain number of the largest/smallest eigenvalues. It is a Krylovbased method, and you can add a preconditioner to it. Also the stopping criterion and the total number of iterations are valid. In other words, a setup for LOBPCG could look like:
opts.solver_par.solver = Magma_LOBPCG; // choose an LOBPCG solver
opts.solver_par.num_eigenvalues = 4; // number of eigenvalues you want to compute
opts.solver_par.maxiter = 1000; // max number of iterations
opts.solver_par.rtol = 1e10; // stopping criterion  relative accuracy of first eigenvalue
opts.precond_par.solver = Magma_ILU; // preconditioner
opts.precond_par.levels = 0; // ILU(0)  no fillin
opts.precond_par.trisolver = Magma_CUSOLVE; //exact triangular solves
Maybe this helps?
Hartwig
Re: Solver parameters for eigenvalue routines
Thanks for the help so far, I'm still having some trouble and need some clarification.
The code I have so far (which I've pasted at the bottom of this reply), which is based on examples provided with the installation, has a seg fault which I'd like some help resolving. I've located the error at the point where I assign values to the 'val' vector at line 43 (the nnz elements of my matrx). I've tried several things but maybe there's something obvious I'm doing wrong (I'm not used to working in C).
Also, I have a couple of more general questions:
1. Is calling the routine magma_z_solver equivalent to calling the routine corresponding to the solver parameters?
2. zlobpcg doesn't seem to have output parameters in its documentation, does it return the eigenvalues as opposed to its error code?
3. Does zlobpcg not calculate eigenvectors at all?
4. I'm used to seeing routines that handle Hermitian or symmetric matrices take a parameter that specifies whether the input matrix is lower or upper triangular, to save memory by only having to create one half of the matrix. Does zlobpcg not do this or does it just do one or the other as a default?
Compiled with the makefile:
Sorry if I'm being a bit demanding but thanks for the help so far.
The code I have so far (which I've pasted at the bottom of this reply), which is based on examples provided with the installation, has a seg fault which I'd like some help resolving. I've located the error at the point where I assign values to the 'val' vector at line 43 (the nnz elements of my matrx). I've tried several things but maybe there's something obvious I'm doing wrong (I'm not used to working in C).
Also, I have a couple of more general questions:
1. Is calling the routine magma_z_solver equivalent to calling the routine corresponding to the solver parameters?
2. zlobpcg doesn't seem to have output parameters in its documentation, does it return the eigenvalues as opposed to its error code?
3. Does zlobpcg not calculate eigenvectors at all?
4. I'm used to seeing routines that handle Hermitian or symmetric matrices take a parameter that specifies whether the input matrix is lower or upper triangular, to save memory by only having to create one half of the matrix. Does zlobpcg not do this or does it just do one or the other as a default?
Code: Select all
//Compiles with the corresponding makefile.
//Takes in the 3 vectors of CSR format, creates MAGMA format sparse matrix, returns its spectrum.
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#include "magma_v2.h"
#include "magmasparse.h"
int main()
{
//Setup
int n = 3; //Dimension of the nxn matrix
magmaDoubleComplex *sol;
//Testing with matrix {{4,1i,2+i},{1+i,1,1},{2i,1,3}} that has spectrum {6.08807,1.61431,1.52624}
int *col; //Column indices of NNZ in A
col = (int*) calloc(n*n, sizeof(int));
col[0] = 0; col[1]=1; col[2] = 2;
col[3] = 0; col[4]=1; col[5] = 2;
col[6] = 0; col[7]=1; col[8] = 2;
int *row; //Row pointer of NNZ in A
row = (int*) calloc(n+1, sizeof(int));
row[0] = 0; row[1] = 3; row[2] = 6; row[3] = 9;
//NNZ of the matrix A:
magmaDoubleComplex *val;
magma_zmalloc(&val, 9);
if (val == NULL)
{
fprintf( stderr, "malloc failed\n" );
return 0;
}
else
{
fprintf( stderr, "malloc succeeded\n" );
}
val[0] = MAGMA_Z_MAKE(4,0);
// val[1] = MAGMA_Z_MAKE(1,1); val[2] = MAGMA_Z_MAKE(2,1);
// val[3] = MAGMA_Z_MAKE(1,1); val[4] = MAGMA_Z_MAKE(1,0); val[5] = MAGMA_Z_MAKE(1,0);
// val[6] = MAGMA_Z_MAKE(2,1); val[7] = MAGMA_Z_MAKE(1,0); val[8] = MAGMA_Z_MAKE(3,0);
printf("val created successfully\n");
//Initialize MAGMA
magma_init();
magma_zopts opts;
magma_queue_t queue;
magma_queue_create(0, &queue);
magma_z_matrix A={Magma_CSR}, dA={Magma_CSR};
magma_z_matrix b={Magma_CSR}, db={Magma_CSR};
magma_z_matrix x={Magma_CSR}, dx={Magma_CSR};
//Pass the system to MAGMA
magma_zcsrset(n,n,row,col,val,&A,queue);
// Choose a solver, preconditioner, etc.  see documentation for options.
opts.solver_par.solver = Magma_LOBPCG; // choose an LOBPCG solver
opts.solver_par.num_eigenvalues = 3; // number of eigenvalues you want to compute
opts.solver_par.maxiter = 1000; // max number of iterations
opts.solver_par.rtol = 1e10; // stopping criterion  relative accuracy of first eigenvalue
opts.precond_par.solver = Magma_ILU; // preconditioner
opts.precond_par.levels = 0; // ILU(0)  no fillin
opts.precond_par.trisolver = Magma_CUSOLVE; //exact triangular solves
//Copy the system to device
magma_zmtransfer(A, &dA, Magma_CPU, Magma_DEV, queue);
//Generate the preconditioner
magma_z_precondsetup(dA, db, &opts.solver_par, &opts.precond_par, queue);
//Calling zlobpcg to find eigenvalues
magma_z_solver(dA, db, &dx, &opts, queue);
// magma_zlobpcg(A, &opts, &opts, queue);
//Then copy the solution back to the host
magma_zmfree( &x, queue );
magma_zmtransfer( dx, &x, Magma_CPU, Magma_DEV, queue );
//and back to the application code
magma_zvget( x, &n, &n, &sol, queue );
//Free the allocated memory...
magma_zmfree( &dx, queue );
magma_zmfree( &db, queue );
magma_zmfree( &dA, queue );
//and finalize MAGMA.
magma_queue_destroy( queue );
magma_finalize();
//Output
int i;
for (i = 0; i < n; ++i) {
printf("%.4f\n", sol[i]);
}
printf("\n\n\n\n");
return 0;
}
Code: Select all
default:
gcc Wall DADD_ I/usr/local/magma/include I/usr/local/magma/sparse/include I/usr/local/cuda/include c o Magma.o Magma.c
gcc Wall o Magma Magma.o L/usr/local/magma/lib lmagma_sparse lmagma L/usr/local/cuda/lib64 lcublas lcudart lcusparse L/usr/local/openblas/lib lopenblas

 Posts: 90
 Joined: Tue Sep 02, 2014 5:44 pm
Re: Solver parameters for eigenvalue routines
Olafur,
1. Is calling the routine magma_z_solver equivalent to calling the routine corresponding to the solver parameters?
> Yes, the magma_z_solver is a wrapper for calling the specific routine. You are free to use the routine directly, if you prefer.
2. zlobpcg doesn't seem to have output parameters in its documentation, does it return the eigenvalues as opposed to its error code?
> zlobpcg returns eigenvalues and eigenvectors. They are stored in magma_z_solver_par>eigenvalues and magma_z_solver_par>eigenvectors. The return value, indeed, is the error code only.
3. Does zlobpcg not calculate eigenvectors at all?
> see above. A good start might be to use testing_dsolver solver LOBPCG testmatrix.mtx to see how the library logic works.
4. I'm used to seeing routines that handle Hermitian or symmetric matrices take a parameter that specifies whether the input matrix is lower or upper triangular, to save memory by only having to create one half of the matrix. Does zlobpcg not do this or does it just do one or the other as a default?
> There is no feature exploiting symmetry. In fact, I am not aware of any GPU library doing this.
Thanks, Hartwig
1. Is calling the routine magma_z_solver equivalent to calling the routine corresponding to the solver parameters?
> Yes, the magma_z_solver is a wrapper for calling the specific routine. You are free to use the routine directly, if you prefer.
2. zlobpcg doesn't seem to have output parameters in its documentation, does it return the eigenvalues as opposed to its error code?
> zlobpcg returns eigenvalues and eigenvectors. They are stored in magma_z_solver_par>eigenvalues and magma_z_solver_par>eigenvectors. The return value, indeed, is the error code only.
3. Does zlobpcg not calculate eigenvectors at all?
> see above. A good start might be to use testing_dsolver solver LOBPCG testmatrix.mtx to see how the library logic works.
4. I'm used to seeing routines that handle Hermitian or symmetric matrices take a parameter that specifies whether the input matrix is lower or upper triangular, to save memory by only having to create one half of the matrix. Does zlobpcg not do this or does it just do one or the other as a default?
> There is no feature exploiting symmetry. In fact, I am not aware of any GPU library doing this.
Thanks, Hartwig
Re: Solver parameters for eigenvalue routines
I've tried calling zlobpcg directly, but I get a 'cusparse: not initialized' error. I've tried running the conjugategradient example provided in CUDA for testing cuSparse and that works, and simply including cusparse.h does not work
Code:
Makefile:
Code:
Code: Select all
//Compiles with the corresponding makefile.
//Takes in the 3 vectors of CSR format, creates MAGMA format sparse matrix, returns its spectrum.
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#include "magma_v2.h"
#include "magmasparse.h"
int main()
{
//Setup
int n = 3; //Dimension of the nxn matrix
magmaDoubleComplex *sol;
//Testing with matrix {{4,1i,2+i},{1+i,1,1},{2i,1,3}} that has spectrum {6.08807,1.61431,1.52624}
int *col; //Column indices of NNZ in A
col = (int*) calloc(n*n, sizeof(int));
col[0] = 0; col[1]=1; col[2] = 2;
col[3] = 0; col[4]=1; col[5] = 2;
col[6] = 0; col[7]=1; col[8] = 2;
int *row; //Row pointer of NNZ in A
row = (int*) calloc(n+1, sizeof(int));
row[0] = 0; row[1] = 3; row[2] = 6; row[3] = 9;
//NNZ of the matrix A:
magmaDoubleComplex *val;
val = malloc(n*n*sizeof(magmaDoubleComplex)); //Allocates the struct on the stack
if (val == NULL)
{
fprintf( stderr, "malloc failed\n" );
return 0;
}
else
{
fprintf( stderr, "malloc succeeded\n" );
}
val[0] = MAGMA_Z_MAKE(12,0); val[1] = MAGMA_Z_MAKE(1,1); val[2] = MAGMA_Z_MAKE(2,1);
val[3] = MAGMA_Z_MAKE(1,1); val[4] = MAGMA_Z_MAKE(1,0); val[5] = MAGMA_Z_MAKE(1,0);
val[6] = MAGMA_Z_MAKE(2,1); val[7] = MAGMA_Z_MAKE(1,0); val[8] = MAGMA_Z_MAKE(3,0);
printf("val created successfully\n");
//Initialize MAGMA
magma_init();
magma_zopts opts;
magma_queue_t queue;
magma_queue_create(0, &queue);
magma_z_matrix A={Magma_CSR}, dA={Magma_CSR};
magma_z_matrix b={Magma_CSR}, db={Magma_CSR};
magma_z_matrix x={Magma_CSR}, dx={Magma_CSR};
//Pass the system to MAGMA
magma_zcsrset(n,n,row,col,val,&A,queue);
printf("system successfully passed to MAGMA\n");
// Choose a solver, preconditioner, etc.  see documentation for options.
opts.solver_par.solver = Magma_LOBPCG; // choose an LOBPCG solver
opts.solver_par.num_eigenvalues = 3; // number of eigenvalues you want to compute
opts.solver_par.maxiter = 1000; // max number of iterations
opts.solver_par.rtol = 1e10; // stopping criterion  relative accuracy of first eigenvalue
opts.precond_par.solver = Magma_ILU; // preconditioner
opts.precond_par.levels = 0; // ILU(0)  no fillin
opts.precond_par.trisolver = Magma_CUSOLVE; //exact triangular solves
// Initialize the solver.
magma_zsolverinfo_init( &opts.solver_par, &opts.precond_par, queue );
//Copy the system to device
magma_zmtransfer(A, &dA, Magma_CPU, Magma_DEV, queue);
printf("system successfully copied to device\n");
//Generate the preconditioner
magma_z_precondsetup(dA, db, &opts.solver_par, &opts.precond_par, queue);
//Calling zlobpcg to find eigenvalues
magma_int_t err;
// magma_z_solver(dA, db, &dx, &opts, queue);
err = magma_zlobpcg(dA, &opts.solver_par, &opts.precond_par, queue);
printf("%s\n", magma_strerror(err));
printf("zlobpcg successfully run\n");
//Then copy the solution back to the host
printf("%f \n", (void*) opts.solver_par.eigenvalues);
magma_zmfree( &x, queue );
magma_zmtransfer( dx, &x, Magma_CPU, Magma_DEV, queue );
printf("system successfully copied back to host\n");
//and back to the application code
magma_zvget( x, &n, &n, &sol, queue );
//Free the allocated memory...
magma_zmfree( &dx, queue );
magma_zmfree( &db, queue );
magma_zmfree( &dA, queue );
//and finalize MAGMA.
magma_queue_destroy( queue );
magma_finalize();
//Output
int i;
for (i = 0; i < n; ++i) {
printf("%.4f\n", sol[i]);
}
printf("\n\n\n\n");
return 0;
}
Code: Select all
default:
gcc g Wall DADD_ I/usr/local/magma/include I/usr/local/magma/sparse/include I/usr/local/cuda/include c o Magma.o Magma.c
gcc g Wall o Magma Magma.o L/usr/local/magma/lib lmagma_sparse lmagma L/usr/local/cuda/lib64 lcublas lcudart lcusparse L/usr/local/openblas/lib lopenblas

 Posts: 90
 Joined: Tue Sep 02, 2014 5:44 pm
Re: Solver parameters for eigenvalue routines
We are looking into this, don't have a solution yet.
Please allow us a few days. Thanks, Hartwig
Please allow us a few days. Thanks, Hartwig