Improved ARPACK wrapper doc

src/shogun/mathematics/arpack.cpp

@@ -85,18 +85,22 @@ void arpack_dsaupd(double* matrix, int n, int nev, const char* which,
 	// shift-invert mode init
 	if (mode==3)
 	{
+		// subtract shift from main diagonal
 		if (shift!=0.0)
 		{
 			for (int i=0; i<n; i++)
 				matrix[i*n+i] -= shift;
 		}
 
+		// compute factorization according to pos value
 		if (pos)
 		{
+			SG_SDEBUG("ARPACK: Using Cholesky factorization");
 			clapack_dpotrf(CblasColMajor,CblasUpper,n,matrix,n);
 		}
 		else
 		{
+			SG_SDEBUG("ARPACK: Using LUP factorization");
 			ipiv = SG_MALLOC(int, n);
 			clapack_dgetrf(CblasColMajor,n,n,matrix,n,ipiv);
 		}
@@ -114,6 +118,7 @@ void arpack_dsaupd(double* matrix, int n, int nev, const char* which,
 		if ((ido==1)||(ido==-1))
 		{
 			if (mode==1)
+			// compute (workd+ipntr[1]-1) = A*(workd+ipntr[0]-1)
 			{
 				cblas_dsymv(CblasColMajor,CblasUpper,
 				            n,1.0,matrix,n,
@@ -122,6 +127,7 @@ void arpack_dsaupd(double* matrix, int n, int nev, const char* which,
 			}
 			if (mode==3)
 			{
+			// solve system of eqs A*(workd+ipntr[0]-1) = (workd+ipntr[1]-1)
 				for (i=0; i<n; i++)
 					tmp[i] = (workd+ipntr[0]-1)[i];
 
@@ -168,17 +174,20 @@ void arpack_dsaupd(double* matrix, int n, int nev, const char* which,
 		// specify that eigenvectors to be computed too		
 		int rvec = 1;
 
+		// call dseupd_ routine
 		dseupd_(&rvec, all_, select, d, v, &ldv, &sigma, bmat,
 		        &n, which_, &nev, &tol, resid, &ncv, v, &ldv,
 		        iparam, ipntr, workd, workl, &lworkl, &ierr);
 
+		// check for errors
 		if (ierr!=0)
 		{
 			SG_SWARNING("DSEUPD failed with status=%d", ierr);
 			status = -1;
 		}
 		else
 		{
+			// store eigenpairs to specified arrays
 			for (i=0; i<nev; i++)
 			{	
 				eigenvalues[i] = d[i];

src/shogun/mathematics/arpack.h

@@ -44,25 +44,28 @@ namespace shogun
  * e.g. ncv and workaround arrays sizes.
  *
  * @param matrix symmetric real matrix of size n*n
+ *        (contains LU or Cholesky factors on exit if mode==3)
  * @param n size of matrix
  * @param nev number of eigenpairs to compute (nev<=n)
  * @param which eigenvalue finding strategy. Possible values:
  *        - "LM": nev Largest Magnitude eigenvalues
  *        - "SM": nev Smallest Magnitude eigenvalues
  *        - "LA": nev Largest Algebraic eigenvalues
+ *                (if mode==3 eigenvalues from the right of shift) 
  *        - "SA": nev Smallest Algebraic eigenvalues
+ *                (if mode==3 eigenvalues from the left of shift)
  *        - "BE": half of nev from each end of the spectrum, i.e. nev%2
  *                smallest and nev%2 largest eigenvalues. If nev is odd,
  *                one more largest eigenvalue will be computed
  * @param mode shift-mode of IRLM. Possible values:
  *        - 1: regular mode
- *        - 2: shift-invert mode
- * @param pos true if positive definite (when Cholesky is used)
- * @param shift shift for shift-invert (2) mode of IRLM. In this mode
+ *        - 3: shift-invert mode
+ * @param pos true if matrix is positive definite (Cholesky factorization)
+ * @param shift shift for shift-invert (3) mode of IRLM. In this mode
  *        routine will compute eigenvalues near provided shift
  * @param eigenvalues array of size nev to hold computed eigenvalues
  * @param eigenvectors array of size nev*n to hold computed eigenvectors
- * @param status on output sets -1 if failed
+ * @param status on output -1 if computation failed
  */
 void arpack_dsaupd(double* matrix, int n, int nev, const char* which, int mode, bool pos,
                    double shift, double* eigenvalues, double* eigenvectors, int& status);