Franck HOUSSEN
aa970cefa6
The idea is:
- cmake generates arpackdef.h from arpackdef.cmake.h.in
- autoheader generates arpackdef.h from arpackdef.autotools.h[.in]
- arpack includes arpackdef.h when/where needed:
- #define provide c_int/a_int according to architecture.
- all f90 who need to include "arpackdef.h" must be moved to F90.
Note: by convention, F90 are preprocessed (f90 are not).
- MPI does not support ILP64: integer*4 must be imposed in all
calls involving MPI (f90 example/test code).
To enable ILP64 users to compile/link, arpackdef.h is added in
the arpack installation (make install).
103 lines
2.6 KiB
C
103 lines
2.6 KiB
C
#include "arpackdef.h"
|
|
|
|
#include <stdio.h>
|
|
#include <stdlib.h>
|
|
#include <math.h>
|
|
|
|
// This test calls fortran from C the old-fashion cumbersome way.
|
|
// Note: icb_arpack_c tests the same kind of things using ICB.
|
|
|
|
#ifdef INCLUDE_FCMANGLE
|
|
#include "FCMangle.h"
|
|
#endif
|
|
|
|
/* test program to solve for the 9 largest eigenvalues of
|
|
* A*x = lambda*x where A is the diagonal matrix
|
|
* with entries 1000, 999, ... , 2, 1 on the diagonal.
|
|
* We're using the non symmetric routines dnaupd and dneupd.
|
|
* This is not efficient since the problem is
|
|
* symmetric but is done to exhibit the bug.
|
|
*/
|
|
|
|
extern void snaupd(a_int *, char *, a_int *, char *, a_int *,
|
|
float *, float *, a_int *, float *,
|
|
a_int *, a_int *, a_int *, float *,
|
|
float *, a_int *, a_int *);
|
|
|
|
extern void sneupd( a_int*, char*, a_int *, float *, float *, float *, a_int*, float *,
|
|
float *, float *, char *, a_int *, char *, a_int *, float *, float *, a_int *,
|
|
float *, a_int *, a_int *, a_int *, float *, float *, a_int *, a_int * );
|
|
|
|
|
|
void matVec(float * x, float * y) {
|
|
int i;
|
|
for ( i = 0; i < 1000; ++i)
|
|
y[i] = ((float) (i+1))*x[i];
|
|
};
|
|
|
|
int main() {
|
|
a_int ido = 0;
|
|
char bmat[] = "I";
|
|
a_int N = 1000;
|
|
char which[] = "LM";
|
|
a_int nev = 9;
|
|
float tol = 0;
|
|
float resid[N];
|
|
a_int ncv = 2*nev+1;
|
|
float V[ncv*N];
|
|
a_int ldv = N;
|
|
a_int iparam[11];
|
|
a_int ipntr[14];
|
|
float workd[3*N];
|
|
a_int rvec = 1;
|
|
char howmny[] = "A";
|
|
float* dr = (float*) malloc((nev+1)*sizeof(float));
|
|
float* di = (float*) malloc((nev+1)*sizeof(float));
|
|
a_int select[3*ncv];
|
|
float z[(N+1)*(nev+1)];
|
|
a_int ldz = N+1;
|
|
float sigmar=0;
|
|
float sigmai=0;
|
|
float workev[3*ncv];
|
|
int k;
|
|
for (k=0; k < 3*N; ++k )
|
|
workd[k] = 0;
|
|
float workl[3*(ncv*ncv) + 6*ncv];
|
|
for (k=0; k < 3*(ncv*ncv) + 6*ncv; ++k )
|
|
workl[k] = 0;
|
|
a_int lworkl = 3*(ncv*ncv) + 6*ncv;
|
|
a_int info = 0;
|
|
|
|
iparam[0] = 1;
|
|
iparam[2] = 10*N;
|
|
iparam[3] = 1;
|
|
iparam[6] = 1;
|
|
|
|
snaupd(&ido, bmat, &N, which, &nev, &tol, resid, &ncv, V, &ldv, iparam, ipntr,
|
|
workd, workl, &lworkl, &info);
|
|
|
|
while(ido == -1 || ido == 1) {
|
|
|
|
matVec(&(workd[ipntr[0]-1]), &(workd[ipntr[1]-1]));
|
|
|
|
snaupd(&ido, bmat, &N, which, &nev, &tol, resid, &ncv, V, &ldv, iparam, ipntr,
|
|
workd, workl, &lworkl, &info);
|
|
}
|
|
|
|
sneupd( &rvec, howmny, select, dr,di, z, &ldz, &sigmar, &sigmai,workev,
|
|
bmat, &N, which, &nev, &tol, resid, &ncv, V, &ldv, iparam, ipntr,
|
|
workd, workl, &lworkl, &info);
|
|
int i;
|
|
for (i = 0; i < nev; ++i) {
|
|
printf("%f\n", dr[i]);
|
|
if(fabs(dr[i] - (float)(1000-i))>1e-2){
|
|
free(dr);
|
|
free(di);
|
|
exit(EXIT_FAILURE);
|
|
}
|
|
}
|
|
free(dr);
|
|
free(di);
|
|
return 0;
|
|
}
|