ports/147487: [PATCH] strange behaviour of math/arpack linking with lapack library

Eijiro Shibusawa ej-sib at ice.uec.ac.jp
Sat Jun 5 03:40:05 UTC 2010


>Number:         147487
>Category:       ports
>Synopsis:       [PATCH] strange behaviour of math/arpack linking with lapack library
>Confidential:   no
>Severity:       non-critical
>Priority:       low
>Responsible:    freebsd-ports-bugs
>State:          open
>Quarter:        
>Keywords:       
>Date-Required:
>Class:          sw-bug
>Submitter-Id:   current-users
>Arrival-Date:   Sat Jun 05 03:40:05 UTC 2010
>Closed-Date:
>Last-Modified:
>Originator:     Eijiro Shibusawa
>Release:        FreeBSD 8.0-RELEASE-p3 i386
>Organization:
N/A
>Environment:
FreeBSD myhost 8.0-RELEASE-p3 FreeBSD 8.0-RELEASE-p3 #0: Wed May 26 05:45:12 UTC 2010     root at i386-builder.daemonology.net:/usr/obj/usr/src/sys/GENERIC  i386
>Description:
When some code is linked with both arpack and lapack library, the arpack routine dsaupd exits with error code -9999, which indicates the routine could not build an Arnoldi factorization.
However the value of iparam[4] indicates all the requested  eigenvalues are successfully converged.
This behaviour seemingly caused by symbol duplication problem described in following website.

http://mathema.tician.de/node/373

After applying the patch provided at the site, I found the code successfully finds all the requested eigen values.
>How-To-Repeat:
# cd /usr/ports/math/arpack 

# make install WITH_ATLAS=YES

$ gcc44 -v
Using built-in specs.
Target: i386-portbld-freebsd8.0
Configured with: ./../gcc-4.4-20090915/configure --disable-nls --with-system-zlib --with-libiconv-prefix=/usr/local --with-gmp=/usr/local --program-suffix=44 --libdir=/usr/local/lib/gcc44 --libexecdir=/usr/local/libexec/gcc44 --with-gxx-include-dir=/usr/local/lib/gcc44/include/c++/ --disable-rpath --prefix=/usr/local --mandir=/usr/local/man --infodir=/usr/local/info/gcc44 --build=i386-portbld-freebsd8.0
Thread model: posix
gcc version 4.4.2 20090915 (prerelease) (GCC) 

$ gcc44 arpack.c -oa.out -Wall -I/usr/local/include -L/usr/local/lib -L. -llapack -lcblas -lf77blas -latlas -lm -lgfortran -larpack

$ ./a.out
Error with dsaupd, info = -9999 
Check documentation in dsaupd
3/3 eigen value(s) are successfully converged.

================================ arpack.c =======================================
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#include <cblas.h>

int dsaupd_(long *ido, char *bmat, long *n, char *which,
	    long *nev, double *tol, double *resid, long *ncv,
	    double *v, long *ldv, long *iparam, long *ipntr,
	    double *workd, double *workl, long *lworkl,
	    long *info);

int dseupd_(long *rvec, char *howmny, long *select,
	    double *d, double *z, long *ldz,
	    double *sigma, char *bmat, long *n,
	    char *which, long *nev,
	    double *tol, double *resid,
	    long *ncv, double *v, long *ldv,
	    long *iparam, long *ipntr,
	    double *workd, double *workl, long *lworkl, long *info);


/*********************************************************************************/
/*  wrapper for arpack routine                                                   */
/*********************************************************************************/
long dsaupd_lm_wrap(long n, const double *H, double *Z,
		   double *s, long nev){

  long dbuf_len;
  double *dbuf;
  long ibuf_len;
  long *ibuf;

  /* dsaupd */
  long ido = 0;
  char bmat[2] = "I"; /* standerd eigen problem */
  char which[3] = "LM"; /* Largest Magnitude */
  double tol = 0.0;
  long ldv;
  double *resid, *v, *d;
  double *workd, *workl;
  long ncv, lworkl, info, *iparam, *ipntr;
  /* dsaupd */
  long rvec = 1;
  char howmny[2] = "A";
  long *select, ierr;
  double sigma;

  /* number of the lanczos vector */
  ncv = nev*2;
  if(ncv > n){
    ncv = n;
  }

  ldv = n;
  lworkl = ncv*ncv + 8*ncv;
  dbuf_len = n + n*ncv + 3*n + lworkl + n;
  dbuf = (double *)malloc(dbuf_len*sizeof(double));
  resid = dbuf;
  v = resid + n;
  workd = v + n*ncv;
  workl = workd + 3*n;
  d = workl + lworkl;

  ibuf_len = 11 + 11 + ncv;
  ibuf = (long *)malloc(ibuf_len*sizeof(long));
  iparam = ibuf;
  ipntr = iparam + 11;
  select = ipntr + 11;

  iparam[0] = 1;
  iparam[2] = 3*n; /* max iteration */
  iparam[6] = 1;   /* standerd eigen problem */

  /* initialize */
  info = 0;

  do {
    /* reverse communication */
    dsaupd_(&ido, bmat, &n, which, &nev, &tol, resid,
	    &ncv, v, &ldv, iparam, ipntr,
	    workd, workl, &lworkl, &info);

    if ((ido == 1) || (ido == -1)){
      /* y = H*x */
      cblas_dsymv(CblasColMajor, CblasLower,
		  n, 1.0, H, n,
		  workd+(ipntr[0]-1), 1, 0.0, workd+(ipntr[1]-1), 1);
      /* z = B*x = x */
      if(ido == 1){
	memcpy(workd+(ipntr[2]-1), workd+(ipntr[0]-1), n*sizeof(double));
      }
    }
    else if (ido != 99){
      fprintf(stderr, "Error with (dsaupd_lm_wrap), ido = %ld \n", ido);
      fprintf(stderr, "Check documentation in dsaupd\n\n");
    }

  }while ((ido == 1) ||
	  (ido == -1));

  if (info < 0) {
    fprintf(stderr, "Error with dsaupd, info = %ld \n", info);
    fprintf(stderr, "Check documentation in dsaupd\n");
    if(info == -9999){
      fprintf(stderr, "%ld/%ld eigen value(s) are successfully converged. \n", iparam[4], nev);
    }
    fprintf(stderr, "\n");
    exit(1);
  }
  else {
    dseupd_(&rvec, howmny, select, s, Z, &ldv, &sigma, bmat,
	    &n, which, &nev, &tol, resid, &ncv, v, &ldv,
	    iparam, ipntr, workd, workl, &lworkl, &ierr);
    if (ierr!=0) {
      fprintf(stderr, "Error with dseupd, info = %ld\n", ierr);
      fprintf(stderr, "Check the documentation of dseupd.\n\n");
    } else if (info==1) {
      fprintf(stderr, "Maximum number of iterations reached.\n\n");
    } else if (info==3) {
      fprintf(stderr, "No shifts could be applied during implicit\n");
      fprintf(stderr, "Arnoldi update, try increasing NCV.\n\n");
    }
  }

  free(dbuf);
  free(ibuf);


  return 0;
}

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

  int i;
  long nR, rank, buf_len;
  double *R, *Z, *s, *buf;

  nR = 8;
  rank = 3;
  buf_len = nR*nR + rank*nR + rank;
  buf = (double *)malloc(buf_len*sizeof(double));
  R = buf;
  Z = R + nR*nR;
  s = Z + rank*nR;

  /* The Rosser Matrix (8x8, Column Major) */
  R[0] =  611;
  R[1] =  196; R[9]  =  899;
  R[2] = -192; R[10] =  113; R[18] = 899;
  R[3] =  407; R[11] = -192; R[19] = 196; R[27] =  611;
  R[4] =   -8; R[12] =  -71; R[20] =  61; R[28] =    8; R[36] =  411;
  R[5] =  -52; R[13] =  -43; R[21] =  49; R[29] =   44; R[37] = -599; R[45] = 411;
  R[6] =  -49; R[14] =   -8; R[22] =   8; R[30] =   59; R[38] =  208; R[46] = 208;
  R[7] =   29; R[15] =  -44; R[23] =  52; R[31] =  -23; R[39] =  208; R[47] = 208;

  R[54] =   99;
  R[55] = -911; R[63] = 99;

  /* Solve eigen problem of the Rosser matrix using arpack */
  dsaupd_lm_wrap(nR, R, Z, s, rank);

  /* Display eigen values */
  for(i=0; i< rank; i++){
    printf("%lf ", s[i]);
  }
  printf("\n");

  free(buf);

  return 0;
}

>Fix:


>Release-Note:
>Audit-Trail:
>Unformatted:



More information about the freebsd-ports-bugs mailing list