[Bf-blender-cvs] SVN commit: /data/svn/bf-blender [38767] branches/soc-2011-avocado/blender/ intern/autoseam: Crash fix on the modularization of eigen value calculation code.

shuvro sarker shuvro05 at gmail.com
Wed Jul 27 21:37:34 CEST 2011


Revision: 38767
          http://projects.blender.org/scm/viewvc.php?view=rev&root=bf-blender&revision=38767
Author:   shuvro
Date:     2011-07-27 19:37:33 +0000 (Wed, 27 Jul 2011)
Log Message:
-----------
Crash fix on the modularization of eigen value calculation code. Eigen calculation related codes are moved from AutoseamAdjacency

Modified Paths:
--------------
    branches/soc-2011-avocado/blender/intern/autoseam/AutoseamAdjacency.cpp
    branches/soc-2011-avocado/blender/intern/autoseam/AutoseamAdjacency.h
    branches/soc-2011-avocado/blender/intern/autoseam/EigenSolver.cpp
    branches/soc-2011-avocado/blender/intern/autoseam/EigenSolver.h
    branches/soc-2011-avocado/blender/intern/autoseam/EigenSolverArpack.cpp
    branches/soc-2011-avocado/blender/intern/autoseam/EigenSolverArpack.h

Modified: branches/soc-2011-avocado/blender/intern/autoseam/AutoseamAdjacency.cpp
===================================================================
--- branches/soc-2011-avocado/blender/intern/autoseam/AutoseamAdjacency.cpp	2011-07-27 19:08:18 UTC (rev 38766)
+++ branches/soc-2011-avocado/blender/intern/autoseam/AutoseamAdjacency.cpp	2011-07-27 19:37:33 UTC (rev 38767)
@@ -45,16 +45,9 @@
 
 #define THRESHOLD_ZERO 0.0005
 #define EXTRA_MEMORY_SIZE 10
-//#define NUM_EIGEN_VAL  25
 
-#define minimum(a,b)  (a <= b) ? a : b;
 
 
-extern "C"{
-    extern  int dsaupd_(integer *, char *, integer *, char *, integer *, doublereal *, doublereal *, integer *, doublereal *, integer *, integer *, integer *, doublereal *, doublereal *, integer *, integer *, ftnlen, ftnlen);
-    extern int dseupd_(logical *, char *, logical *, doublereal *, doublereal *, integer *, doublereal *, char *, integer *, char *, integer *, doublereal *, doublereal *, integer *, doublereal *, integer *, integer *, integer *, doublereal *, doublereal *, integer *, integer *, ftnlen, ftnlen, ftnlen);
-}
-
 AutoseamAdjacency::AutoseamAdjacency(int dimension)
 {
 	m_adjacency = MatrixXd::Zero(dimension, dimension);
@@ -64,30 +57,6 @@
 	//threshold_value = min_value;
 }
 
-int AutoseamAdjacency:: mul_matrix_vector(int dimension, double *vec, double * result)
-{
-    int i;
-    Eigen::VectorXd output(dimension);
-    
-    
-    //printf("the dimesion passed here : %d\n The output values are: \n", dimension);
-    
-    for(i = 0; i < dimension; i++){
-        //printf("%lf ", *(vec + i));
-        mul_vec(i) = *(vec + i);
-    }
-    
-    output = m_adjacency * mul_vec;
-    
-    //pack the values into result
-    for(i = 0; i < dimension; i++){
-        //printf("%lf ", output(i));
-        *(result + i) = output(i);
-    }
-    
-    return 0;
-
-}
 void AutoseamAdjacency::set(int row, int col, float value)
 {
 	m_adjacency(row, col) = value;
@@ -151,245 +120,7 @@
 	}
 }
 
-//int AutoseamAdjacency::is_adjacent(int index1, int index2)
-//{
-//    //if the indexes are not adjacent it returns -1
-//    return m_adjacency(index1, index1) ? 1 : -1;
-//}
 
-
-/**
- * The function will calculate eigen values using the functions
- * of arpack.
- */
-double * AutoseamAdjacency:: calculate_eigen_arpack()
-{
-    #ifdef WITH_ARPACK
-    
-        
-        integer  n;
-        integer  ido, ncv, nev;
-        doublereal tol;
-        char bmat[1];
-        integer info;
-        logical rvec;
-        integer ierr, mode1;
-        doublereal sigma;
-        char which[3] = "SM"; // we are looking for the smallest eigen values
-    
-    
-        //integer nconv;
-        integer ipntr[11];
-        integer iparam[11];
-        integer ldv;
-    
-        integer ishfts, maxitr, lworkl;
-        int vectors_size;
-    
-        // size : number of eigen vectors
-        logical *select = NULL;
-        // LWORKL must be at least NCV**2 + 8*NCV.
-        doublereal *workl = NULL; 
-        // 3*N
-        doublereal *workd = NULL;
-        doublereal *d__ = NULL;  
-        doublereal *v = NULL;
-        doublereal *resid = NULL;
-
-
-    
-        
-    
-        // This is the dimension of the matrix
-        n = matrix_dimension;
-                    
-        /**
-         *The number of eigen values we want to calculate
-         * This should be the minimum number between matrix_dimension and NUM_EIGEN_VAL
-         * The constraints - nev <= maxnev
-         * nev + 1 <= ncv <= maxncv
-         * n >= ncv
-         */
-        
-        
-        ldv = (integer)matrix_dimension;
-        nev = minimum(matrix_dimension - 1, NUM_EIGEN_VAL);
-        num_eigen_values = nev;
-        ncv = minimum(nev + 1, n);
-    
-    
-        *(unsigned char *)bmat = 'I';
-    
-        
-        if (n > ldv) {
-            printf(" ERROR with _SSIMP: N is greater than MAXN \n");
-            //goto L9000;
-            return NULL;
-        }
-        else if (nev > NUM_EIGEN_VAL) {
-            printf(" ERROR with _SSIMP: NEV is greater than MAXNEV \n");
-            //goto L9000;
-            return NULL;
-                                                                            
-        }
-        else if (ncv <= nev || ncv > n) {
-            printf(" ERROR with _SSIMP: NCV is greater than MAXNCV \n");
-            // goto L9000;
-            return NULL;
-        }
-    
-            
-        //vectors_size = matrix_dimension * nev;
-        vectors_size = ldv * ncv;
-        v            = (doublereal*)MEM_callocN(vectors_size * sizeof(doublereal), "eigen_vectors");
-        workl        = (doublereal*)MEM_callocN((ncv*ncv + 8*ncv) * sizeof(doublereal), "workl");
-        //d__          = (doublereal*)MEM_callocN(2*matrix_dimension * sizeof(doublereal), "eigen_values");
-        d__          = (doublereal*)MEM_callocN(2*ncv * sizeof(doublereal), "eigen_values");
-    
-        resid        = (doublereal*)MEM_callocN(matrix_dimension * sizeof(doublereal), "resid");
-        workd        = (doublereal*)MEM_callocN(3 * matrix_dimension * sizeof(doublereal), "workd");
-        
-        select       = (logical*)MEM_callocN(ncv * sizeof(logical), "select");
-        
-        
-
-
-        lworkl = ncv * (ncv + 8);
-        tol = 0.;
-        info = 0;
-        ido = 0;
-    
-        ishfts = 1;
-        maxitr = 350;
-        mode1 = 1;
-    
-        iparam[0] = ishfts;
-        iparam[2] = maxitr;
-        iparam[6] = mode1;
-    
-        // main loop        
-        //L10:
-        while (true) {
-        
-    
-    
-            /* %---------------------------------------------% */
-            /* | Repeatedly call the routine DSAUPD and take | */
-            /* | actions indicated by parameter IDO until    | */
-            /* | either convergence is indicated or maxitr   | */
-            /* | has been exceeded.                          | */
-            /* %---------------------------------------------% */
-    
-            //dsaupd_(&ido, bmat, &n, which, &nev, &tol, resid, &ncv, v, &c__256, iparam, ipntr, workd, workl, &lworkl, &info, (ftnlen)1, (ftnlen)2);
-            dsaupd_(&ido, bmat, &n, which, &nev, &tol, resid, &ncv, v, &ldv, iparam, ipntr, workd, workl, &lworkl, &info, (ftnlen)1, (ftnlen)2);
-
-    
-            if (ido == -1 || ido == 1) {
-                //perform matrix-vector multiplication
-                //av_(&nx, &workd[ipntr[0] - 1], &workd[ipntr[1] - 1]);
-                mul_matrix_vector(n, &workd[ipntr[0] - 1], &workd[ipntr[1] - 1]);
-                //goto L10;
-            }
-            else{
-                break;
-            }
-            
-        }
-    
-    
-        if (info < 0) {
-        
-            // error occured
-            printf(" Error with _saupd, info = %d\n", info);
-            printf(" Check documentation in _saupd \n");
-        } 
-        else {
-        
-            // no major error
-            rvec = TRUE_;
-        
-            dseupd_(&rvec, "All", select, d__, v, &ldv, &sigma, bmat, &n, 
-                which, &nev, &tol, resid, &ncv, v, &ldv, iparam, ipntr, 
-                workd, workl, &lworkl, &ierr, (ftnlen)3, (ftnlen)1, (ftnlen)2);
-        
-        
-            /*         %----------------------------------------------% */
-            /*         | Eigenvalues are returned in the first column | */
-            /*         | of the two dimensional array D and the       | */
-            /*         | corresponding eigenvectors are returned in   | */
-            /*         | the first NCONV (=IPARAM(5)) columns of the  | */
-            /*         | two dimensional array V if requested.        | */
-            /*         | Otherwise, an orthogonal basis for the       | */
-            /*         | invariant subspace corresponding to the      | */
-            /*         | eigenvalues in D is returned in V.           | */
-            /*         %----------------------------------------------% */
-        
-            if (ierr != 0) {
-            
-               printf(" Error with _seupd, info = %d\n", ierr);
-               printf(" Check the documentation of _seupd. \n");
-                
-            } 
-            else {
-                
-                // Block to print eigen values and eigen vectors
-                
-                #ifdef EIGEN_DEBUG_ARPACK
-                {
-                    int i,j;
-                    
-                    for(i = 0; i < nev; i++)
-                    {
-                        printf("eigen value : %f\n", d__[i]);
-                        printf("The eigen vector %d :\n", i);
-                        
-                        //prints one eigen vector
-                        for(j = 0; j < matrix_dimension; j++){
-                            printf("%lf ", v[i*ldv + j]);
-                        }
-                    }
-                }
-                #endif
-                //set the eigen values
-                int i;
-                eigen_values = Eigen::VectorXd(nev);
-                for(i = 0; i < nev; i++) eigen_values(i) = d__[i];
-                
-                //memory deallocation
-                
-                MEM_freeN(d__);
-                MEM_freeN(resid);
-                MEM_freeN(workd);
-                MEM_freeN(workl);
-                MEM_freeN(select);
-
-
-                
-                // return the eigen vectors
-                return v;
-                                
-            
-            }
-        
-        
-        }
-    
-    
-    MEM_freeN(d__);
-    MEM_freeN(resid);
-    MEM_freeN(workd);
-    MEM_freeN(workl);
-    MEM_freeN(select);
-    
-    *v = NULL;
-    
-    return v;
-    
-    
-    #endif
-    
-}
-
 bool AutoseamAdjacency::generate_seam()
 {
     
@@ -407,14 +138,14 @@
 		a(i,i) = -rowsum;
 	}
     
-    eigen_vectors = calculate_eigen_arpack();
+    //eigen_vectors = calculate_eigen_arpack();
     
 
     // Need to modify the errors of the commented code.
     
-    /*EigenSolverArpack solver(matrix_dimension, num_eigen_values);
+    EigenSolverArpack solver(matrix_dimension);
     solver.matrix =  a;
-    eigen_vectors = solver.calculate_eigen_space();*/
+    eigen_vectors = solver.calculate_eigen_space();
 	
     
     //dssimp_();
@@ -429,12 +160,16 @@
 		MatrixXd aplus;
 		MatrixXd aminus;
 	
-		while ( (f < num_eigen_values)) {
+		//while ( (f < num_eigen_values)) {
+        while ( (f < solver.num_eigen_vectors)) {

@@ Diff output truncated at 10240 characters. @@



More information about the Bf-blender-cvs mailing list