diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..9eca6c8 --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +*.a +*.o diff --git a/EXAMPLES/C/main_all.c b/EXAMPLES/C/main_all.c index fe669bc..affa6e9 100644 --- a/EXAMPLES/C/main_all.c +++ b/EXAMPLES/C/main_all.c @@ -12,30 +12,30 @@ #include "mpi.h" #include "pmrrr.h" -static int read_tri_mat(char*, double**, double**); -static void print_vector(char*, double*, char*, int); -static void print_matrix(char*, double*, int, int, int); +static PMRRR_Int read_tri_mat(char*, double**, double**); +static void print_vector(char*, double*, char*, PMRRR_Int); +static void print_matrix(char*, double*, PMRRR_Int, PMRRR_Int, PMRRR_Int); int main(int argc, char **argv) { /* Input parameter to 'pmrrr' */ - int n; /* Matrix size */ - int il, iu; - int tryRAC = 1; /* Try high rel. accuracy */ - double *D, *E; /* Diagonal and off-diagonal elements */ - double vl, vu; - - double *W; /* eigenvalues */ - int nz; /* # local eigenvectors */ - int offset; - double *Z; /* eigenvectors; stored by colums */ - int ldz; - int *Zsupp; /* eigenvector support */ + PMRRR_Int n; /* Matrix size */ + PMRRR_Int il, iu; + PMRRR_Int tryRAC = 1; /* Try high rel. accuracy */ + double *D, *E; /* Diagonal and off-diagonal elements */ + double vl, vu; + + double *W; /* eigenvalues */ + PMRRR_Int nz; /* # local eigenvectors */ + PMRRR_Int offset; + double *Z; /* eigenvectors; stored by colums */ + PMRRR_Int ldz; + PMRRR_Int *Zsupp; /* eigenvector support */ /* Others */ - int pid, nproc, info, status; - int i; + int pid, nproc, status; + PMRRR_Int info, i; MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &status); MPI_Comm_rank(MPI_COMM_WORLD, &pid); @@ -55,9 +55,9 @@ int main(int argc, char **argv) /* Allocate memory */ W = (double *) malloc( n * sizeof(double) ); - Zsupp = (int *) malloc( 2*n * sizeof(int) ); + Zsupp = (PMRRR_Int *) malloc( 2*n * sizeof(PMRRR_Int) ); - nz = (int) ceil(n / (double) nproc); + nz = (PMRRR_Int) ceil(n / (double) nproc); ldz = n; Z = (double *) malloc((size_t) n * nz * sizeof(double) ); @@ -101,9 +101,9 @@ int main(int argc, char **argv) /* * Reads the triadiagonal matrix from a file. */ -static int read_tri_mat(char *filename, double **Dp, double **Ep) +static PMRRR_Int read_tri_mat(char *filename, double **Dp, double **Ep) { - int i, n; + PMRRR_Int i, n; FILE *filedes; filedes = fopen(filename, "r"); @@ -129,9 +129,9 @@ static int read_tri_mat(char *filename, double **Dp, double **Ep) -static void print_vector(char *pre, double *v, char *post, int n) +static void print_vector(char *pre, double *v, char *post, PMRRR_Int n) { - int i; + PMRRR_Int i; printf("\n%s\n", pre); for (i=0; i (b) ? (a) : (b) ) #define fmin(a,b) ( (a) < (b) ? (a) : (b) ) +#ifndef __cplusplus typedef int bool; #define false 0 #define true 1 +#endif /* __cplusplus */ #endif #else /* __STDC__ not defined */ @@ -70,16 +74,20 @@ typedef int bool; /* in case compiler does not support type qualifiers * see Harbison and Steele p. 89*/ +#ifndef __cplusplus #define const /*nothing*/ #define volatile /*nothing*/ #define restrict /*nothing*/ +#endif /* __cplusplus */ #define fmax(a,b) ( (a) > (b) ? (a) : (b) ) #define fmin(a,b) ( (a) < (b) ? (a) : (b) ) +#ifndef __cplusplus typedef int bool; #define false 0 #define true 1 +#endif /* __cplusplus */ /* need also remove all // style comments */ /* inline ... */ @@ -92,5 +100,13 @@ typedef int bool; #define imin(a,b) ( (a) < (b) ? (a) : (b) ) #define iceil(a,b) ( (((a)%(b))!=0) ? (((a)/(b))+1) : ((a)/(b)) ) +#ifdef MPRRR_INT64 +typedef int64_t PMRRR_Int; +#define PMRRR_MPI_INT_TYPE MPI_LONG_LONG_INT +#else +typedef int PMRRR_Int; +#define PMRRR_MPI_INT_TYPE MPI_INT +#endif /* MPRRR_INT64 */ + /* End of header file */ #endif diff --git a/INCLUDE/plarre.h b/INCLUDE/plarre.h index 66a4707..e19f855 100644 --- a/INCLUDE/plarre.h +++ b/INCLUDE/plarre.h @@ -51,8 +51,8 @@ * Note: this implementation is here is not really optimized in * terms of performance and memory usage. */ -int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, - val_t *Wstruct, tol_t *tolstruct, int *nzp, int *offsetp); +PMRRR_Int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, + val_t *Wstruct, tol_t *tolstruct, PMRRR_Int *nzp, PMRRR_Int *offsetp); /* Perturb the initial root representation by "1 + eps*RAND_FACTOR*rand"; * default: 8.0 */ diff --git a/INCLUDE/plarrv.h b/INCLUDE/plarrv.h index 8d921aa..aea14d5 100644 --- a/INCLUDE/plarrv.h +++ b/INCLUDE/plarrv.h @@ -47,9 +47,9 @@ /* * Computation of eigenvectors of a symmetric tridiagonal */ -int plarrv(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, - vec_t *Zstruct, tol_t *tolstruct, int *nzp, - int *myfirstp); +PMRRR_Int plarrv(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, + vec_t *Zstruct, tol_t *tolstruct, PMRRR_Int *nzp, + PMRRR_Int *myfirstp); #define COMM_COMPLETE 0 #define COMM_INCOMPLETE 1 diff --git a/INCLUDE/pmrrr.h b/INCLUDE/pmrrr.h index d0c652d..e613cca 100644 --- a/INCLUDE/pmrrr.h +++ b/INCLUDE/pmrrr.h @@ -55,10 +55,10 @@ * * Function prototype: */ -int pmrrr(char *jobz, char *range, int *n, double *D, - double *E, double *vl, double *vu, int *il, int *iu, - int *tryrac, MPI_Comm comm, int *nz, int *offset, - double *W, double *Z, int *ldz, int *Zsupp); +PMRRR_Int pmrrr(char *jobz, char *range, PMRRR_Int *n, double *D, + double *E, double *vl, double *vu, PMRRR_Int *il, PMRRR_Int *iu, + PMRRR_Int *tryrac, MPI_Comm comm, PMRRR_Int *nz, PMRRR_Int *offset, + double *W, double *Z, PMRRR_Int *ldz, PMRRR_Int *Zsupp); /* Arguments: * ---------- @@ -167,7 +167,7 @@ int pmrrr(char *jobz, char *range, int *n, double *D, * EXAMPLE CALL: * ------------- * char *jobz, *range; - * int n, il, iu, tryRAC=0, nz, offset, ldz, *Zsupp; + * PMRRR_Int n, il, iu, tryRAC=0, nz, offset, ldz, *Zsupp; * double *D, *E, *W, *Z, vl, vu; * * // allocate space for D, E, W, Z @@ -229,7 +229,7 @@ int pmrrr(char *jobz, char *range, int *n, double *D, * all computed eigenvalues (iu-il+1) in W; this routine is designed * to be called right after 'pmrrr'. */ -int PMR_comm_eigvals(MPI_Comm comm, int *nz, int *ifirst, double *W); +PMRRR_Int PMR_comm_eigvals(MPI_Comm comm, PMRRR_Int *nz, PMRRR_Int *ifirst, double *W); /* Arguments: * ---------- * @@ -262,43 +262,45 @@ int PMR_comm_eigvals(MPI_Comm comm, int *nz, int *ifirst, double *W); /* LAPACK and BLAS function prototypes * Note: type specifier 'extern' does not matter in declaration * so here used to mark routines from LAPACK and BLAS libraries */ -extern void pmrrr_dscal(int*, double*, double*, int*); +extern void pmrrr_dscal(PMRRR_Int*, double*, double*, PMRRR_Int*); -extern double odnst_(char*, int*, double*, double*); -extern void odrrr_(int*, double*, double*, int*); -extern void odrra_(int*, double*, double*, double*, double*, - double*, int*, int*, int*); -extern void odrrc_(char*, int*, double*, double*, double*, double*, - double*, int*, int*, int*, int*); -extern void odrrd_(char*, char*, int*, double*, double*, int*, - int*, double*, double*, double*, double*, - double*, double*, int*, int*, int*, double*, - double*, double*, double*, int*, int*, double*, - int*, int*); -extern void odrrb_(int*, double*, double*, int*, int*, double*, - double*, int*, double*, double*, double*, double*, - int*, double*, double*, int*, int*); -extern void odrrk_(int*, int*, double*, double*, double*, double*, - double*, double*, double*, double*, int*); -extern void odebz_(int*, int*, int*, int*, int*, int*, double*, +extern double odnst_(char*, PMRRR_Int*, double*, double*); +extern void odrrr_(PMRRR_Int*, double*, double*, PMRRR_Int*); +extern void odrra_(PMRRR_Int*, double*, double*, double*, double*, + double*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*); +extern void odrrc_(char*, PMRRR_Int*, double*, double*, double*, double*, + double*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*); +extern void odrrd_(char*, char*, PMRRR_Int*, double*, double*, PMRRR_Int*, + PMRRR_Int*, double*, double*, double*, double*, + double*, double*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, double*, + double*, double*, double*, PMRRR_Int*, PMRRR_Int*, double*, + PMRRR_Int*, PMRRR_Int*); +extern void odrrb_(PMRRR_Int*, double*, double*, PMRRR_Int*, PMRRR_Int*, double*, + double*, PMRRR_Int*, double*, double*, double*, double*, + PMRRR_Int*, double*, double*, PMRRR_Int*, PMRRR_Int*); +extern void odrrk_(PMRRR_Int*, PMRRR_Int*, double*, double*, double*, double*, + double*, double*, double*, double*, PMRRR_Int*); +extern void odebz_(PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, double*, double*, double*, double*, double*, double*, - int*, double*, double*, int*, int*, double*, - int*, int*); -extern void odrnv_(int*, int*, int*, double*); -extern void odrrf_(int*, double*, double*, double*, int*, int*, + PMRRR_Int*, double*, double*, PMRRR_Int*, PMRRR_Int*, double*, + PMRRR_Int*, PMRRR_Int*); +extern void odrnv_(PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, double*); +extern void odrrf_(PMRRR_Int*, double*, double*, double*, PMRRR_Int*, PMRRR_Int*, double*, double*, double*, double*, double*, double*, double*, double*, double*, double*, - double*, int*); -extern void odr1v_(int*, int*, int*, double*, double*, double*, + double*, PMRRR_Int*); +extern void odr1v_(PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, double*, double*, double*, double*, double*, double*, double*, double*, - bool*, int*, double*, double*, int*, int*, + bool*, PMRRR_Int*, double*, double*, PMRRR_Int*, PMRRR_Int*, double*, double*, double*, double*); -extern void odrrj_(int*, double*, double*, int*, int*, double*, - int*, double*, double*, double*, int*, double*, - double*, int*); -extern void odstmr_(char*, char*, int*, double*, double*, double*, - double*, int*, int*, int*, double*, double*, - int*, int*, int*, int*, double*, int*, int*, - int*, int*); +extern void odrrj_(PMRRR_Int*, double*, double*, PMRRR_Int*, PMRRR_Int*, double*, + PMRRR_Int*, double*, double*, double*, PMRRR_Int*, double*, + double*, PMRRR_Int*); +extern void odstmr_(char*, char*, PMRRR_Int*, double*, double*, double*, + double*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, double*, double*, + PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, double*, PMRRR_Int*, PMRRR_Int*, + PMRRR_Int*, PMRRR_Int*); + +extern PMRRR_Int odscl_(PMRRR_Int *n, double *da, double *dx, PMRRR_Int *incx); #endif /* End of header file */ diff --git a/INCLUDE/process_task.h b/INCLUDE/process_task.h index 5add95c..c8bc69a 100644 --- a/INCLUDE/process_task.h +++ b/INCLUDE/process_task.h @@ -48,23 +48,23 @@ #include "queue.h" /* Function prototypes */ -int PMR_process_c_task(cluster_t *cl, int tid, proc_t *procinfo, +PMRRR_Int PMR_process_c_task(cluster_t *cl, PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, workQ_t *workQ, - counter_t *num_left, double *work, int *iwork); + counter_t *num_left, double *work, PMRRR_Int *iwork); -int PMR_process_s_task(singleton_t *sng, int tid, proc_t *procinfo, +PMRRR_Int PMR_process_s_task(singleton_t *sng, PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, counter_t *num_left, - double *work, int *iwork); + double *work, PMRRR_Int *iwork); -int PMR_process_r_task(refine_t *rf, proc_t *procinfo, val_t *Wstruct, - tol_t *tolstruct, double *work, int *iwork); +PMRRR_Int PMR_process_r_task(refine_t *rf, proc_t *procinfo, val_t *Wstruct, + tol_t *tolstruct, double *work, PMRRR_Int *iwork); -void PMR_process_r_queue(int tid, proc_t *procinfo, val_t *Wstruct, +void PMR_process_r_queue(PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, workQ_t *workQ, counter_t *num_left, double *work, - int *iwork); + PMRRR_Int *iwork); #endif /* End of header file */ diff --git a/INCLUDE/queue.h b/INCLUDE/queue.h index e9ac9bc..b542421 100644 --- a/INCLUDE/queue.h +++ b/INCLUDE/queue.h @@ -47,13 +47,13 @@ typedef struct task_aux task_t; struct task_aux { void *data; /* ptr to data, has to be casted */ - int flag; /* flag specifying the task */ + PMRRR_Int flag; /* flag specifying the task */ task_t *next; /* ptr to next task; NULL if non-existing; */ task_t *prev; /* ptr to prev. task; NULL if non-existing; */ }; typedef struct { - int num_tasks; + PMRRR_Int num_tasks; task_t *head; task_t *back; #ifdef NOSPINLOCKS @@ -66,11 +66,11 @@ typedef struct { /* functionality of the queue */ queue_t *PMR_create_empty_queue (void); -int PMR_insert_task_at_front (queue_t *queue, task_t *task); -int PMR_insert_task_at_back (queue_t *queue, task_t *task); +PMRRR_Int PMR_insert_task_at_front (queue_t *queue, task_t *task); +PMRRR_Int PMR_insert_task_at_back (queue_t *queue, task_t *task); task_t *PMR_remove_task_at_front(queue_t *queue); task_t *PMR_remove_task_at_back (queue_t *queue); -int PMR_get_num_tasks(queue_t *queue); +PMRRR_Int PMR_get_num_tasks(queue_t *queue); void PMR_destroy_queue(queue_t *queue); #endif diff --git a/INCLUDE/rrr.h b/INCLUDE/rrr.h index c93e6ac..7a3543c 100644 --- a/INCLUDE/rrr.h +++ b/INCLUDE/rrr.h @@ -49,26 +49,26 @@ typedef struct { double *restrict L; double *restrict DL; double *restrict DLL; - int size; - int depth; + PMRRR_Int size; + PMRRR_Int depth; bool parent_processed; bool copied_parent_rrr; - int ndepend; + PMRRR_Int ndepend; pthread_mutex_t mutex; } rrr_t; rrr_t *PMR_create_rrr(double *restrict D, double *restrict L, double *restrict DL, double *restrict DLL, - int size, int depth); + PMRRR_Int size, PMRRR_Int depth); rrr_t *PMR_reset_rrr (rrr_t *restrict RRR, double *restrict D, double *restrict L, double *restrict DL, - double *restrict DLL, int size, int depth); + double *restrict DLL, PMRRR_Int size, PMRRR_Int depth); -int PMR_increment_rrr_dependencies(rrr_t *RRR); -int PMR_set_parent_processed_flag (rrr_t *RRR); -int PMR_set_copied_parent_rrr_flag(rrr_t *RRR, bool val); -int PMR_try_destroy_rrr(rrr_t *RRR); +PMRRR_Int PMR_increment_rrr_dependencies(rrr_t *RRR); +PMRRR_Int PMR_set_parent_processed_flag (rrr_t *RRR); +PMRRR_Int PMR_set_copied_parent_rrr_flag(rrr_t *RRR, bool val); +PMRRR_Int PMR_try_destroy_rrr(rrr_t *RRR); #endif diff --git a/INCLUDE/structs.h b/INCLUDE/structs.h index e40ef35..f3201c7 100644 --- a/INCLUDE/structs.h +++ b/INCLUDE/structs.h @@ -47,44 +47,44 @@ #include "queue.h" typedef struct { - int n; + PMRRR_Int n; double *restrict D; double *restrict E; - int nsplit; - int *restrict isplit ; + PMRRR_Int nsplit; + PMRRR_Int *restrict isplit ; double spdiam; } in_t; typedef struct { - int n; + PMRRR_Int n; double *vl; double *vu; - int *il; - int *iu; + PMRRR_Int *il; + PMRRR_Int *iu; double *restrict W; double *restrict Werr; double *restrict Wgap; - int *restrict Windex; - int *restrict iblock; - int *restrict iproc; + PMRRR_Int *restrict Windex; + PMRRR_Int *restrict iblock; + PMRRR_Int *restrict iproc; double *restrict Wshifted; double *restrict gersch; } val_t; typedef struct { - int ldz; - int nz; + PMRRR_Int ldz; + PMRRR_Int nz; double *restrict Z; - int *restrict Zsupp; - int *restrict Zindex; + PMRRR_Int *restrict Zsupp; + PMRRR_Int *restrict Zindex; } vec_t; typedef struct { - int pid; - int nproc; + PMRRR_Int pid; + PMRRR_Int nproc; MPI_Comm comm; - int nthreads; - int thread_support; + PMRRR_Int nthreads; + PMRRR_Int thread_support; } proc_t; typedef struct { @@ -95,7 +95,7 @@ typedef struct { } tol_t; typedef struct { - int num_messages; + PMRRR_Int num_messages; MPI_Request *requests; MPI_Status *stats; } comm_t; @@ -108,41 +108,41 @@ typedef struct { typedef struct { double lambda; - int local_ind; - int block_ind; - int ind; + PMRRR_Int local_ind; + PMRRR_Int block_ind; + PMRRR_Int ind; } sort_struct_t; typedef struct { - int n; + PMRRR_Int n; double *D; double *E; double *E2; - int il; - int iu; - int my_il; - int my_iu; - int nsplit; - int *isplit; + PMRRR_Int il; + PMRRR_Int iu; + PMRRR_Int my_il; + PMRRR_Int my_iu; + PMRRR_Int nsplit; + PMRRR_Int *isplit; double bsrtol; double pivmin; double *gersch; double *W; double *Werr; - int *Windex; - int *iblock; + PMRRR_Int *Windex; + PMRRR_Int *iblock; } auxarg1_t; typedef struct { - int bl_size; + PMRRR_Int bl_size; double *D; double *DE2; - int rf_begin; - int rf_end; + PMRRR_Int rf_begin; + PMRRR_Int rf_end; double *W; double *Werr; double *Wgap; - int *Windex; + PMRRR_Int *Windex; double rtol1; double rtol2; double pivmin; @@ -150,7 +150,7 @@ typedef struct { } auxarg2_t; typedef struct { - int tid; + PMRRR_Int tid; proc_t *procinfo; val_t *Wstruct; vec_t *Zstruct; diff --git a/INCLUDE/tasks.h b/INCLUDE/tasks.h index ec8f9df..2800266 100644 --- a/INCLUDE/tasks.h +++ b/INCLUDE/tasks.h @@ -52,28 +52,28 @@ #define REFINE_TASK_FLAG 2 typedef struct { - int begin; - int end; - int depth; - int bl_begin; - int bl_end; + PMRRR_Int begin; + PMRRR_Int end; + PMRRR_Int depth; + PMRRR_Int bl_begin; + PMRRR_Int bl_end; double bl_spdiam; double lgap; rrr_t *RRR; } singleton_t; typedef struct { - int begin; - int end; - int depth; - int bl_begin; /* In priciple not needed since info */ - int bl_end; /* also contained in iblock+isplit */ + PMRRR_Int begin; + PMRRR_Int end; + PMRRR_Int depth; + PMRRR_Int bl_begin; /* In priciple not needed since info */ + PMRRR_Int bl_end; /* also contained in iblock+isplit */ double bl_spdiam; double lgap; - int proc_W_begin; - int proc_W_end; - int left_pid; - int right_pid; + PMRRR_Int proc_W_begin; + PMRRR_Int proc_W_end; + PMRRR_Int left_pid; + PMRRR_Int right_pid; rrr_t *RRR; bool wait_until_refined; comm_t *messages; @@ -81,33 +81,33 @@ typedef struct { typedef struct { - int begin; - int end; + PMRRR_Int begin; + PMRRR_Int end; double *D; double *DLL; - int p; - int q; - int bl_size; + PMRRR_Int p; + PMRRR_Int q; + PMRRR_Int bl_size; double bl_spdiam; - int producer_tid; // not longer needed + PMRRR_Int producer_tid; // not longer needed sem_t *sem; /* since semt_t is a handle could also store it instead of pointer to it, but pointer is all that is needed */ } refine_t; -task_t *PMR_create_s_task(int first, int last, int depth, - int bl_begin, int bl_end, double spdiam, +task_t *PMR_create_s_task(PMRRR_Int first, PMRRR_Int last, PMRRR_Int depth, + PMRRR_Int bl_begin, PMRRR_Int bl_end, double spdiam, double lgap, rrr_t *RRR); -task_t *PMR_create_c_task(int first, int last, int depth, - int bl_begin, int bl_end, double spdiam, - double lgap, int proc_W_begin, - int proc_W_end, int left_pid, int right_pid, +task_t *PMR_create_c_task(PMRRR_Int first, PMRRR_Int last, PMRRR_Int depth, + PMRRR_Int bl_begin, PMRRR_Int bl_end, double spdiam, + double lgap, PMRRR_Int proc_W_begin, + PMRRR_Int proc_W_end, PMRRR_Int left_pid, PMRRR_Int right_pid, rrr_t *RRR); -task_t *PMR_create_r_task(int begin, int end, double *D, - double *DLL, int p, int q, int bl_size, - double bl_spdiam, int tid, sem_t *sem); +task_t *PMR_create_r_task(PMRRR_Int begin, PMRRR_Int end, double *D, + double *DLL, PMRRR_Int p, PMRRR_Int q, PMRRR_Int bl_size, + double bl_spdiam, PMRRR_Int tid, sem_t *sem); #endif diff --git a/README b/README index 5fc0ece..1446891 100644 --- a/README +++ b/README @@ -40,10 +40,12 @@ For more information please see 'INCLUDE/pmrrr.h'. ###################################################################### # C function prototype: # ###################################################################### -# int pmrrr(char *jobz, char *range, int *n, double *D, # -# double *E, double *vl, double *vu, int *il, int *iu, # -# int *tryrac, MPI_Comm comm, int *nz, int *offset, # -# double *W, double *Z, int *ldz, int *Zsupp); # +# PMRRR_Int pmrrr(char *jobz, char *range, PMRRR_Int *n, # +# double *D, double *E, double *vl, double *vu, # +# PMRRR_Int *il, PMRRR_Int *iu, # +# PMRRR_Int *tryrac, MPI_Comm comm, # +# PMRRR_Int *nz, PMRRR_Int *offset, # +# double *W, double *Z, PMRRR_Int *ldz, PMRRR_Int *Zsupp); # # # # Arguments: # # ---------- # @@ -172,6 +174,18 @@ Note: Do not forget to tell mpiexec/mpirun about the environment variable PMR_NUM_THREADS if you use it. This is typically done via the option '-env' or '-x'. +The 'PMRRR_Int' type is either 32-bit or 64-bit signed integer. +64-bit is used, if '-DPMRRR_INT64' compilation flag is set. +64-bit integers allow you to solve problems with more than 4 Gb +blocks per node. + +Please note that in Fortran, you must use 'INTEGER*8' for all the +pmrrr-related integers when the library is compiled with '-DPMRRR_INT64'. +In C/C++, you can use the 'PMRRR_Int' alias defined via 'pmrrr.h', but +you must also compile your program with or without '-DPMRRR_INT64', +depending on how the library was compiled. Maybe in the future, auto-generated +headers will be introduced to make this easier. + Some additinal comments: ------------------------ @@ -189,8 +203,8 @@ given below. For more information please see 'pmrrr.h'. ###################################################################### # C function prototype: # ###################################################################### -# int PMR_comm_eigvals(MPI_Comm comm, int *nz, int *offset, # -# double *W); # +# PMRRR_Int PMR_comm_eigvals(MPI_Comm comm, # +# PMRRR_Int *nz, PMRRR_Int *offset, double *W); # ###################################################################### ###################################################################### diff --git a/SRC/BLAS/odcpy.c b/SRC/BLAS/odcpy.c index d9a394f..3ce1009 100644 --- a/SRC/BLAS/odcpy.c +++ b/SRC/BLAS/odcpy.c @@ -6,16 +6,17 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odcpy_(int *n, double *dx, int *incx, - double *dy, int *incy) +PMRRR_Int odcpy_(PMRRR_Int *n, double *dx, PMRRR_Int *incx, + double *dy, PMRRR_Int *incy) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; /* Local variables */ - int i__, m, ix, iy, mp1; + PMRRR_Int i__, m, ix, iy, mp1; /* .. Scalar Arguments .. */ /* .. */ diff --git a/SRC/BLAS/odscl.c b/SRC/BLAS/odscl.c index f0e8a66..56c06c6 100644 --- a/SRC/BLAS/odscl.c +++ b/SRC/BLAS/odscl.c @@ -6,16 +6,17 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odscl_(int *n, double *da, double *dx, - int *incx) +PMRRR_Int odscl_(PMRRR_Int *n, double *da, double *dx, + PMRRR_Int *incx) { /* System generated locals */ - int i__1, i__2; + PMRRR_Int i__1, i__2; /* Local variables */ - int i__, m, mp1, nincx; + PMRRR_Int i__, m, mp1, nincx; /* .. Scalar Arguments .. */ /* .. */ diff --git a/SRC/BLAS/odswap.c b/SRC/BLAS/odswap.c index 688005c..c221609 100644 --- a/SRC/BLAS/odswap.c +++ b/SRC/BLAS/odswap.c @@ -6,16 +6,17 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odswap_(int *n, double *dx, int *incx, - double *dy, int *incy) +PMRRR_Int odswap_(PMRRR_Int *n, double *dx, PMRRR_Int *incx, + double *dy, PMRRR_Int *incy) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; /* Local variables */ - int i__, m, ix, iy, mp1; + PMRRR_Int i__, m, ix, iy, mp1; double dtemp; /* .. Scalar Arguments .. */ diff --git a/SRC/LAPACK/ode2.c b/SRC/LAPACK/ode2.c index 17c09cf..bb2bb57 100644 --- a/SRC/LAPACK/ode2.c +++ b/SRC/LAPACK/ode2.c @@ -6,9 +6,10 @@ #include #include #include +#include "global.h" /* Subroutine */ -int ode2_(double *a, double *b, double *c__, +PMRRR_Int ode2_(double *a, double *b, double *c__, double *rt1, double *rt2) { /* System generated locals */ diff --git a/SRC/LAPACK/odebz.c b/SRC/LAPACK/odebz.c index f80373d..fb83b72 100644 --- a/SRC/LAPACK/odebz.c +++ b/SRC/LAPACK/odebz.c @@ -6,27 +6,28 @@ #include #include #include +#include "global.h" #define imax(a,b) ( (a) > (b) ? (a) : (b) ) #define imin(a,b) ( (a) < (b) ? (a) : (b) ) /* Subroutine */ -int odebz_(int *ijob, int *nitmax, int *n, - int *mmax, int *minp, int *nbmin, double *abstol, +PMRRR_Int odebz_(PMRRR_Int *ijob, PMRRR_Int *nitmax, PMRRR_Int *n, + PMRRR_Int *mmax, PMRRR_Int *minp, PMRRR_Int *nbmin, double *abstol, double *reltol, double *pivmin, double *d__, double * - e, double *e2, int *nval, double *ab, double *c__, - int *mout, int *nab, double *work, int *iwork, - int *info) + e, double *e2, PMRRR_Int *nval, double *ab, double *c__, + PMRRR_Int *mout, PMRRR_Int *nab, double *work, PMRRR_Int *iwork, + PMRRR_Int *info) { /* System generated locals */ - int nab_dim1, nab_offset, ab_dim1, ab_offset, i__1, i__2, i__3, i__4, + PMRRR_Int nab_dim1, nab_offset, ab_dim1, ab_offset, i__1, i__2, i__3, i__4, i__5, i__6; double d__1, d__2, d__3, d__4; /* Local variables */ - int j, kf, ji, kl, jp, jit; + PMRRR_Int j, kf, ji, kl, jp, jit; double tmp1, tmp2; - int itmp1, itmp2, kfnew, klnew; + PMRRR_Int itmp1, itmp2, kfnew, klnew; /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odev2.c b/SRC/LAPACK/odev2.c index ba67ed5..ce1e7d2 100644 --- a/SRC/LAPACK/odev2.c +++ b/SRC/LAPACK/odev2.c @@ -6,9 +6,10 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odev2_(double *a, double *b, double *c__, +PMRRR_Int odev2_(double *a, double *b, double *c__, double *rt1, double *rt2, double *cs1, double *sn1) { /* System generated locals */ @@ -19,7 +20,7 @@ int odev2_(double *a, double *b, double *c__, /* Local variables */ double ab, df, cs, ct, tb, sm, tn, rt, adf, acs; - int sgn1, sgn2; + PMRRR_Int sgn1, sgn2; double acmn, acmx; diff --git a/SRC/LAPACK/odnan.c b/SRC/LAPACK/odnan.c index 2d3aa54..387be15 100644 --- a/SRC/LAPACK/odnan.c +++ b/SRC/LAPACK/odnan.c @@ -6,14 +6,15 @@ #include #include #include +#include "global.h" -int odnan_(double *din) +PMRRR_Int odnan_(double *din) { /* System generated locals */ - int ret_val; + PMRRR_Int ret_val; /* Local variables */ - extern int odsnan_(double *, double *); + extern PMRRR_Int odsnan_(double *, double *); /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odneg.c b/SRC/LAPACK/odneg.c index f5c63cd..7ebd5cb 100644 --- a/SRC/LAPACK/odneg.c +++ b/SRC/LAPACK/odneg.c @@ -6,26 +6,27 @@ #include #include #include +#include "global.h" #define imax(a,b) ( (a) > (b) ? (a) : (b) ) #define imin(a,b) ( (a) < (b) ? (a) : (b) ) -int odneg_(int *n, double *d__, double *lld, double * - sigma, double *pivmin, int *r__) +PMRRR_Int odneg_(PMRRR_Int *n, double *d__, double *lld, double * + sigma, double *pivmin, PMRRR_Int *r__) { /* System generated locals */ - int ret_val, i__1, i__2, i__3, i__4; + PMRRR_Int ret_val, i__1, i__2, i__3, i__4; /* Local variables */ - int j; + PMRRR_Int j; double p, t; - int bj; + PMRRR_Int bj; double tmp; - int neg1, neg2; + PMRRR_Int neg1, neg2; double bsav, gamma, dplus; - extern int odnan_(double *); - int negcnt; - int sawnan; + extern PMRRR_Int odnan_(double *); + PMRRR_Int negcnt; + PMRRR_Int sawnan; double dminus; diff --git a/SRC/LAPACK/odnst.c b/SRC/LAPACK/odnst.c index f18dce8..8dbb6b6 100644 --- a/SRC/LAPACK/odnst.c +++ b/SRC/LAPACK/odnst.c @@ -6,25 +6,26 @@ #include #include #include +#include "global.h" /* Table of constant values */ -static int c__1 = 1; +static PMRRR_Int c__1 = 1; -double odnst_(char *norm, int *n, double *d__, double *e) +double odnst_(char *norm, PMRRR_Int *n, double *d__, double *e) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double ret_val, d__1, d__2, d__3, d__4, d__5; /* Builtin functions */ // double sqrt(double); /* Local variables */ - int i__; + PMRRR_Int i__; double sum, scale; - extern int olsame_(char *, char *); + extern PMRRR_Int olsame_(char *, char *); double anorm; - extern /* Subroutine */ int odssq_(int *, double *, int *, + extern /* Subroutine */ PMRRR_Int odssq_(PMRRR_Int *, double *, PMRRR_Int *, double *, double *); @@ -97,6 +98,8 @@ double odnst_(char *norm, int *n, double *d__, double *e) --e; --d__; + anorm = 0.0; + /* Function Body */ if (*n <= 0) { anorm = 0.; diff --git a/SRC/LAPACK/odr1v.c b/SRC/LAPACK/odr1v.c index dde8858..c2b5603 100644 --- a/SRC/LAPACK/odr1v.c +++ b/SRC/LAPACK/odr1v.c @@ -6,37 +6,38 @@ #include #include #include +#include "global.h" #define TRUE_ (1) #define FALSE_ (0) /* Subroutine */ -int odr1v_(int *n, int *b1, int *bn, double +PMRRR_Int odr1v_(PMRRR_Int *n, PMRRR_Int *b1, PMRRR_Int *bn, double *lambda, double *d__, double *l, double *ld, double * - lld, double *pivmin, double *gaptol, double *z__, int - *wantnc, int *negcnt, double *ztz, double *mingma, - int *r__, int *isuppz, double *nrminv, double *resid, + lld, double *pivmin, double *gaptol, double *z__, PMRRR_Int + *wantnc, PMRRR_Int *negcnt, double *ztz, double *mingma, + PMRRR_Int *r__, PMRRR_Int *isuppz, double *nrminv, double *resid, double *rqcorr, double *work) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2, d__3; /* Builtin functions */ // double sqrt(double); /* Local variables */ - int i__; + PMRRR_Int i__; double s; - int r1, r2; + PMRRR_Int r1, r2; double eps, tmp; - int neg1, neg2, indp, inds; + PMRRR_Int neg1, neg2, indp, inds; double dplus; // extern double odmch_(char *); - extern int odnan_(double *); - int indlpl, indumn; + extern PMRRR_Int odnan_(double *); + PMRRR_Int indlpl, indumn; double dminus; - int sawnan1, sawnan2; + PMRRR_Int sawnan1, sawnan2; /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odrnv.c b/SRC/LAPACK/odrnv.c index cff0093..f29ddde 100644 --- a/SRC/LAPACK/odrnv.c +++ b/SRC/LAPACK/odrnv.c @@ -6,24 +6,25 @@ #include #include #include +#include "global.h" #define imin(a,b) ( (a) < (b) ? (a) : (b) ) /* Subroutine */ -int odrnv_(int *idist, int *iseed, int *n, +PMRRR_Int odrnv_(PMRRR_Int *idist, PMRRR_Int *iseed, PMRRR_Int *n, double *x) { /* System generated locals */ - int i__1, i__2, i__3; + PMRRR_Int i__1, i__2, i__3; /* Builtin functions */ // double log(double), sqrt(double), cos(double); /* Local variables */ - int i__; + PMRRR_Int i__; double u[128]; - int il, iv, il2; - extern int odruv_(int *, int *, double *); + PMRRR_Int il, iv, il2; + extern PMRRR_Int odruv_(PMRRR_Int *, PMRRR_Int *, double *); /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odrra.c b/SRC/LAPACK/odrra.c index 8eb8b05..a824205 100644 --- a/SRC/LAPACK/odrra.c +++ b/SRC/LAPACK/odrra.c @@ -6,21 +6,22 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odrra_(int *n, double *d__, double *e, - double *e2, double *spltol, double *tnrm, int *nsplit, - int *isplit, int *info) +PMRRR_Int odrra_(PMRRR_Int *n, double *d__, double *e, + double *e2, double *spltol, double *tnrm, PMRRR_Int *nsplit, + PMRRR_Int *isplit, PMRRR_Int *info) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Builtin functions */ // double sqrt(double); /* Local variables */ - int i__; + PMRRR_Int i__; double tmp1, eabs; diff --git a/SRC/LAPACK/odrrb.c b/SRC/LAPACK/odrrb.c index 9403123..ccc3371 100644 --- a/SRC/LAPACK/odrrb.c +++ b/SRC/LAPACK/odrrb.c @@ -6,31 +6,32 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odrrb_(int *n, double *d__, double *lld, - int *ifirst, int *ilast, double *rtol1, double *rtol2, - int *offset, double *w, double *wgap, double *werr, - double *work, int *iwork, double *pivmin, double * - spdiam, int *twist, int *info) +PMRRR_Int odrrb_(PMRRR_Int *n, double *d__, double *lld, + PMRRR_Int *ifirst, PMRRR_Int *ilast, double *rtol1, double *rtol2, + PMRRR_Int *offset, double *w, double *wgap, double *werr, + double *work, PMRRR_Int *iwork, double *pivmin, double * + spdiam, PMRRR_Int *twist, PMRRR_Int *info) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Builtin functions */ // double log(double); /* Local variables */ - int i__, k, r__, i1, ii, ip; + PMRRR_Int i__, k, r__, i1, ii, ip; double gap, mid, tmp, back, lgap, rgap, left; - int iter, nint, prev, next; + PMRRR_Int iter, nint, prev, next; double cvrgd, right, width; - extern int odneg_(int *, double *, double *, double * -, double *, int *); - int negcnt; + extern PMRRR_Int odneg_(PMRRR_Int *, double *, double *, double * +, double *, PMRRR_Int *); + PMRRR_Int negcnt; double mnwdth; - int olnint, maxitr; + PMRRR_Int olnint, maxitr; /* -- LAPACK auxiliary routine (version 3.2) -- */ /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ @@ -158,7 +159,7 @@ int odrrb_(int *n, double *d__, double *lld, /* Function Body */ *info = 0; - maxitr = (int) ((log(*spdiam + *pivmin) - log(*pivmin)) / log(2.)) + + maxitr = (PMRRR_Int) ((log(*spdiam + *pivmin) - log(*pivmin)) / log(2.)) + 2; mnwdth = *pivmin * 2.; @@ -169,7 +170,7 @@ int odrrb_(int *n, double *d__, double *lld, /* Initialize unconverged intervals in [ WORK(2*I-1), WORK(2*I) ]. */ /* The Sturm Count, Count( WORK(2*I-1) ) is arranged to be I-1, while */ -/* Count( WORK(2*I) ) is stored in IWORK( 2*I ). The int IWORK( 2*I-1 ) */ +/* Count( WORK(2*I) ) is stored in IWORK( 2*I ). The PMRRR_Int IWORK( 2*I-1 ) */ /* for an unconverged interval is set to the index of the next unconverged */ /* interval, and is -1 or 0 for a converged interval. Thus a linked */ /* list of unconverged intervals is set up. */ diff --git a/SRC/LAPACK/odrrc.c b/SRC/LAPACK/odrrc.c index f0c92ca..157ac40 100644 --- a/SRC/LAPACK/odrrc.c +++ b/SRC/LAPACK/odrrc.c @@ -6,22 +6,23 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odrrc_(char *jobt, int *n, double *vl, +PMRRR_Int odrrc_(char *jobt, PMRRR_Int *n, double *vl, double *vu, double *d__, double *e, double *pivmin, - int *eigcnt, int *lcnt, int *rcnt, int *info) + PMRRR_Int *eigcnt, PMRRR_Int *lcnt, PMRRR_Int *rcnt, PMRRR_Int *info) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1; /* Local variables */ - int i__; + PMRRR_Int i__; double sl, su, tmp, tmp2; - int matt; - extern int olsame_(char *, char *); + PMRRR_Int matt; + extern PMRRR_Int olsame_(char *, char *); double lpivot, rpivot; diff --git a/SRC/LAPACK/odrrd.c b/SRC/LAPACK/odrrd.c index 3484312..0129db1 100644 --- a/SRC/LAPACK/odrrd.c +++ b/SRC/LAPACK/odrrd.c @@ -6,6 +6,7 @@ #include #include #include +#include "global.h" #define imax(a,b) ( (a) > (b) ? (a) : (b) ) #define imin(a,b) ( (a) < (b) ? (a) : (b) ) @@ -13,57 +14,56 @@ #define FALSE_ (0) /* Table of constant values */ -static int c__1 = 1; -static int c_n1 = -1; -static int c__3 = 3; -static int c__2 = 2; -static int c__0 = 0; +static PMRRR_Int c__1 = 1; +static PMRRR_Int c__3 = 3; +static PMRRR_Int c__2 = 2; +static PMRRR_Int c__0 = 0; /* Subroutine */ -int odrrd_(char *range, char *order, int *n, double *vl, - double *vu, int *il, int *iu, double *gers, +PMRRR_Int odrrd_(char *range, char *order, PMRRR_Int *n, double *vl, + double *vu, PMRRR_Int *il, PMRRR_Int *iu, double *gers, double *reltol, double *d__, double *e, double *e2, - double *pivmin, int *nsplit, int *isplit, int *m, + double *pivmin, PMRRR_Int *nsplit, PMRRR_Int *isplit, PMRRR_Int *m, double *w, double *werr, double *wl, double *wu, - int *iblock, int *indexw, double *work, int *iwork, - int *info) + PMRRR_Int *iblock, PMRRR_Int *indexw, double *work, PMRRR_Int *iwork, + PMRRR_Int *info) { /* System generated locals */ - int i__1, i__2, i__3; + PMRRR_Int i__1, i__2, i__3; double d__1, d__2; /* Builtin functions */ // double log(double); /* Local variables */ - int i__, j, ib, ie, je, nb; + PMRRR_Int i__, j, ib, ie, je, nb; double gl; - int im, in; + PMRRR_Int im, in; double gu; - int iw, jee; + PMRRR_Int iw, jee; double eps; - int nwl; + PMRRR_Int nwl; double wlu, wul; - int nwu; + PMRRR_Int nwu; double tmp1, tmp2; - int iend, jblk, ioff, iout, itmp1, itmp2, jdisc; - extern int olsame_(char *, char *); - int iinfo; + PMRRR_Int iend, jblk, ioff, iout, itmp1, itmp2, jdisc; + extern PMRRR_Int olsame_(char *, char *); + PMRRR_Int iinfo; double atoli; - int iwoff, itmax; + PMRRR_Int iwoff, itmax; double wkill, rtoli, uflow, tnorm; // extern double dlamch_(char *); - int ibegin; - extern /* Subroutine */ int odebz_(int *, int *, int *, - int *, int *, int *, double *, double *, - double *, double *, double *, double *, int *, - double *, double *, int *, int *, double *, - int *, int *); - int irange, idiscl, idumma[1]; - /* extern int ilaenv_(int *, char *, char *, int *, int *, */ - /* int *, int *); */ - int idiscu; - int ncnvrg, toofew; + PMRRR_Int ibegin; + extern /* Subroutine */ PMRRR_Int odebz_(PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, + PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, double *, double *, + double *, double *, double *, double *, PMRRR_Int *, + double *, double *, PMRRR_Int *, PMRRR_Int *, double *, + PMRRR_Int *, PMRRR_Int *); + PMRRR_Int irange, idiscl, idumma[1]; + /* extern PMRRR_Int ilaenv_(PMRRR_Int *, char *, char *, PMRRR_Int *, PMRRR_Int *, */ + /* PMRRR_Int *, PMRRR_Int *); */ + PMRRR_Int idiscu; + PMRRR_Int ncnvrg, toofew; /* -- LAPACK auxiliary routine (version 3.2.1) -- */ @@ -291,6 +291,8 @@ int odrrd_(char *range, char *order, int *n, double *vl, /* Function Body */ *info = 0; + wlu = wul = 0.0; + /* Decode RANGE */ if (olsame_(range, "A")) { @@ -378,8 +380,8 @@ int odrrd_(char *range, char *order, int *n, double *vl, /* Computing MAX */ d__1 = fabs(gl), d__2 = fabs(gu); tnorm = fmax(d__1,d__2); - gl = gl - tnorm * 2. * eps * *n - *pivmin * 4.; - gu = gu + tnorm * 2. * eps * *n + *pivmin * 4.; + gl = gl - tnorm * 2. * eps * (double)(*n) - *pivmin * 4.; + gu = gu + tnorm * 2. * eps * (double)(*n) + *pivmin * 4.; /* [JAN/28/2009] remove the line below since SPDIAM variable not use */ /* SPDIAM = GU - GL */ /* Input arguments for ODEBZ: */ @@ -397,7 +399,7 @@ int odrrd_(char *range, char *order, int *n, double *vl, /* RANGE='I': Compute an interval containing eigenvalues */ /* IL through IU. The initial interval [GL,GU] from the global */ /* Gerschgorin bounds GL and GU is refined by ODEBZ. */ - itmax = (int) ((log(tnorm + *pivmin) - log(*pivmin)) / log(2.)) + + itmax = (PMRRR_Int) ((log(tnorm + *pivmin) - log(*pivmin)) / log(2.)) + 2; work[*n + 1] = gl; work[*n + 2] = gl; @@ -551,8 +553,8 @@ int odrrd_(char *range, char *order, int *n, double *vl, /* SPDIAM = GU - GL */ /* GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN */ /* GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN */ - gl = gl - tnorm * 2. * eps * in - *pivmin * 2.; - gu = gu + tnorm * 2. * eps * in + *pivmin * 2.; + gl = gl - tnorm * 2. * eps * (double)in - *pivmin * 2.; + gu = gu + tnorm * 2. * eps * (double)in + *pivmin * 2.; if (irange > 1) { if (gu < *wl) { @@ -584,7 +586,7 @@ int odrrd_(char *range, char *order, int *n, double *vl, nwu += iwork[in + 1]; iwoff = *m - iwork[1]; /* Compute Eigenvalues */ - itmax = (int) ((log(gu - gl + *pivmin) - log(*pivmin)) / log( + itmax = (PMRRR_Int) ((log(gu - gl + *pivmin) - log(*pivmin)) / log( 2.)) + 2; odebz_(&c__2, &itmax, &in, &in, &c__1, &nb, &atoli, &rtoli, pivmin, &d__[ibegin], &e[ibegin], &e2[ibegin], idumma, & diff --git a/SRC/LAPACK/odrre.c b/SRC/LAPACK/odrre.c index 10b56a1..4fa0f2d 100644 --- a/SRC/LAPACK/odrre.c +++ b/SRC/LAPACK/odrre.c @@ -6,81 +6,82 @@ #include #include #include +#include "global.h" /* Table of constant values */ -static int c__1 = 1; -static int c__2 = 2; +static PMRRR_Int c__1 = 1; +static PMRRR_Int c__2 = 2; #define TRUE_ (1) #define FALSE_ (0) /* Subroutine */ -int odrre_(char *range, int *n, double *vl, - double *vu, int *il, int *iu, double *d__, double *e, +PMRRR_Int odrre_(char *range, PMRRR_Int *n, double *vl, + double *vu, PMRRR_Int *il, PMRRR_Int *iu, double *d__, double *e, double *e2, double *rtol1, double *rtol2, double *spltol, - int *nsplit, int *isplit, int *m, double *w, - double *werr, double *wgap, int *iblock, int *indexw, - double *gers, double *pivmin, double *work, int * - iwork, int *info) + PMRRR_Int *nsplit, PMRRR_Int *isplit, PMRRR_Int *m, double *w, + double *werr, double *wgap, PMRRR_Int *iblock, PMRRR_Int *indexw, + double *gers, double *pivmin, double *work, PMRRR_Int * + iwork, PMRRR_Int *info) { /* System generated locals */ - int i__1, i__2; + PMRRR_Int i__1, i__2; double d__1, d__2, d__3; /* Builtin functions */ // double sqrt(double), log(double); /* Local variables */ - int i__, j; + PMRRR_Int i__, j; double s1, s2; - int mb; + PMRRR_Int mb; double gl; - int in, mm; + PMRRR_Int in, mm; double gu; - int cnt; + PMRRR_Int cnt; double eps, tau, tmp, rtl; - int cnt1, cnt2; + PMRRR_Int cnt1, cnt2; double tmp1, eabs; - int iend, jblk; + PMRRR_Int iend, jblk; double eold; - int indl; + PMRRR_Int indl; double dmax__, emax; - int wend, idum, indu; + PMRRR_Int wend, idum, indu; double rtol; - int iseed[4]; + PMRRR_Int iseed[4]; double avgap, sigma; - extern int olsame_(char *, char *); - int iinfo; - extern /* Subroutine */ int odcpy_(int *, double *, int *, - double *, int *); + extern PMRRR_Int olsame_(char *, char *); + PMRRR_Int iinfo; + extern /* Subroutine */ PMRRR_Int odcpy_(PMRRR_Int *, double *, PMRRR_Int *, + double *, PMRRR_Int *); long double norep; - extern /* Subroutine */ int odsq2_(int *, double *, int *); + extern /* Subroutine */ PMRRR_Int odsq2_(PMRRR_Int *, double *, PMRRR_Int *); // extern double odmch_(char *); - int ibegin; + PMRRR_Int ibegin; long double forceb; - int irange; + PMRRR_Int irange; double sgndef; - extern /* Subroutine */ int odrra_(int *, double *, double *, - double *, double *, double *, int *, int *, - int *), odrrb_(int *, double *, double *, - int *, int *, double *, double *, int *, - double *, double *, double *, double *, int *, - double *, double *, int *, int *), odrrc_(char * -, int *, double *, double *, double *, double - *, double *, int *, int *, int *, int *); - int wbegin; - extern /* Subroutine */ int odrrd_(char *, char *, int *, double - *, double *, int *, int *, double *, double *, - double *, double *, double *, double *, int * -, int *, int *, double *, double *, double *, - double *, int *, int *, double *, int *, - int *); + extern /* Subroutine */ PMRRR_Int odrra_(PMRRR_Int *, double *, double *, + double *, double *, double *, PMRRR_Int *, PMRRR_Int *, + PMRRR_Int *), odrrb_(PMRRR_Int *, double *, double *, + PMRRR_Int *, PMRRR_Int *, double *, double *, PMRRR_Int *, + double *, double *, double *, double *, PMRRR_Int *, + double *, double *, PMRRR_Int *, PMRRR_Int *), odrrc_(char * +, PMRRR_Int *, double *, double *, double *, double + *, double *, PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, PMRRR_Int *); + PMRRR_Int wbegin; + extern /* Subroutine */ PMRRR_Int odrrd_(char *, char *, PMRRR_Int *, double + *, double *, PMRRR_Int *, PMRRR_Int *, double *, double *, + double *, double *, double *, double *, PMRRR_Int * +, PMRRR_Int *, PMRRR_Int *, double *, double *, double *, + double *, PMRRR_Int *, PMRRR_Int *, double *, PMRRR_Int *, + PMRRR_Int *); double safmin, spdiam; - extern /* Subroutine */ int odrrk_(int *, int *, double *, + extern /* Subroutine */ PMRRR_Int odrrk_(PMRRR_Int *, PMRRR_Int *, double *, double *, double *, double *, double *, - double *, double *, double *, int *); + double *, double *, double *, PMRRR_Int *); long double usedqd; double clwdth, isleft; - extern /* Subroutine */ int odrnv_(int *, int *, int *, + extern /* Subroutine */ PMRRR_Int odrnv_(PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, double *); double isrght, bsrtol, dpivot; @@ -282,6 +283,9 @@ int odrre_(char *range, int *n, double *vl, /* Function Body */ *info = 0; + mb = 0; + wend = 0; + irange = 0; /* Decode RANGE */ @@ -460,7 +464,7 @@ int odrre_(char *range, int *n, double *vl, goto L170; } else { /* Decide whether dqds or bisection is more efficient */ - usedqd = (double) mb > in * .5 && ! forceb; + usedqd = (double) mb > (double) in * .5 && ! forceb; wend = wbegin + mb - 1; /* Calculate gaps for the current block */ /* In later stages, when representations for individual */ @@ -595,7 +599,7 @@ int odrre_(char *range, int *n, double *vl, if (usedqd) { /* The initial SIGMA was to the outer end of the spectrum */ /* the matrix is definite and we need not retreat. */ - tau = spdiam * eps * *n + *pivmin * 2.; + tau = spdiam * eps * (double)(*n) + *pivmin * 2.; tau = fmax(tau, 2 * eps * fabs(sigma)); } else { if (mb > 1) { @@ -668,9 +672,9 @@ int odrre_(char *range, int *n, double *vl, if (idum == 5) { if (sgndef == 1.) { /* The fudged Gerschgorin shift should succeed */ - sigma = gl - spdiam * 2. * eps * *n - *pivmin * 4.; + sigma = gl - spdiam * 2. * eps * (double)(*n) - *pivmin * 4.; } else { - sigma = gu + spdiam * 2. * eps * *n + *pivmin * 4.; + sigma = gu + spdiam * 2. * eps * (double)(*n) + *pivmin * 4.; } } else { sigma -= sgndef * tau; diff --git a/SRC/LAPACK/odrrf.c b/SRC/LAPACK/odrrf.c index b22f9fe..fe98f6f 100644 --- a/SRC/LAPACK/odrrf.c +++ b/SRC/LAPACK/odrrf.c @@ -6,46 +6,47 @@ #include #include #include +#include "global.h" /* Table of constant values */ #define TRUE_ (1) #define FALSE_ (0) -static int c__1 = 1; +static PMRRR_Int c__1 = 1; /* Subroutine */ -int odrrf_(int *n, double *d__, double *l, - double *ld, int *clstrt, int *clend, double *w, +PMRRR_Int odrrf_(PMRRR_Int *n, double *d__, double *l, + double *ld, PMRRR_Int *clstrt, PMRRR_Int *clend, double *w, double *wgap, double *werr, double *spdiam, double * clgapl, double *clgapr, double *pivmin, double *sigma, - double *dplus, double *lplus, double *work, int *info) + double *dplus, double *lplus, double *work, PMRRR_Int *info) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2, d__3; /* Builtin functions */ // double sqrt(double); /* Local variables */ - int i__; + PMRRR_Int i__; double s, bestshift, smlgrowth, eps, tmp, max1, max2, rrr1, rrr2, znm2, growthbound, fail, fact, oldp; - int indx; + PMRRR_Int indx; double prod; - int ktry; + PMRRR_Int ktry; double fail2, avgap, ldmax, rdmax; - int shift; - extern /* Subroutine */ int odcpy_(int *, double *, int *, - double *, int *); - int dorrr1; + PMRRR_Int shift; + extern /* Subroutine */ PMRRR_Int odcpy_(PMRRR_Int *, double *, PMRRR_Int *, + double *, PMRRR_Int *); + PMRRR_Int dorrr1; // extern double odmch_(char *); double ldelta; - int nofail; + PMRRR_Int nofail; double mingap, lsigma, rdelta; - extern int odnan_(double *); - int forcer; + extern PMRRR_Int odnan_(double *); + PMRRR_Int forcer; double rsigma, clwdth; - int sawnan1, sawnan2, tryrrr1; + PMRRR_Int sawnan1, sawnan2, tryrrr1; /* -- LAPACK auxiliary routine (version 3.2) -- */ @@ -177,6 +178,8 @@ int odrrf_(int *n, double *d__, double *l, /* representations despite large element growth or signal INFO=1 */ nofail = TRUE_; + indx = 0; + /* Compute the average gap length of the cluster */ clwdth = (d__1 = w[*clend] - w[*clstrt], fabs(d__1)) + werr[*clend] + werr[ *clstrt]; diff --git a/SRC/LAPACK/odrrj.c b/SRC/LAPACK/odrrj.c index 89647b1..a082974 100644 --- a/SRC/LAPACK/odrrj.c +++ b/SRC/LAPACK/odrrj.c @@ -6,30 +6,31 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odrrj_(int *n, double *d__, double *e2, - int *ifirst, int *ilast, double *rtol, int *offset, - double *w, double *werr, double *work, int *iwork, - double *pivmin, double *spdiam, int *info) +PMRRR_Int odrrj_(PMRRR_Int *n, double *d__, double *e2, + PMRRR_Int *ifirst, PMRRR_Int *ilast, double *rtol, PMRRR_Int *offset, + double *w, double *werr, double *work, PMRRR_Int *iwork, + double *pivmin, double *spdiam, PMRRR_Int *info) { /* System generated locals */ - int i__1, i__2; + PMRRR_Int i__1, i__2; double d__1, d__2; /* Builtin functions */ // double log(double); /* Local variables */ - int i__, j, k, p; + PMRRR_Int i__, j, k, p; double s; - int i1, i2, ii; + PMRRR_Int i1, i2, ii; double fac, mid; - int cnt; + PMRRR_Int cnt; double tmp, left; - int iter, nint, prev, next, savi1; + PMRRR_Int iter, nint, prev, next, savi1; double right, width, dplus; - int olnint, maxitr; + PMRRR_Int olnint, maxitr; /* -- LAPACK auxiliary routine (version 3.2) -- */ @@ -137,12 +138,12 @@ int odrrj_(int *n, double *d__, double *e2, /* Function Body */ *info = 0; - maxitr = (int) ((log(*spdiam + *pivmin) - log(*pivmin)) / log(2.)) + + maxitr = (PMRRR_Int) ((log(*spdiam + *pivmin) - log(*pivmin)) / log(2.)) + 2; /* Initialize unconverged intervals in [ WORK(2*I-1), WORK(2*I) ]. */ /* The Sturm Count, Count( WORK(2*I-1) ) is arranged to be I-1, while */ -/* Count( WORK(2*I) ) is stored in IWORK( 2*I ). The int IWORK( 2*I-1 ) */ +/* Count( WORK(2*I) ) is stored in IWORK( 2*I ). The PMRRR_Int IWORK( 2*I-1 ) */ /* for an unconverged interval is set to the index of the next unconverged */ /* interval, and is -1 or 0 for a converged interval. Thus a linked */ /* list of unconverged intervals is set up. */ diff --git a/SRC/LAPACK/odrrk.c b/SRC/LAPACK/odrrk.c index ac73f85..b58cf5e 100644 --- a/SRC/LAPACK/odrrk.c +++ b/SRC/LAPACK/odrrk.c @@ -6,26 +6,27 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odrrk_(int *n, int *iw, double *gl, +PMRRR_Int odrrk_(PMRRR_Int *n, PMRRR_Int *iw, double *gl, double *gu, double *d__, double *e2, double *pivmin, - double *reltol, double *w, double *werr, int *info) + double *reltol, double *w, double *werr, PMRRR_Int *info) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Builtin functions */ // double log(double); /* Local variables */ - int i__, it; + PMRRR_Int i__, it; double mid, eps, tmp1, tmp2, left, atoli, right; - int itmax; + PMRRR_Int itmax; double rtoli, tnorm; // extern double odmch_(char *); - int negcnt; + PMRRR_Int negcnt; /* -- LAPACK auxiliary routine (version 3.2) -- */ @@ -122,10 +123,10 @@ int odrrk_(int *n, int *iw, double *gl, tnorm = fmax(d__1,d__2); rtoli = *reltol; atoli = *pivmin * 4.; - itmax = (int) ((log(tnorm + *pivmin) - log(*pivmin)) / log(2.)) + 2; + itmax = (PMRRR_Int) ((log(tnorm + *pivmin) - log(*pivmin)) / log(2.)) + 2; *info = -1; - left = *gl - tnorm * 2. * eps * *n - *pivmin * 4.; - right = *gu + tnorm * 2. * eps * *n + *pivmin * 4.; + left = *gl - tnorm * 2. * eps * (double)(*n) - *pivmin * 4.; + right = *gu + tnorm * 2. * eps * (double)(*n) + *pivmin * 4.; it = 0; L10: diff --git a/SRC/LAPACK/odrrr.c b/SRC/LAPACK/odrrr.c index 0384ff4..ab2f9c3 100644 --- a/SRC/LAPACK/odrrr.c +++ b/SRC/LAPACK/odrrr.c @@ -6,27 +6,28 @@ #include #include #include +#include "global.h" #define TRUE_ (1) #define FALSE_ (0) /* Subroutine */ -int odrrr_(int *n, double *d__, double *e, - int *info) +PMRRR_Int odrrr_(PMRRR_Int *n, double *d__, double *e, + PMRRR_Int *info) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1; /* Builtin functions */ // double sqrt(double); /* Local variables */ - int i__; + PMRRR_Int i__; double eps, tmp, tmp2, rmin; // extern double odmch_(char *); double offdig, safmin; - int yesrel; + PMRRR_Int yesrel; double smlnum, offdig2; diff --git a/SRC/LAPACK/odrrv.c b/SRC/LAPACK/odrrv.c index cdbf88f..dcaa756 100644 --- a/SRC/LAPACK/odrrv.c +++ b/SRC/LAPACK/odrrv.c @@ -6,6 +6,7 @@ #include #include #include +#include "global.h" #define imax(a,b) ( (a) > (b) ? (a) : (b) ) #define imin(a,b) ( (a) < (b) ? (a) : (b) ) @@ -14,95 +15,95 @@ /* Table of constant values */ static double c_b5 = 0.; -static int c__1 = 1; -static int c__2 = 2; +static PMRRR_Int c__1 = 1; +static PMRRR_Int c__2 = 2; -/* Subroutine */ int odrrv_(int *n, double *vl, double *vu, - double *d__, double *l, double *pivmin, int *isplit, - int *m, int *dol, int *dou, double *minrgp, +/* Subroutine */ PMRRR_Int odrrv_(PMRRR_Int *n, double *vl, double *vu, + double *d__, double *l, double *pivmin, PMRRR_Int *isplit, + PMRRR_Int *m, PMRRR_Int *dol, PMRRR_Int *dou, double *minrgp, double *rtol1, double *rtol2, double *w, double *werr, - double *wgap, int *iblock, int *indexw, double *gers, - double *z__, int *ldz, int *isuppz, double *work, - int *iwork, int *info) + double *wgap, PMRRR_Int *iblock, PMRRR_Int *indexw, double *gers, + double *z__, PMRRR_Int *ldz, PMRRR_Int *isuppz, double *work, + PMRRR_Int *iwork, PMRRR_Int *info) { /* System generated locals */ - int z_dim1, z_offset, i__1, i__2, i__3, i__4, i__5; + PMRRR_Int z_dim1, z_offset, i__1, i__2, i__3, i__4, i__5; double d__1, d__2; - int L__1; + PMRRR_Int L__1; /* Builtin functions */ // double log(double); /* Local variables */ - int minwsize, i__, j, k, p, q, miniwsize, ii; + PMRRR_Int minwsize, i__, j, k, p, q, miniwsize, ii; double gl; - int im, in; + PMRRR_Int im, in; double gu, gap, eps, tau, tol, tmp; - int zto; + PMRRR_Int zto; double ztz; - int iend, jblk; + PMRRR_Int iend, jblk; double lgap; - int done; + PMRRR_Int done; double rgap, left; - int wend, iter; + PMRRR_Int wend, iter; double bstw; - int itmp1; - extern /* Subroutine */ int odscl_(int *, double *, double *, - int *); - int indld; + PMRRR_Int itmp1; + extern /* Subroutine */ PMRRR_Int odscl_(PMRRR_Int *, double *, double *, + PMRRR_Int *); + PMRRR_Int indld; double fudge; - int idone; + PMRRR_Int idone; double sigma; - int iinfo, iindr; + PMRRR_Int iinfo, iindr; double resid; - int eskip; + PMRRR_Int eskip; double right; - extern /* Subroutine */ int odcpy_(int *, double *, int *, - double *, int *); - int nclus, zfrom; + extern /* Subroutine */ PMRRR_Int odcpy_(PMRRR_Int *, double *, PMRRR_Int *, + double *, PMRRR_Int *); + PMRRR_Int nclus, zfrom; double rqtol; - int iindc1, iindc2; - extern /* Subroutine */ int odr1v_(int *, int *, int *, + PMRRR_Int iindc1, iindc2; + extern /* Subroutine */ PMRRR_Int odr1v_(PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, double *, double *, double *, double *, - double *, double *, double *, double *, int *, - int *, double *, double *, int *, int *, + double *, double *, double *, double *, PMRRR_Int *, + PMRRR_Int *, double *, double *, PMRRR_Int *, PMRRR_Int *, double *, double *, double *, double *); - int stp2ii; + PMRRR_Int stp2ii; double lambda; // extern double odmch_(char *); - int ibegin, indeig; - int needbs; - int indlld; + PMRRR_Int ibegin, indeig; + PMRRR_Int needbs; + PMRRR_Int indlld; double sgndef, mingma; - extern /* Subroutine */ int odrrb_(int *, double *, double *, - int *, int *, double *, double *, int *, - double *, double *, double *, double *, int *, - double *, double *, int *, int *); - int oldien, oldncl, wbegin; + extern /* Subroutine */ PMRRR_Int odrrb_(PMRRR_Int *, double *, double *, + PMRRR_Int *, PMRRR_Int *, double *, double *, PMRRR_Int *, + double *, double *, double *, double *, PMRRR_Int *, + double *, double *, PMRRR_Int *, PMRRR_Int *); + PMRRR_Int oldien, oldncl, wbegin; double spdiam; - int negcnt; - extern /* Subroutine */ int odrrf_(int *, double *, double *, - double *, int *, int *, double *, double *, + PMRRR_Int negcnt; + extern /* Subroutine */ PMRRR_Int odrrf_(PMRRR_Int *, double *, double *, + double *, PMRRR_Int *, PMRRR_Int *, double *, double *, double *, double *, double *, double *, double *, double *, double *, double *, - double *, int *); - int oldcls; + double *, PMRRR_Int *); + PMRRR_Int oldcls; double savgap; - int ndepth; + PMRRR_Int ndepth; double ssigma; - extern /* Subroutine */ int odset_(char *, int *, int *, - double *, double *, double *, int *); - int usedbs; - int iindwk, offset; + extern /* Subroutine */ PMRRR_Int odset_(char *, PMRRR_Int *, PMRRR_Int *, + double *, double *, double *, PMRRR_Int *); + PMRRR_Int usedbs; + PMRRR_Int iindwk, offset; double gaptol; - int newcls, oldfst, indwrk, windex, oldlst; - int usedrq; - int newfst, newftt, parity, windmn, windpl, isupmn, newlst, zusedl; + PMRRR_Int newcls, oldfst, indwrk, windex, oldlst; + PMRRR_Int usedrq; + PMRRR_Int newfst, newftt, parity, windmn, windpl, isupmn, newlst, zusedl; double bstres; - int newsiz, zusedu, zusedw; + PMRRR_Int newsiz, zusedu, zusedw; double nrminv, rqcorr; - int tryrqc; - int isupmx; + PMRRR_Int tryrqc; + PMRRR_Int isupmx; /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odruv.c b/SRC/LAPACK/odruv.c index b1e231c..1535155 100644 --- a/SRC/LAPACK/odruv.c +++ b/SRC/LAPACK/odruv.c @@ -6,16 +6,17 @@ #include #include #include +#include "global.h" #define imax(a,b) ( (a) > (b) ? (a) : (b) ) #define imin(a,b) ( (a) < (b) ? (a) : (b) ) /* Subroutine */ -int odruv_(int *iseed, int *n, double *x) +PMRRR_Int odruv_(PMRRR_Int *iseed, PMRRR_Int *n, double *x) { /* Initialized data */ - static int mm[512] /* was [128][4] */ = { 494,2637,255,2008,1253, + static PMRRR_Int mm[512] /* was [128][4] */ = { 494,2637,255,2008,1253, 3344,4084,1739,3143,3468,688,1657,1238,3166,1292,3422,1270,2016, 154,2862,697,1706,491,931,1444,444,3577,3944,2184,1661,3482,657, 3023,3618,1267,1828,164,3798,3087,2400,2870,3876,1905,1593,1797, @@ -56,10 +57,10 @@ int odruv_(int *iseed, int *n, double *x) 3537,517,3017,2141,1537 }; /* System generated locals */ - int i__1; + PMRRR_Int i__1; /* Local variables */ - int i__, i1, i2, i3, i4, it1, it2, it3, it4; + PMRRR_Int i__, i1, i2, i3, i4, it1, it2, it3, it4; /* -- LAPACK auxiliary routine (version 3.2) -- */ @@ -103,7 +104,7 @@ int odruv_(int *iseed, int *n, double *x) /* 2**b: an exhaustive analysis for b = 32 and a partial analysis for */ /* b = 48', Math. Comp. 189, pp 331-344, 1990). */ -/* 48-bit ints are stored in 4 int array elements with 12 bits */ +/* 48-bit ints are stored in 4 PMRRR_Int array elements with 12 bits */ /* per element. Hence the routine is portable across machines with */ /* ints of 32 bits or more. */ @@ -131,6 +132,8 @@ int odruv_(int *iseed, int *n, double *x) i3 = iseed[3]; i4 = iseed[4]; + it1 = it2 = it3 = it4 = 0; + i__1 = imin(*n,128); for (i__ = 1; i__ <= i__1; ++i__) { @@ -152,7 +155,7 @@ int odruv_(int *iseed, int *n, double *x) 127] + i4 * mm[i__ - 1]; it1 %= 4096; -/* Convert 48-bit int to a real number in the interval (0,1) */ +/* Convert 48-bit PMRRR_Int to a real number in the interval (0,1) */ x[i__] = ((double) it1 + ((double) it2 + ((double) it3 + ( double) it4 * 2.44140625e-4) * 2.44140625e-4) * @@ -160,7 +163,7 @@ int odruv_(int *iseed, int *n, double *x) if (x[i__] == 1.) { /* If a real number has n bits of precision, and the first */ -/* n bits of the 48-bit int above happen to be all 1 (which */ +/* n bits of the 48-bit PMRRR_Int above happen to be all 1 (which */ /* will occur about once every 2**n calls), then X( I ) will */ /* be rounded to exactly 1.0. */ /* Since X( I ) is not supposed to return exactly 0.0 or 1.0, */ diff --git a/SRC/LAPACK/odset.c b/SRC/LAPACK/odset.c index 62cbacf..79011e7 100644 --- a/SRC/LAPACK/odset.c +++ b/SRC/LAPACK/odset.c @@ -7,18 +7,19 @@ #include #include #include +#include "global.h" #define imin(a,b) ( (a) < (b) ? (a) : (b) ) -/* Subroutine */ int odset_(char *uplo, int *m, int *n, double *alpha, -double *beta, double *a, int *lda) +/* Subroutine */ PMRRR_Int odset_(char *uplo, PMRRR_Int *m, PMRRR_Int *n, double *alpha, +double *beta, double *a, PMRRR_Int *lda) { /* System generated locals */ - int a_dim1, a_offset, i__1, i__2, i__3; + PMRRR_Int a_dim1, a_offset, i__1, i__2, i__3; /* Local variables */ - int i__, j; - extern int olsame_(char *, char *); + PMRRR_Int i__, j; + extern PMRRR_Int olsame_(char *, char *); /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odsnan.c b/SRC/LAPACK/odsnan.c index 530e705..8afc8ce 100644 --- a/SRC/LAPACK/odsnan.c +++ b/SRC/LAPACK/odsnan.c @@ -6,11 +6,12 @@ #include #include #include +#include "global.h" -int odsnan_(double *din1, double *din2) +PMRRR_Int odsnan_(double *din1, double *din2) { /* System generated locals */ - int ret_val; + PMRRR_Int ret_val; /* -- LAPACK auxiliary routine (version 3.2) -- */ /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ diff --git a/SRC/LAPACK/odsq2.c b/SRC/LAPACK/odsq2.c index 88b0116..81d7491 100644 --- a/SRC/LAPACK/odsq2.c +++ b/SRC/LAPACK/odsq2.c @@ -6,20 +6,17 @@ #include #include #include +#include "global.h" /* Table of constant values */ -static int c__1 = 1; -static int c__2 = 2; -static int c__10 = 10; -static int c__3 = 3; -static int c__4 = 4; -static int c__11 = 11; +static PMRRR_Int c__1 = 1; +static PMRRR_Int c__2 = 2; /* Subroutine */ -int odsq2_(int *n, double *z__, int *info) +PMRRR_Int odsq2_(PMRRR_Int *n, double *z__, PMRRR_Int *info) { /* System generated locals */ - int i__1, i__2, i__3; + PMRRR_Int i__1, i__2, i__3; double d__1, d__2; /* Builtin functions */ @@ -27,38 +24,38 @@ int odsq2_(int *n, double *z__, int *info) /* Local variables */ double d__, e, g; - int k; + PMRRR_Int k; double s, t; - int i0, i4, n0; + PMRRR_Int i0, i4, n0; double dn; - int pp; + PMRRR_Int pp; double dn1, dn2, dee, eps, tau, tol; - int ipn4; + PMRRR_Int ipn4; double tol2; - int ieee; - int nbig; + PMRRR_Int ieee; + PMRRR_Int nbig; double dmin__, emin, emax; - int kmin, ndiv, iter; + PMRRR_Int kmin, ndiv, iter; double qmin, temp, qmax, zmax; - int splt; + PMRRR_Int splt; double dmin1, dmin2; - int nfail; + PMRRR_Int nfail; double desig, trace, sigma; - int iinfo, ttype; - extern /* Subroutine */ int odsq3_(int *, int *, double *, - int *, double *, double *, double *, double *, - int *, int *, int *, int *, int *, + PMRRR_Int iinfo, ttype; + extern /* Subroutine */ PMRRR_Int odsq3_(PMRRR_Int *, PMRRR_Int *, double *, + PMRRR_Int *, double *, double *, double *, double *, + PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, double *, double *, double *, double *, double *, double *, double *); // extern double odmch_(char *); double deemin; - int iwhila, iwhilb; + PMRRR_Int iwhila, iwhilb; double oldemn, safmin; - extern /* Subroutine */ int oerbla_(char *, int *); - /* extern int oienv_(int *, char *, char *, int *, int *, */ - /* int *, int *); */ - extern /* Subroutine */ int odsrt_(char *, int *, double *, - int *); + extern /* Subroutine */ PMRRR_Int oerbla_(char *, PMRRR_Int *); + /* extern PMRRR_Int oienv_(PMRRR_Int *, char *, char *, PMRRR_Int *, PMRRR_Int *, */ + /* PMRRR_Int *, PMRRR_Int *); */ + extern /* Subroutine */ PMRRR_Int odsrt_(char *, PMRRR_Int *, double *, + PMRRR_Int *); /* -- LAPACK routine (version 3.2) -- */ @@ -90,7 +87,7 @@ int odsq2_(int *n, double *z__, int *info) /* Z(1,3,5,,..). The tridiagonal is L*U or, if you prefer, the */ /* symmetric tridiagonal to which it is similar. */ -/* Note : ODSQ2 defines a int variable, IEEE, which is true */ +/* Note : ODSQ2 defines a PMRRR_Int variable, IEEE, which is true */ /* on machines which follow ieee-754 floating-point standard in their */ /* handling of infinities and NaNs, and false otherwise. This variable */ /* is passed to ODSQ3. */ @@ -585,7 +582,7 @@ int odsq2_(int *n, double *z__, int *info) /* Computing 2nd power */ i__1 = *n; z__[(*n << 1) + 4] = (double) ndiv / (double) (i__1 * i__1); - z__[(*n << 1) + 5] = nfail * 100. / (double) iter; + z__[(*n << 1) + 5] = (double) nfail * 100. / (double) iter; return 0; /* End of ODSQ2 */ diff --git a/SRC/LAPACK/odsq3.c b/SRC/LAPACK/odsq3.c index 2332f3b..31903a7 100644 --- a/SRC/LAPACK/odsq3.c +++ b/SRC/LAPACK/odsq3.c @@ -6,16 +6,17 @@ #include #include #include +#include "global.h" -/* Subroutine */ int odsq3_(int *i0, int *n0, double *z__, - int *pp, double *dmin__, double *sigma, double *desig, - double *qmax, int *nfail, int *iter, int *ndiv, - int *ieee, int *ttype, double *dmin1, double *dmin2, +/* Subroutine */ PMRRR_Int odsq3_(PMRRR_Int *i0, PMRRR_Int *n0, double *z__, + PMRRR_Int *pp, double *dmin__, double *sigma, double *desig, + double *qmax, PMRRR_Int *nfail, PMRRR_Int *iter, PMRRR_Int *ndiv, + PMRRR_Int *ieee, PMRRR_Int *ttype, double *dmin1, double *dmin2, double *dn, double *dn1, double *dn2, double *g, double *tau) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Builtin functions */ @@ -23,21 +24,21 @@ /* Local variables */ double s, t; - int j4, nn; + PMRRR_Int j4, nn; double eps, tol; - int n0in, ipn4; + PMRRR_Int n0in, ipn4; double tol2, temp; - extern /* Subroutine */ int odsq4_(int *, int *, double *, - int *, int *, double *, double *, double *, - double *, double *, double *, double *, int *, - double *), odsq5_(int *, int *, double *, - int *, double *, double *, double *, double *, - double *, double *, double *, int *), odsq6_( - int *, int *, double *, int *, double *, + extern /* Subroutine */ PMRRR_Int odsq4_(PMRRR_Int *, PMRRR_Int *, double *, + PMRRR_Int *, PMRRR_Int *, double *, double *, double *, + double *, double *, double *, double *, PMRRR_Int *, + double *), odsq5_(PMRRR_Int *, PMRRR_Int *, double *, + PMRRR_Int *, double *, double *, double *, double *, + double *, double *, double *, PMRRR_Int *), odsq6_( + PMRRR_Int *, PMRRR_Int *, double *, PMRRR_Int *, double *, double *, double *, double *, double *, double *); // extern double odmch_(char *); - extern int odnan_(double *); + extern PMRRR_Int odnan_(double *); /* -- LAPACK routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odsq4.c b/SRC/LAPACK/odsq4.c index 390d92b..678a1da 100644 --- a/SRC/LAPACK/odsq4.c +++ b/SRC/LAPACK/odsq4.c @@ -6,15 +6,16 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odsq4_(int *i0, int *n0, double *z__, - int *pp, int *n0in, double *dmin__, double *dmin1, +PMRRR_Int odsq4_(PMRRR_Int *i0, PMRRR_Int *n0, double *z__, + PMRRR_Int *pp, PMRRR_Int *n0in, double *dmin__, double *dmin1, double *dmin2, double *dn, double *dn1, double *dn2, - double *tau, int *ttype, double *g) + double *tau, PMRRR_Int *ttype, double *g) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Builtin functions */ @@ -22,7 +23,7 @@ int odsq4_(int *i0, int *n0, double *z__, /* Local variables */ double s, a2, b1, b2; - int i4, nn, np; + PMRRR_Int i4, nn, np; double gam, gap1, gap2; @@ -110,6 +111,8 @@ int odsq4_(int *i0, int *n0, double *z__, /* Parameter adjustments */ --z__; + s = 0.0; + /* Function Body */ if (*dmin__ <= 0.) { *tau = -(*dmin__); diff --git a/SRC/LAPACK/odsq5.c b/SRC/LAPACK/odsq5.c index ed15030..2de2bfa 100644 --- a/SRC/LAPACK/odsq5.c +++ b/SRC/LAPACK/odsq5.c @@ -6,20 +6,21 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odsq5_(int *i0, int *n0, double *z__, - int *pp, double *tau, double *dmin__, double *dmin1, +PMRRR_Int odsq5_(PMRRR_Int *i0, PMRRR_Int *n0, double *z__, + PMRRR_Int *pp, double *tau, double *dmin__, double *dmin1, double *dmin2, double *dn, double *dnm1, double *dnm2, - int *ieee) + PMRRR_Int *ieee) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Local variables */ double d__; - int j4, j4p2; + PMRRR_Int j4, j4p2; double emin, temp; diff --git a/SRC/LAPACK/odsq6.c b/SRC/LAPACK/odsq6.c index 969c65a..f8d75eb 100644 --- a/SRC/LAPACK/odsq6.c +++ b/SRC/LAPACK/odsq6.c @@ -6,19 +6,20 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odsq6_(int *i0, int *n0, double *z__, - int *pp, double *dmin__, double *dmin1, double *dmin2, +PMRRR_Int odsq6_(PMRRR_Int *i0, PMRRR_Int *n0, double *z__, + PMRRR_Int *pp, double *dmin__, double *dmin1, double *dmin2, double *dn, double *dnm1, double *dnm2) { /* System generated locals */ - int i__1; + PMRRR_Int i__1; double d__1, d__2; /* Local variables */ double d__; - int j4, j4p2; + PMRRR_Int j4, j4p2; double emin, temp; // extern double odmch_(char *); double safmin; diff --git a/SRC/LAPACK/odsrt.c b/SRC/LAPACK/odsrt.c index 7111cde..259dfb5 100644 --- a/SRC/LAPACK/odsrt.c +++ b/SRC/LAPACK/odsrt.c @@ -6,26 +6,27 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odsrt_(char *id, int *n, double *d__, int * +PMRRR_Int odsrt_(char *id, PMRRR_Int *n, double *d__, PMRRR_Int * info) { /* System generated locals */ - int i__1, i__2; + PMRRR_Int i__1, i__2; /* Local variables */ - int i__, j; + PMRRR_Int i__, j; double d1, d2, d3; - int dir; + PMRRR_Int dir; double tmp; - int endd; - extern int olsame_(char *, char *); - int stack[64] /* was [2][32] */; + PMRRR_Int endd; + extern PMRRR_Int olsame_(char *, char *); + PMRRR_Int stack[64] /* was [2][32] */; double dmnmx; - int start; - extern /* Subroutine */ int oerbla_(char *, int *); - int stkpnt; + PMRRR_Int start; + extern /* Subroutine */ PMRRR_Int oerbla_(char *, PMRRR_Int *); + PMRRR_Int stkpnt; /* -- LAPACK routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odssq.c b/SRC/LAPACK/odssq.c index 54302c5..248bac1 100644 --- a/SRC/LAPACK/odssq.c +++ b/SRC/LAPACK/odssq.c @@ -6,17 +6,18 @@ #include #include #include +#include "global.h" /* Subroutine */ -int odssq_(int *n, double *x, int *incx, +PMRRR_Int odssq_(PMRRR_Int *n, double *x, PMRRR_Int *incx, double *scale, double *sumsq) { /* System generated locals */ - int i__1, i__2; + PMRRR_Int i__1, i__2; double d__1; /* Local variables */ - int ix; + PMRRR_Int ix; double absxi; /* -- LAPACK auxiliary routine (version 3.2) -- */ diff --git a/SRC/LAPACK/odstmr.c b/SRC/LAPACK/odstmr.c index fa90440..6fd539c 100644 --- a/SRC/LAPACK/odstmr.c +++ b/SRC/LAPACK/odstmr.c @@ -6,9 +6,10 @@ #include #include #include +#include "global.h" /* Table of constant values */ -static int c__1 = 1; +static PMRRR_Int c__1 = 1; static double c_b18 = .001; #define TRUE_ (1) #define FALSE_ (0) @@ -16,87 +17,87 @@ static double c_b18 = .001; /* Subroutine */ -int odstmr_(char *jobz, char *range, int *n, double *d__, - double *e, double *vl, double *vu, int *il, - int *iu, int *m, double *w, double *z__, int *ldz, - int *nzc, int *isuppz, int *tryrac, double *work, - int *lwork, int *iwork, int *liwork, int *info) +PMRRR_Int odstmr_(char *jobz, char *range, PMRRR_Int *n, double *d__, + double *e, double *vl, double *vu, PMRRR_Int *il, + PMRRR_Int *iu, PMRRR_Int *m, double *w, double *z__, PMRRR_Int *ldz, + PMRRR_Int *nzc, PMRRR_Int *isuppz, PMRRR_Int *tryrac, double *work, + PMRRR_Int *lwork, PMRRR_Int *iwork, PMRRR_Int *liwork, PMRRR_Int *info) { /* System generated locals */ - int z_dim1, z_offset, i__1, i__2; + PMRRR_Int z_dim1, z_offset, i__1, i__2; double d__1, d__2; /* Builtin functions */ // double sqrt(double); /* Local variables */ - int i__, j; + PMRRR_Int i__, j; double r1, r2; - int jj; + PMRRR_Int jj; double cs; - int in; + PMRRR_Int in; double sn, wl, wu; - int iil, iiu; + PMRRR_Int iil, iiu; double eps, tmp; - int indd, iend, jblk, wend; + PMRRR_Int indd, iend, jblk, wend; double rmin, rmax; - int itmp; + PMRRR_Int itmp; double tnrm; - extern /* Subroutine */ int ode2_(double *, double *, double + extern /* Subroutine */ PMRRR_Int ode2_(double *, double *, double *, double *, double *); - int inde2, itmp2; + PMRRR_Int inde2, itmp2; double rtol1, rtol2; - extern /* Subroutine */ int odscl_(int *, double *, double *, - int *); + extern /* Subroutine */ PMRRR_Int odscl_(PMRRR_Int *, double *, double *, + PMRRR_Int *); double scale; - int indgp; - extern int olsame_(char *, char *); - int iinfo, iindw, ilast; - extern /* Subroutine */ int odcpy_(int *, double *, int *, - double *, int *), odswap_(int *, double *, int - *, double *, int *); - int lwmin; - int wantz; - extern /* Subroutine */ int odev2_(double *, double *, + PMRRR_Int indgp; + extern PMRRR_Int olsame_(char *, char *); + PMRRR_Int iinfo, iindw, ilast; + extern /* Subroutine */ PMRRR_Int odcpy_(PMRRR_Int *, double *, PMRRR_Int *, + double *, PMRRR_Int *), odswap_(PMRRR_Int *, double *, PMRRR_Int + *, double *, PMRRR_Int *); + PMRRR_Int lwmin; + PMRRR_Int wantz; + extern /* Subroutine */ PMRRR_Int odev2_(double *, double *, double *, double *, double *, double *, double *); // extern double odmch_(char *); - int alleig; - int ibegin; - int indeig; - int iindbl; - int valeig; - extern /* Subroutine */ int odrrc_(char *, int *, double *, - double *, double *, double *, double *, int *, - int *, int *, int *), odrre_(char *, - int *, double *, double *, int *, int *, + PMRRR_Int alleig; + PMRRR_Int ibegin; + PMRRR_Int indeig; + PMRRR_Int iindbl; + PMRRR_Int valeig; + extern /* Subroutine */ PMRRR_Int odrrc_(char *, PMRRR_Int *, double *, + double *, double *, double *, double *, PMRRR_Int *, + PMRRR_Int *, PMRRR_Int *, PMRRR_Int *), odrre_(char *, + PMRRR_Int *, double *, double *, PMRRR_Int *, PMRRR_Int *, double *, double *, double *, double *, - double *, double *, int *, int *, int *, - double *, double *, double *, int *, int *, - double *, double *, double *, int *, int *); - int wbegin; + double *, double *, PMRRR_Int *, PMRRR_Int *, PMRRR_Int *, + double *, double *, double *, PMRRR_Int *, PMRRR_Int *, + double *, double *, double *, PMRRR_Int *, PMRRR_Int *); + PMRRR_Int wbegin; double safmin; - extern /* Subroutine */ int odrrj_(int *, double *, double *, - int *, int *, double *, int *, double *, - double *, double *, int *, double *, double *, - int *), oerbla_(char *, int *); + extern /* Subroutine */ PMRRR_Int odrrj_(PMRRR_Int *, double *, double *, + PMRRR_Int *, PMRRR_Int *, double *, PMRRR_Int *, double *, + double *, double *, PMRRR_Int *, double *, double *, + PMRRR_Int *), oerbla_(char *, PMRRR_Int *); double bignum; - int inderr, iindwk, indgrs, offset; - extern double odnst_(char *, int *, double *, double *); - extern /* Subroutine */ int odrrr_(int *, double *, double *, - int *), odrrv_(int *, double *, double *, - double *, double *, double *, int *, int *, - int *, int *, double *, double *, double *, - double *, double *, double *, int *, int *, - double *, double *, int *, int *, double *, - int *, int *), odsrt_(char *, int *, double *, - int *); + PMRRR_Int inderr, iindwk, indgrs, offset; + extern double odnst_(char *, PMRRR_Int *, double *, double *); + extern /* Subroutine */ PMRRR_Int odrrr_(PMRRR_Int *, double *, double *, + PMRRR_Int *), odrrv_(PMRRR_Int *, double *, double *, + double *, double *, double *, PMRRR_Int *, PMRRR_Int *, + PMRRR_Int *, PMRRR_Int *, double *, double *, double *, + double *, double *, double *, PMRRR_Int *, PMRRR_Int *, + double *, double *, PMRRR_Int *, PMRRR_Int *, double *, + PMRRR_Int *, PMRRR_Int *), odsrt_(char *, PMRRR_Int *, double *, + PMRRR_Int *); double thresh; - int iinspl, ifirst, indwrk, liwmin, nzcmin; + PMRRR_Int iinspl, ifirst, indwrk, liwmin, nzcmin; double pivmin; - int nsplit; + PMRRR_Int nsplit; double smlnum; - int lquery, zquery; + PMRRR_Int lquery, zquery; /* -- LAPACK computational routine (version 3.2) -- */ diff --git a/SRC/LAPACK/oerbla.c b/SRC/LAPACK/oerbla.c index 744a637..38ad3c8 100644 --- a/SRC/LAPACK/oerbla.c +++ b/SRC/LAPACK/oerbla.c @@ -6,12 +6,11 @@ #include #include #include +#include "global.h" -/* Table of constant values */ -static int c__1 = 1; /* Subroutine */ -int oerbla_(char *srname, int *info) +PMRRR_Int oerbla_(char *srname, PMRRR_Int *info) { @@ -49,7 +48,7 @@ int oerbla_(char *srname, int *info) /* .. Executable Statements .. */ printf("** On entry to %6s, parameter number %2i had an illegal value\n", - srname, *info); + srname, (int)(*info)); /* End of OERBLA */ diff --git a/SRC/LAPACK/olsame.c b/SRC/LAPACK/olsame.c index da8e30f..640acd6 100644 --- a/SRC/LAPACK/olsame.c +++ b/SRC/LAPACK/olsame.c @@ -6,16 +6,17 @@ #include #include #include +#include "global.h" -int olsame_(char *ca, char *cb) +PMRRR_Int olsame_(char *ca, char *cb) { /* System generated locals */ - int ret_val; + PMRRR_Int ret_val; /* Local variables */ - static int inta, intb, zcode; + static PMRRR_Int inta, intb, zcode; /* -- LAPACK auxiliary routine (version 2.0) -- diff --git a/SRC/counter.c b/SRC/counter.c index 230c27d..54e4255 100644 --- a/SRC/counter.c +++ b/SRC/counter.c @@ -46,9 +46,9 @@ -counter_t *PMR_create_counter(int init_value) +counter_t *PMR_create_counter(PMRRR_Int init_value) { - int info; + PMRRR_Int info; counter_t *counter; counter = (counter_t *) malloc( sizeof(counter_t) ); @@ -78,9 +78,9 @@ void PMR_destroy_counter(counter_t *counter) -int PMR_get_counter_value(counter_t *counter) +PMRRR_Int PMR_get_counter_value(counter_t *counter) { - int value, info; + PMRRR_Int value, info; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&counter->lock); @@ -103,9 +103,9 @@ int PMR_get_counter_value(counter_t *counter) -inline int PMR_set_counter_value(counter_t *counter, int value) +inline PMRRR_Int PMR_set_counter_value(counter_t *counter, PMRRR_Int value) { - int info; + PMRRR_Int info; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&counter->lock); @@ -128,9 +128,9 @@ inline int PMR_set_counter_value(counter_t *counter, int value) -int PMR_decrement_counter(counter_t *counter, int amount) +PMRRR_Int PMR_decrement_counter(counter_t *counter, PMRRR_Int amount) { - int value, info; + PMRRR_Int value, info; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&counter->lock); @@ -155,9 +155,9 @@ int PMR_decrement_counter(counter_t *counter, int amount) -int PMR_increment_counter(counter_t *counter, int amount) +PMRRR_Int PMR_increment_counter(counter_t *counter, PMRRR_Int amount) { - int value, info; + PMRRR_Int value, info; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&counter->lock); diff --git a/SRC/odscal.c b/SRC/odscal.c index 7a5873a..e1eb9c9 100644 --- a/SRC/odscal.c +++ b/SRC/odscal.c @@ -2,11 +2,11 @@ #include "global.h" /* non-optimized, non-threaded DSCAL replacement */ -void pmrrr_dscal(int *n, double *alpha, double *restrict x, int *incx) +void pmrrr_dscal(PMRRR_Int *n, double *alpha, double *restrict x, PMRRR_Int *incx) { - int i; - int stride = *incx; - int size = *n; + PMRRR_Int i; + PMRRR_Int stride = *incx; + PMRRR_Int size = *n; double s = *alpha; for (i=0; ipid; - int nproc = procinfo->nproc; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int nproc = procinfo->nproc; bool wantZ = (jobz[0] == 'V' || jobz[0] == 'v'); bool cntval = (jobz[0] == 'C' || jobz[0] == 'c'); - int n = Dstruct->n; + PMRRR_Int n = Dstruct->n; double *restrict D = Dstruct->D; double *restrict E = Dstruct->E; - int *restrict isplit = Dstruct->isplit; + PMRRR_Int *restrict isplit = Dstruct->isplit; double *vl = Wstruct->vl; double *vu = Wstruct->vu; - int *il = Wstruct->il; - int *iu = Wstruct->iu; + PMRRR_Int *il = Wstruct->il; + PMRRR_Int *iu = Wstruct->iu; double *restrict W = Wstruct->W; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; - int *restrict Windex = Wstruct->Windex; - int *restrict iblock = Wstruct->iblock; + PMRRR_Int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict iblock = Wstruct->iblock; double *restrict gersch = Wstruct->gersch; /* constants */ - int IZERO = 0, IONE = 1; + PMRRR_Int IZERO = 0, IONE = 1; double DZERO = 0.0; /* work space */ double *E2; double *work; - int *iwork; + PMRRR_Int *iwork; /* compute geschgorin disks and spectral diameter */ double gl, gu, eold, emax, eabs; /* compute splitting points */ - int bl_begin, bl_end, bl_size; + PMRRR_Int bl_begin, bl_end, bl_size; /* distribute work among processes */ - int ifirst, ilast, ifirst_tmp, ilast_tmp; - int chunk, isize, iil, iiu; + PMRRR_Int ifirst, ilast, ifirst_tmp, ilast_tmp; + PMRRR_Int chunk, isize, iil, iiu; /* gather results */ int *rcount, *rdispl; /* others */ - int info, i, j, jbl, idummy; + PMRRR_Int info, i, j, jbl, idummy; double tmp1, dummy; bool sorted; enum range_enum {allrng=1, valrng=2, indrng=3} irange; double intervals[2]; - int negcounts[2]; + PMRRR_Int negcounts[2]; double sigma; if (range[0] == 'A' || range[0] == 'a') { @@ -176,15 +174,15 @@ int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, } /* allocate work space */ - E2 = (double *) malloc( n * sizeof(double) ); + E2 = (double *) malloc( (size_t)n * sizeof(double) ); assert(E2 != NULL); - work = (double *) malloc( 4*n * sizeof(double) ); + work = (double *) malloc( 4*(size_t)n * sizeof(double) ); assert(work != NULL); - iwork = (int *) malloc( 3*n * sizeof(int) ); + iwork = (PMRRR_Int *) malloc( 3*(size_t)n * sizeof(PMRRR_Int) ); assert(iwork != NULL); - rcount = (int *) malloc( nproc * sizeof(int) ); + rcount = (int *) malloc( (size_t)nproc * sizeof(int) ); assert(rcount != NULL); - rdispl = (int *) malloc( nproc * sizeof(int) ); + rdispl = (int *) malloc( (size_t)nproc * sizeof(int) ); assert(rdispl != NULL); /* Compute square of off-diagonal elements */ @@ -255,7 +253,8 @@ int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, /* loop over unreduced blocks */ - bl_begin = 0; + bl_begin = 0; + isize = ifirst = ilast = 0; for (jbl=0; jblnsplit; jbl++) { @@ -295,8 +294,8 @@ int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, *offsetp = ifirst - iil; *nzp = isize; } - rcount[i] = ilast_tmp - ifirst_tmp + 1; - rdispl[i] = ifirst_tmp - iil; + rcount[i] = (int)(ilast_tmp - ifirst_tmp + 1); + rdispl[i] = (int)(ifirst_tmp - iil); ifirst_tmp = ilast_tmp + 1; ifirst_tmp = imin(ifirst_tmp, iiu + 1); } @@ -329,17 +328,17 @@ int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, assert(info == 0); } - memcpy(work, &W[bl_begin], isize * sizeof(double) ); - MPI_Allgatherv(work, isize, MPI_DOUBLE, &W[bl_begin], rcount, rdispl, + memcpy(work, &W[bl_begin], (size_t)isize * sizeof(double) ); + MPI_Allgatherv(work, (int)isize, MPI_DOUBLE, &W[bl_begin], rcount, rdispl, MPI_DOUBLE, procinfo->comm); - memcpy(work, &Werr[bl_begin], isize * sizeof(double) ); - MPI_Allgatherv(work, isize, MPI_DOUBLE, &Werr[bl_begin], rcount, rdispl, + memcpy(work, &Werr[bl_begin], (size_t)isize * sizeof(double) ); + MPI_Allgatherv(work, (int)isize, MPI_DOUBLE, &Werr[bl_begin], rcount, rdispl, MPI_DOUBLE, procinfo->comm); - memcpy(iwork, &Windex[bl_begin], isize * sizeof(int) ); - MPI_Allgatherv(iwork, isize, MPI_INT, &Windex[bl_begin], rcount, rdispl, - MPI_INT, procinfo->comm); + memcpy(iwork, &Windex[bl_begin], (size_t)isize * sizeof(PMRRR_Int) ); + MPI_Allgatherv(iwork, (int)isize, PMRRR_MPI_INT_TYPE, &Windex[bl_begin], rcount, rdispl, + PMRRR_MPI_INT_TYPE, procinfo->comm); /* Ensure that within block eigenvalues sorted */ sorted = false; @@ -395,7 +394,7 @@ int plarre(proc_t *procinfo, char *jobz, char *range, in_t *Dstruct, * Free's on allocated memory of plarre routine */ static -void clean_up_plarre(double *E2, double *work, int *iwork, +void clean_up_plarre(double *E2, double *work, PMRRR_Int *iwork, int *rcount, int *rdispl) { free(E2); @@ -409,16 +408,15 @@ void clean_up_plarre(double *E2, double *work, int *iwork, static -int eigval_approx_proc(proc_t *procinfo, int ifirst, int ilast, - int n, double *D, double *E, double *E2, - int *Windex, int *iblock, double *gersch, tol_t *tolstruct, +PMRRR_Int eigval_approx_proc(proc_t *procinfo, PMRRR_Int ifirst, PMRRR_Int ilast, + PMRRR_Int n, double *D, double *E, double *E2, + PMRRR_Int *Windex, PMRRR_Int *iblock, double *gersch, tol_t *tolstruct, double *W, double *Werr, double *Wgap, double *work, - int *iwork) + PMRRR_Int *iwork) { /* Input parameter */ - int pid = procinfo->pid; - int isize = ilast-ifirst+1; - double pivmin = tolstruct->pivmin; + PMRRR_Int isize = ilast-ifirst+1; + double pivmin = tolstruct->pivmin; /* double gl, gu, wl, wu; */ double wl, wu; @@ -427,23 +425,23 @@ int eigval_approx_proc(proc_t *procinfo, int ifirst, int ilast, double bsrtol; /* /\* Multithreading *\/ */ - int nthreads; - int max_nthreads = procinfo->nthreads; - int iifirst, iilast, chunk; + PMRRR_Int nthreads; + PMRRR_Int max_nthreads = procinfo->nthreads; + PMRRR_Int iifirst, iilast, chunk; pthread_t *threads; pthread_attr_t attr; auxarg1_t *auxarg1; void *status; /* Others */ - int nsplit, *isplit; - int info, m, i, j; + PMRRR_Int nsplit, *isplit; + PMRRR_Int info, m, i, j; double dummy; /* Allocate workspace */ - isplit = (int *) malloc( n * sizeof(int) ); + isplit = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(isplit != NULL); - threads = (pthread_t *) malloc( max_nthreads * sizeof(pthread_t) ); + threads = (pthread_t *) malloc( (size_t)max_nthreads * sizeof(pthread_t) ); assert(threads != NULL); /* This is an unreduced block */ @@ -546,36 +544,34 @@ int eigval_approx_proc(proc_t *procinfo, int ifirst, int ilast, static -int eigval_root_proc(proc_t *procinfo, int ifirst, int ilast, - int n, double *D, double *E, double *E2, - int *Windex, int *iblock, double *gersch, tol_t *tolstruct, +PMRRR_Int eigval_root_proc(proc_t *procinfo, PMRRR_Int ifirst, PMRRR_Int ilast, + PMRRR_Int n, double *D, double *E, double *E2, + PMRRR_Int *Windex, PMRRR_Int *iblock, double *gersch, tol_t *tolstruct, double *W, double *Werr, double *Wgap, double *work, - int *iwork) + PMRRR_Int *iwork) { /* Input parameter */ - int pid = procinfo->pid; - /* int isize = ilast-ifirst+1; */ - double pivmin = tolstruct->pivmin; + double pivmin = tolstruct->pivmin; /* Tolerances */ double rtl; /* Create random vector to perturb rrr, same seed */ - int two_n = 2*n; - int iseed[4] = {1,1,1,1}; + PMRRR_Int two_n = 2*n; + PMRRR_Int iseed[4] = {1,1,1,1}; double *randvec; double isleft, isright, spdiam; double sigma, s1, s2; - int sgndef, cnt, negcnt_lft, negcnt_rgt; + PMRRR_Int sgndef, cnt, negcnt_lft, negcnt_rgt; double tau; - int jtry, off_L, off_invD; + PMRRR_Int jtry, off_L, off_invD; double Dpivot, Dmax; bool noREP; - int info, i, j; - int IONE = 1, ITWO = 2; + PMRRR_Int info, i, j; + PMRRR_Int IONE = 1, ITWO = 2; double tmp, tmp1, tmp2; double gl, gu; @@ -583,7 +579,7 @@ int eigval_root_proc(proc_t *procinfo, int ifirst, int ilast, rtl = sqrt(DBL_EPSILON); /* Allocate workspace */ - randvec = (double *) malloc( 2*n * sizeof(double) ); + randvec = (double *) malloc( 2*(size_t)n * sizeof(double) ); assert(randvec != NULL); /* create random vector to perturb rrr and broadcast it */ @@ -644,7 +640,7 @@ int eigval_root_proc(proc_t *procinfo, int ifirst, int ilast, /* define increment to perturb initial shift to find RRR * with not too much element growth */ - tau = spdiam*DBL_EPSILON*n + 2.0*pivmin; + tau = spdiam*DBL_EPSILON*(double)n + 2.0*pivmin; /* try to find initial RRR of block: @@ -686,17 +682,17 @@ int eigval_root_proc(proc_t *procinfo, int ifirst, int ilast, * so we should not end here */ if (jtry == MAX_TRY_RRR-2) { if (sgndef == ONE) { /* floating point comparison okay here */ - sigma = gl - FUDGE_FACTOR*spdiam*DBL_EPSILON*n + sigma = gl - FUDGE_FACTOR*spdiam*DBL_EPSILON*(double)n - FUDGE_FACTOR*2.0*pivmin; } else { - sigma = gu + FUDGE_FACTOR*spdiam*DBL_EPSILON*n + sigma = gu + FUDGE_FACTOR*spdiam*DBL_EPSILON*(double)n + FUDGE_FACTOR*2.0*pivmin; } } else if (jtry == MAX_TRY_RRR-1) { fprintf(stderr,"No initial representation could be found.\n"); exit(3); } else { - sigma -= sgndef*tau; + sigma -= (double)sgndef*tau; tau *= 2.0; continue; } @@ -707,8 +703,8 @@ int eigval_root_proc(proc_t *procinfo, int ifirst, int ilast, /* end trying to find initial RRR of block */ /* save initial RRR and corresponding shift */ - memcpy(D, &work[0], n * sizeof(double) ); - memcpy(E, &work[n], (n-1) * sizeof(double) ); + memcpy(D, &work[0], (size_t)n * sizeof(double) ); + memcpy(E, &work[n], (size_t)(n-1) * sizeof(double) ); E[n-1] = sigma; /* work[0:4*n-1] can now be used again for anything */ @@ -736,48 +732,46 @@ int eigval_root_proc(proc_t *procinfo, int ifirst, int ilast, static -int eigval_refine_proc(proc_t *procinfo, int ifirst, int ilast, - int n, double *D, double *E, double *E2, - int *Windex, int *iblock, double *gersch, tol_t *tolstruct, +PMRRR_Int eigval_refine_proc(proc_t *procinfo, PMRRR_Int ifirst, PMRRR_Int ilast, + PMRRR_Int n, double *D, double *E, double *E2, + PMRRR_Int *Windex, PMRRR_Int *iblock, double *gersch, tol_t *tolstruct, double *W, double *Werr, double *Wgap, double *work, - int *iwork) + PMRRR_Int *iwork) { /* Input parameter */ - int pid = procinfo->pid; - int isize = ilast-ifirst+1; - double pivmin = tolstruct->pivmin; + PMRRR_Int isize = ilast-ifirst+1; + double pivmin = tolstruct->pivmin; /* double gl, gu, wl, wu; */ double gl, gu; /* Multithreading */ - int nthreads; - int max_nthreads = procinfo->nthreads; - int iifirst, iilast, chunk; + PMRRR_Int nthreads; + PMRRR_Int max_nthreads = procinfo->nthreads; + PMRRR_Int chunk; pthread_t *threads; pthread_attr_t attr; auxarg2_t *auxarg2; void *status; /* Others */ - int nsplit, *isplit; + PMRRR_Int *isplit; double spdiam; - int i_low, i_upp; + PMRRR_Int i_low, i_upp; double sigma; - int off_DE2, offset; - int rf_begin, rf_end; + PMRRR_Int off_DE2, offset; + PMRRR_Int rf_begin, rf_end; - int info, i; + PMRRR_Int info, i; /* Allocate space */ - threads = (pthread_t *) malloc( max_nthreads * sizeof(pthread_t) ); + threads = (pthread_t *) malloc( (size_t)max_nthreads * sizeof(pthread_t) ); assert(threads != NULL); - isplit = (int *) malloc( n * sizeof(int) ); + isplit = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(isplit != NULL); /* This is an unreduced block */ - nsplit = 1; isplit[0] = n; /* Prepare multi-threading */ @@ -902,19 +896,19 @@ static void *eigval_subset_thread_a(void *argin) { /* from input argument */ - int n, il, iu, my_il, my_iu; + PMRRR_Int n, il, iu, my_il, my_iu; double *D, *E, *E2, *gersch; double bsrtol, pivmin; - int nsplit, *isplit; + PMRRR_Int nsplit, *isplit; /* others */ - int info; + PMRRR_Int info; double dummy1, dummy2; - int num_vals; + PMRRR_Int num_vals; double *W_tmp, *Werr_tmp, *W, *Werr; - int *iblock_tmp, *Windex_tmp, *iblock, *Windex; + PMRRR_Int *iblock_tmp, *Windex_tmp, *iblock, *Windex; double *work; - int *iwork; + PMRRR_Int *iwork; retrieve_auxarg1((auxarg1_t *) argin, &n, &D, &E, &E2, &il, &iu, &my_il, &my_iu, &nsplit, @@ -922,22 +916,22 @@ void *eigval_subset_thread_a(void *argin) &W, &Werr, &Windex, &iblock); /* allocate memory */ - W_tmp = (double *) malloc( n * sizeof(double) ); + W_tmp = (double *) malloc( (size_t)n * sizeof(double) ); assert(W_tmp != NULL); - Werr_tmp = (double *) malloc( n * sizeof(double) ); + Werr_tmp = (double *) malloc( (size_t)n * sizeof(double) ); assert(Werr_tmp != NULL); - Windex_tmp = (int *) malloc( n * sizeof(int) ); + Windex_tmp = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(Windex_tmp != NULL); - iblock_tmp = (int *) malloc( n * sizeof(int) ); + iblock_tmp = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(iblock_tmp != NULL); - work = (double *) malloc( 4*n * sizeof(double) ); + work = (double *) malloc( 4*(size_t)n * sizeof(double) ); assert (work != NULL); - iwork = (int *) malloc( 3*n * sizeof(int) ); + iwork = (PMRRR_Int *) malloc( 3*(size_t)n * sizeof(PMRRR_Int) ); assert (iwork != NULL); /* compute eigenvalues 'my_il' to 'my_iu', put into temporary arrays */ @@ -949,10 +943,10 @@ void *eigval_subset_thread_a(void *argin) assert(info == 0); /* copy computed values in W, Werr, Windex, iblock (which are work space) */ - memcpy(&W[my_il-il], W_tmp, num_vals * sizeof(double) ); - memcpy(&Werr[my_il-il], Werr_tmp, num_vals * sizeof(double) ); - memcpy(&Windex[my_il-il], Windex_tmp, num_vals * sizeof(int) ); - memcpy(&iblock[my_il-il], iblock_tmp, num_vals * sizeof(int) ); + memcpy(&W[my_il-il], W_tmp, (size_t)num_vals * sizeof(double) ); + memcpy(&Werr[my_il-il], Werr_tmp, (size_t)num_vals * sizeof(double) ); + memcpy(&Windex[my_il-il], Windex_tmp, (size_t)num_vals * sizeof(PMRRR_Int) ); + memcpy(&iblock[my_il-il], iblock_tmp, (size_t)num_vals * sizeof(PMRRR_Int) ); free(W_tmp); free(Werr_tmp); @@ -968,11 +962,11 @@ void *eigval_subset_thread_a(void *argin) static -auxarg1_t *create_auxarg1(int n, double *D, double *E, double *E2, - int il, int iu, int my_il, int my_iu, - int nsplit, int *isplit, double bsrtol, +auxarg1_t *create_auxarg1(PMRRR_Int n, double *D, double *E, double *E2, + PMRRR_Int il, PMRRR_Int iu, PMRRR_Int my_il, PMRRR_Int my_iu, + PMRRR_Int nsplit, PMRRR_Int *isplit, double bsrtol, double pivmin, double *gersch, double *W, - double *Werr, int *Windex, int *iblock) + double *Werr, PMRRR_Int *Windex, PMRRR_Int *iblock) { auxarg1_t *arg; @@ -1004,12 +998,12 @@ auxarg1_t *create_auxarg1(int n, double *D, double *E, double *E2, static -void retrieve_auxarg1(auxarg1_t *arg, int *n, double **D, double **E, - double **E2, int *il, int *iu, int *my_il, - int *my_iu, int *nsplit, int **isplit, +void retrieve_auxarg1(auxarg1_t *arg, PMRRR_Int *n, double **D, double **E, + double **E2, PMRRR_Int *il, PMRRR_Int *iu, PMRRR_Int *my_il, + PMRRR_Int *my_iu, PMRRR_Int *nsplit, PMRRR_Int **isplit, double *bsrtol, double *pivmin, double **gersch, - double **W, double **Werr, int **Windex, - int **iblock) + double **W, double **Werr, PMRRR_Int **Windex, + PMRRR_Int **iblock) { *n = arg->n; *D = arg->D; @@ -1039,28 +1033,27 @@ static void *eigval_subset_thread_r(void *argin) { /* from input argument */ - int bl_size, rf_begin, rf_end; + PMRRR_Int bl_size, rf_begin, rf_end; double *D, *DE2; double rtol1, rtol2, pivmin; double bl_spdiam; - val_t *Wstruct; /* others */ - int info, offset; + PMRRR_Int info, offset; double *W, *Werr, *Wgap; - int *Windex; + PMRRR_Int *Windex; double *work; - int *iwork; + PMRRR_Int *iwork; retrieve_auxarg2((auxarg2_t *) argin, &bl_size, &D, &DE2, &rf_begin, &rf_end, &W, &Werr, &Wgap, &Windex, &rtol1, &rtol2, &pivmin, &bl_spdiam); /* malloc work space */ - work = (double *) malloc( 2*bl_size * sizeof(double) ); + work = (double *) malloc( 2*(size_t)bl_size * sizeof(double) ); assert(work != NULL); - iwork = (int *) malloc( 2*bl_size * sizeof(int) ); + iwork = (PMRRR_Int *) malloc( 2*(size_t)bl_size * sizeof(PMRRR_Int) ); assert(iwork != NULL); /* special case of only one eigenvalue */ @@ -1087,9 +1080,9 @@ void *eigval_subset_thread_r(void *argin) static -auxarg2_t *create_auxarg2(int bl_size, double *D, double *DE2, - int rf_begin, int rf_end, double *W, double *Werr, - double *Wgap, int *Windex, +auxarg2_t *create_auxarg2(PMRRR_Int bl_size, double *D, double *DE2, + PMRRR_Int rf_begin, PMRRR_Int rf_end, double *W, double *Werr, + double *Wgap, PMRRR_Int *Windex, double rtol1, double rtol2, double pivmin, double bl_spdiam) { @@ -1119,9 +1112,9 @@ auxarg2_t *create_auxarg2(int bl_size, double *D, double *DE2, static -void retrieve_auxarg2(auxarg2_t *arg, int *bl_size, double **D, - double **DE2, int *rf_begin, int *rf_end, - double **W, double **Werr, double **Wgap, int **Windex, +void retrieve_auxarg2(auxarg2_t *arg, PMRRR_Int *bl_size, double **D, + double **DE2, PMRRR_Int *rf_begin, PMRRR_Int *rf_end, + double **W, double **Werr, double **Wgap, PMRRR_Int **Windex, double *rtol1, double *rtol2, double *pivmin, double *bl_spdiam) { @@ -1141,21 +1134,3 @@ void retrieve_auxarg2(auxarg2_t *arg, int *bl_size, double **D, free(arg); } - - - -/* - * Compare function for using qsort() on an array - * of doubles - */ -static -int cmp(const void *a1, const void *a2) -{ - double arg1 = *(double *)a1; - double arg2 = *(double *)a2; - - if (arg1 < arg2) - return(-1); - else - return(1); -} diff --git a/SRC/plarrv.c b/SRC/plarrv.c index 8184ade..4c9c30d 100644 --- a/SRC/plarrv.c +++ b/SRC/plarrv.c @@ -58,19 +58,19 @@ #include "counter.h" -static int assign_to_proc(proc_t *procinfo, in_t *Dstruct, - val_t *Wstruct, vec_t *Zstruct, int *nzp, - int *myfirstp); +static PMRRR_Int assign_to_proc(proc_t *procinfo, in_t *Dstruct, + val_t *Wstruct, vec_t *Zstruct, PMRRR_Int *nzp, + PMRRR_Int *myfirstp); static int cmpa(const void*, const void*); -static int init_workQ(proc_t *procinfo, in_t *Dstruct, - val_t *Wstruct, int *nzp, +static PMRRR_Int init_workQ(proc_t *procinfo, in_t *Dstruct, + val_t *Wstruct, PMRRR_Int *nzp, workQ_t *workQ); static void *empty_workQ(void*); static workQ_t *create_workQ(); static void destroy_workQ(workQ_t*); -static auxarg3_t *create_auxarg3(int, proc_t*, val_t*, vec_t*, +static auxarg3_t *create_auxarg3(PMRRR_Int, proc_t*, val_t*, vec_t*, tol_t*, workQ_t*, counter_t*); -static void retrieve_auxarg3(auxarg3_t*, int*, proc_t**, val_t**, +static void retrieve_auxarg3(auxarg3_t*, PMRRR_Int*, proc_t**, val_t**, vec_t**, tol_t**, workQ_t**, counter_t**); @@ -79,14 +79,13 @@ static void retrieve_auxarg3(auxarg3_t*, int*, proc_t**, val_t**, /* * Computation of eigenvectors of a symmetric tridiagonal */ -int plarrv(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, - vec_t *Zstruct, tol_t *tolstruct, int *nzp, - int *myfirstp) +PMRRR_Int plarrv(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, + vec_t *Zstruct, tol_t *tolstruct, PMRRR_Int *nzp, + PMRRR_Int *myfirstp) { /* Input variables */ - int pid = procinfo->pid; - int nthreads = procinfo->nthreads; - int n = Dstruct->n; + PMRRR_Int nthreads = procinfo->nthreads; + PMRRR_Int n = Dstruct->n; double *W = Wstruct->W; /* Work space */ @@ -100,18 +99,18 @@ int plarrv(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, counter_t *num_left; /* Others */ - int info, i; + PMRRR_Int info, i; workQ_t *workQ; /* Allocate work space and copy eigenvalues */ - Wshifted = (double *) malloc( n * sizeof(double) ); + Wshifted = (double *) malloc( (size_t)n * sizeof(double) ); assert(Wshifted != NULL); - memcpy(Wshifted, W, n*sizeof(double)); + memcpy(Wshifted, W, (size_t)n*sizeof(double)); Wstruct->Wshifted = Wshifted; - threads = (pthread_t *) malloc(nthreads * sizeof(pthread_t)); + threads = (pthread_t *) malloc((size_t)nthreads * sizeof(pthread_t)); assert(threads != NULL); /* Assign eigenvectors to processes */ @@ -168,30 +167,30 @@ int plarrv(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, * Assign the computation of eigenvectors to the processes */ static -int assign_to_proc(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, - vec_t *Zstruct, int *nzp, int *myfirstp) +PMRRR_Int assign_to_proc(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, + vec_t *Zstruct, PMRRR_Int *nzp, PMRRR_Int *myfirstp) { /* From inputs */ - int pid = procinfo->pid; - int nproc = procinfo->nproc; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int nproc = procinfo->nproc; double *restrict L = Dstruct->E; - int *restrict isplit = Dstruct->isplit; - int n = Wstruct->n; - int il = *(Wstruct->il); - int iu = *(Wstruct->iu); + PMRRR_Int *restrict isplit = Dstruct->isplit; + PMRRR_Int n = Wstruct->n; + PMRRR_Int il = *(Wstruct->il); + PMRRR_Int iu = *(Wstruct->iu); double *restrict W = Wstruct->W; - int *restrict Windex = Wstruct->Windex; - int *restrict iblock = Wstruct->iblock; - int *restrict iproc = Wstruct->iproc; - int *restrict Zindex = Zstruct->Zindex; + PMRRR_Int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict iblock = Wstruct->iblock; + PMRRR_Int *restrict iproc = Wstruct->iproc; + PMRRR_Int *restrict Zindex = Zstruct->Zindex; /* Others */ - int i, id, j, k, isize, iblk, ishift, + PMRRR_Int i, id, j, k, isize, iblk, ishift, ibegin, iend, chunk, ind; double sigma; sort_struct_t *array; - array = (sort_struct_t *) malloc(n*sizeof(sort_struct_t)); + array = (sort_struct_t *) malloc((size_t)n*sizeof(sort_struct_t)); for (i=0; iblock_ind == arg2->block_ind) { - return (arg1->local_ind - arg2->local_ind); + return(arg1->local_ind < arg2->local_ind ? -1 : 1); } else { if (arg1->lambda < arg2->lambda) { return(-1); @@ -298,42 +298,42 @@ int cmpa(const void *a1, const void *a2) * into the work queue. */ static -int init_workQ(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, - int *nzp, workQ_t *workQ) +PMRRR_Int init_workQ(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, + PMRRR_Int *nzp, workQ_t *workQ) { /* Input arguments */ - int pid = procinfo->pid; - int nproc = procinfo->nproc; - int nthreads = procinfo->nthreads; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int nproc = procinfo->nproc; + PMRRR_Int nthreads = procinfo->nthreads; double *restrict D = Dstruct->D; double *restrict L = Dstruct->E; - int nsplit = Dstruct->nsplit; - int *restrict isplit = Dstruct->isplit; + PMRRR_Int nsplit = Dstruct->nsplit; + PMRRR_Int *restrict isplit = Dstruct->isplit; double *restrict W = Wstruct->W; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; - int *restrict iproc = Wstruct->iproc; + PMRRR_Int *restrict iproc = Wstruct->iproc; double *restrict Wshifted = Wstruct->Wshifted; double *restrict gersch = Wstruct->gersch; - int nz = *nzp; + PMRRR_Int nz = *nzp; /* Loop over blocks */ - int ibegin, iend, isize, iWbegin, iWend, nbl; + PMRRR_Int ibegin, iend, isize, iWbegin, iWend, nbl; double sigma, gl, gu, avggap, spdiam; double *restrict DL; double *restrict DLL; rrr_t *RRR, *RRR_parent; /* Splitting into singletons and cluster */ - int new_first, new_last, new_size; - int sn_first, sn_last, sn_size; - int cl_first, cl_last, cl_size; + PMRRR_Int new_first, new_last, new_size; + PMRRR_Int sn_first, sn_last, sn_size; + PMRRR_Int cl_first, cl_last, cl_size; bool task_inserted; - int max_size, left_pid, right_pid; + PMRRR_Int max_size, left_pid, right_pid; double lgap; /* Others */ - int i, j, k, l; + PMRRR_Int i, j, k, l; double tmp; task_t *task; @@ -353,7 +353,7 @@ int init_workQ(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, gu = fmax(gu, gersch[2*i + 1]); } spdiam = gu - gl; - avggap = spdiam / (isize-1); + avggap = spdiam / (double)(isize-1); /* Find eigenvalues in block */ nbl = 0; @@ -382,10 +382,10 @@ int init_workQ(proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, /* Compute DL and DLL for later computation of singletons * (freed when all singletons of root are computed) */ - DL = (double *) malloc(isize * sizeof(double)); + DL = (double *) malloc((size_t)isize * sizeof(double)); assert(DL != NULL); - DLL = (double *) malloc(isize * sizeof(double)); + DLL = (double *) malloc((size_t)isize * sizeof(double)); assert(DLL != NULL); for (i=0; in; /* max. needed double precision work space: odr1v */ - work = (double *) malloc(4*n * sizeof(double)); + work = (double *) malloc(4*(size_t)n * sizeof(double)); assert(work != NULL); /* max. needed double precision work space: odrrb */ - iwork = (int *) malloc(2*n * sizeof(int) ); + iwork = (PMRRR_Int *) malloc(2*(size_t)n * sizeof(PMRRR_Int) ); assert(iwork != NULL); @@ -686,7 +687,7 @@ static void destroy_workQ(workQ_t *wq) static auxarg3_t* -create_auxarg3(int tid, proc_t *procinfo, val_t *Wstruct, +create_auxarg3(PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, workQ_t *workQ, counter_t *num_left) { @@ -710,7 +711,7 @@ create_auxarg3(int tid, proc_t *procinfo, val_t *Wstruct, static void -retrieve_auxarg3(auxarg3_t *arg, int *tid, proc_t **procinfo, +retrieve_auxarg3(auxarg3_t *arg, PMRRR_Int *tid, proc_t **procinfo, val_t **Wstruct, vec_t **Zstruct, tol_t **tolstruct, workQ_t **workQ, counter_t **num_left) diff --git a/SRC/pmrrr.c b/SRC/pmrrr.c index 439a348..34ca631 100644 --- a/SRC/pmrrr.c +++ b/SRC/pmrrr.c @@ -59,20 +59,20 @@ #include "structs.h" -static int handle_small_cases(char*, char*, int*, double*, double*, - double*, double*, int*, int*, int*, - MPI_Comm, int*, int*, double*, double*, - int*, int*); +static PMRRR_Int handle_small_cases(char*, char*, PMRRR_Int*, double*, double*, + double*, double*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, + MPI_Comm, PMRRR_Int*, PMRRR_Int*, double*, double*, + PMRRR_Int*, PMRRR_Int*); static int cmp(const void*, const void*); static int cmpb(const void*, const void*); static double scale_matrix(in_t*, val_t*, bool); -static void invscale_eigenvalues(val_t*, double, int); -static int sort_eigenpairs(proc_t*, val_t*, vec_t*); +static void invscale_eigenvalues(val_t*, double, PMRRR_Int); +static PMRRR_Int sort_eigenpairs(proc_t*, val_t*, vec_t*); static void clean_up(MPI_Comm, double*, double*, double*, - int*, int*, int*, int*, int*, proc_t*, + PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, PMRRR_Int*, proc_t*, in_t*, val_t*, vec_t*, tol_t*); -static int refine_to_highrac(proc_t*, char*, double*, double*, - in_t*, int*, val_t*, tol_t*); +static PMRRR_Int refine_to_highrac(proc_t*, char*, double*, double*, + in_t*, PMRRR_Int*, val_t*, tol_t*); @@ -84,14 +84,14 @@ static int refine_to_highrac(proc_t*, char*, double*, double*, * See README or 'pmrrr.h' for details. */ -int pmrrr(char *jobz, char *range, int *np, double *D, - double *E, double *vl, double *vu, int *il, - int *iu, int *tryracp, MPI_Comm comm, int *nzp, - int *offsetp, double *W, double *Z, int *ldz, - int *Zsupp) +PMRRR_Int pmrrr(char *jobz, char *range, PMRRR_Int *np, double *D, + double *E, double *vl, double *vu, PMRRR_Int *il, + PMRRR_Int *iu, PMRRR_Int *tryracp, MPI_Comm comm, PMRRR_Int *nzp, + PMRRR_Int *offsetp, double *W, double *Z, PMRRR_Int *ldz, + PMRRR_Int *Zsupp) { /* Input parameter */ - int n = *np; + PMRRR_Int n = *np; bool onlyW = (jobz[0] == 'N' || jobz[0] == 'n'); bool wantZ = (jobz[0] == 'V' || jobz[0] == 'v'); bool cntval = (jobz[0] == 'C' || jobz[0] == 'c'); @@ -101,7 +101,7 @@ int pmrrr(char *jobz, char *range, int *np, double *D, /* Work space */ double *Werr, *Wgap, *gersch, *Dcopy, *E2copy; - int *iblock, *iproc, *Windex, *Zindex, *isplit; + PMRRR_Int *iblock, *iproc, *Windex, *Zindex, *isplit; /* Structures to store variables */ proc_t *procinfo; @@ -118,8 +118,8 @@ int pmrrr(char *jobz, char *range, int *np, double *D, /* Others */ double scale; - int i, info; - int ifirst, ilast, isize, ifirst_tmp, ilast_tmp, chunk, iil, iiu; + PMRRR_Int i, info; + PMRRR_Int ifirst, ilast, isize, ifirst_tmp, ilast_tmp, chunk, iil, iiu; /* Check input parameters */ if(!(onlyW || wantZ || cntval)) return(1); @@ -160,7 +160,7 @@ int pmrrr(char *jobz, char *range, int *np, double *D, #if defined(MVAPICH2_VERSION) if (nthreads>1) { - int mv2_affinity=1; + PMRRR_Int mv2_affinity=1; char *mv2_string = getenv("MV2_ENABLE_AFFINITY"); if (mv2_string != NULL) { mv2_affinity = atoi(mv2_string); @@ -201,21 +201,21 @@ int pmrrr(char *jobz, char *range, int *np, double *D, } /* Allocate memory */ - Werr = (double *) malloc( n * sizeof(double) ); + Werr = (double *) malloc( (size_t)n * sizeof(double) ); assert(Werr != NULL); - Wgap = (double *) malloc( n * sizeof(double) ); + Wgap = (double *) malloc( (size_t)n * sizeof(double) ); assert(Wgap != NULL); - gersch = (double *) malloc( 2*n*sizeof(double) ); + gersch = (double *) malloc( 2*(size_t)n*sizeof(double) ); assert(gersch != NULL); - iblock = (int *) calloc( n , sizeof(int) ); + iblock = (PMRRR_Int *) calloc( (size_t)n , sizeof(PMRRR_Int) ); assert(iblock != NULL); - iproc = (int *) malloc( n * sizeof(int) ); + iproc = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(iproc != NULL); - Windex = (int *) malloc( n * sizeof(int) ); + Windex = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(Windex != NULL); - isplit = (int *) malloc( n * sizeof(int) ); + isplit = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(isplit != NULL); - Zindex = (int *) malloc( n * sizeof(int) ); + Zindex = (PMRRR_Int *) malloc( (size_t)n * sizeof(PMRRR_Int) ); assert(Zindex != NULL); procinfo = (proc_t *) malloc( sizeof(proc_t) ); assert(procinfo != NULL); @@ -259,6 +259,9 @@ int pmrrr(char *jobz, char *range, int *np, double *D, Zstruct->Zsupp = Zsupp; Zstruct->Zindex = Zindex; + isize = 0; + Dcopy = E2copy = NULL; + /* Scale matrix to allowable range, returns 1.0 if not scaled */ scale = scale_matrix(Dstruct, Wstruct, valeig); @@ -273,10 +276,10 @@ int pmrrr(char *jobz, char *range, int *np, double *D, /* This case is extremely rare in practice */ tolstruct->split = DBL_EPSILON; /* Copy original data needed for refinement later */ - Dcopy = (double *) malloc( n * sizeof(double) ); + Dcopy = (double *) malloc( (size_t)n * sizeof(double) ); assert(Dcopy != NULL); - memcpy(Dcopy, D, n*sizeof(double)); - E2copy = (double *) malloc( n * sizeof(double) ); + memcpy(Dcopy, D, (size_t)n*sizeof(double)); + E2copy = (double *) malloc( (size_t)n * sizeof(double) ); assert(E2copy != NULL); for (i=0; i 0) { - memmove(W, &W[ifirst-1], *nzp * sizeof(double)); + memmove(W, &W[ifirst-1], (size_t)(*nzp) * sizeof(double)); } /* If matrix was scaled, rescale eigenvalues */ @@ -396,8 +399,8 @@ int pmrrr(char *jobz, char *range, int *np, double *D, */ static void clean_up(MPI_Comm comm, double *Werr, double *Wgap, - double *gersch, int *iblock, int *iproc, - int *Windex, int *isplit, int *Zindex, + double *gersch, PMRRR_Int *iblock, PMRRR_Int *iproc, + PMRRR_Int *Windex, PMRRR_Int *isplit, PMRRR_Int *Zindex, proc_t *procinfo, in_t *Dstruct, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct) @@ -425,28 +428,30 @@ void clean_up(MPI_Comm comm, double *Werr, double *Wgap, * Wrapper to call LAPACKs DSTEMR for small matrices */ static -int handle_small_cases(char *jobz, char *range, int *np, double *D, - double *E, double *vlp, double *vup, int *ilp, - int *iup, int *tryracp, MPI_Comm comm, int *nzp, - int *myfirstp, double *W, double *Z, int *ldzp, - int *Zsupp) +PMRRR_Int handle_small_cases(char *jobz, char *range, PMRRR_Int *np, double *D, + double *E, double *vlp, double *vup, PMRRR_Int *ilp, + PMRRR_Int *iup, PMRRR_Int *tryracp, MPI_Comm comm, PMRRR_Int *nzp, + PMRRR_Int *myfirstp, double *W, double *Z, PMRRR_Int *ldzp, + PMRRR_Int *Zsupp) { bool cntval = (jobz[0] == 'C' || jobz[0] == 'c'); bool onlyW = (jobz[0] == 'N' || jobz[0] == 'n'); bool wantZ = (jobz[0] == 'V' || jobz[0] == 'v'); bool indeig = (range[0] == 'I' || range[0] == 'i'); - int n = *np; - int ldz_tmp = *np; - int ldz = *ldzp; + PMRRR_Int n = *np; + PMRRR_Int ldz_tmp = *np; + PMRRR_Int ldz = *ldzp; int nproc, pid; - int m, lwork, *iwork, liwork, info; + PMRRR_Int m, lwork, *iwork, liwork, info; double *Z_tmp, *work, cnt; - int i, itmp, MINUSONE=-1; - int chunk, myfirst, mylast, mysize; + PMRRR_Int i, itmp, MINUSONE=-1; + PMRRR_Int chunk, myfirst, mylast, mysize; MPI_Comm_size(comm, &nproc); MPI_Comm_rank(comm, &pid); + + Z_tmp = NULL; if (onlyW) { lwork = 12*n; @@ -459,15 +464,15 @@ int handle_small_cases(char *jobz, char *range, int *np, double *D, liwork = 10*n; if (indeig) itmp = *iup-*ilp+1; else itmp = n; - Z_tmp = (double *) malloc(n*itmp * sizeof(double)); + Z_tmp = (double *) malloc((size_t)n*(size_t)itmp * sizeof(double)); assert(Z_tmp != NULL); } else { return(1); } - work = (double *) malloc( lwork * sizeof(double)); + work = (double *) malloc( (size_t)lwork * sizeof(double)); assert(work != NULL); - iwork = (int *) malloc( liwork * sizeof(int)); + iwork = (PMRRR_Int *) malloc( (size_t)liwork * sizeof(PMRRR_Int)); assert(iwork != NULL); if (cntval) { @@ -478,7 +483,7 @@ int handle_small_cases(char *jobz, char *range, int *np, double *D, &liwork, &info); assert(info == 0); - *nzp = (int) ceil(cnt/nproc); + *nzp = (PMRRR_Int) ceil(cnt/nproc); free(work); free(iwork); return(0); } @@ -494,16 +499,16 @@ int handle_small_cases(char *jobz, char *range, int *np, double *D, mysize = mylast - myfirst + 1; if (mysize > 0) { - memmove(W, &W[myfirst], mysize*sizeof(double)); + memmove(W, &W[myfirst], (size_t)mysize*sizeof(double)); if (wantZ) { if (ldz == ldz_tmp) { /* copy everything in one chunk */ - memcpy(Z, &Z_tmp[myfirst*ldz_tmp], n*mysize*sizeof(double)); + memcpy(Z, &Z_tmp[myfirst*ldz_tmp], (size_t)n*(size_t)mysize*sizeof(double)); } else { /* copy each vector seperately */ for (i=0; in; + PMRRR_Int n = Dstruct->n; double *restrict D = Dstruct->D; double *restrict E = Dstruct->E; double *vl = Wstruct->vl; @@ -536,7 +541,7 @@ double scale_matrix(in_t *Dstruct, val_t *Wstruct, bool valeig) double scale = 1.0; double T_norm; double smlnum, bignum, rmin, rmax; - int IONE = 1, itmp; + PMRRR_Int IONE = 1, itmp; /* Set some machine dependent constants */ smlnum = DBL_MIN / DBL_EPSILON; @@ -575,13 +580,13 @@ double scale_matrix(in_t *Dstruct, val_t *Wstruct, bool valeig) */ static void invscale_eigenvalues(val_t *Wstruct, double scale, - int size) + PMRRR_Int size) { double *vl = Wstruct->vl; double *vu = Wstruct->vu; double *W = Wstruct->W; double invscale = 1.0 / scale; - int IONE = 1; + PMRRR_Int IONE = 1; if (scale != 1.0) { /* FP cmp okay */ *vl *= invscale; @@ -594,20 +599,19 @@ void invscale_eigenvalues(val_t *Wstruct, double scale, static -int sort_eigenpairs_local(proc_t *procinfo, int m, val_t *Wstruct, vec_t *Zstruct) +PMRRR_Int sort_eigenpairs_local(proc_t *procinfo, PMRRR_Int m, val_t *Wstruct, vec_t *Zstruct) { - int pid = procinfo->pid; - int n = Wstruct->n; + PMRRR_Int n = Wstruct->n; double *restrict W = Wstruct->W; double *restrict work = Wstruct->gersch; - int ldz = Zstruct->ldz; + PMRRR_Int ldz = Zstruct->ldz; double *restrict Z = Zstruct->Z; - int *restrict Zsupp = Zstruct->Zsupp; + PMRRR_Int *restrict Zsupp = Zstruct->Zsupp; bool sorted; - int j; + PMRRR_Int j; double tmp; - int itmp1, itmp2; + PMRRR_Int itmp1, itmp2; /* Make sure that sorted correctly; ineffective implementation, * but usually no or very little swapping should be done here */ @@ -629,9 +633,9 @@ int sort_eigenpairs_local(proc_t *procinfo, int m, val_t *Wstruct, vec_t *Zstruc Zsupp[2*j + 1] = Zsupp[2*(j+1) + 1]; Zsupp[2*(j+1) +1 ] = itmp2; /* swap eigenvector */ - memcpy(work, &Z[j*ldz], n*sizeof(double)); - memcpy(&Z[j*ldz], &Z[(j+1)*ldz], n*sizeof(double)); - memcpy(&Z[(j+1)*ldz], work, n*sizeof(double)); + memcpy(work, &Z[j*ldz], (size_t)n*sizeof(double)); + memcpy(&Z[j*ldz], &Z[(j+1)*ldz], (size_t)n*sizeof(double)); + memcpy(&Z[(j+1)*ldz], work, (size_t)n*sizeof(double)); } } } /* end while */ @@ -643,29 +647,29 @@ int sort_eigenpairs_local(proc_t *procinfo, int m, val_t *Wstruct, vec_t *Zstruc static -int sort_eigenpairs_global(proc_t *procinfo, int m, val_t *Wstruct, +PMRRR_Int sort_eigenpairs_global(proc_t *procinfo, PMRRR_Int m, val_t *Wstruct, vec_t *Zstruct) { - int pid = procinfo->pid; - int nproc = procinfo->nproc; - int n = Wstruct->n; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int nproc = procinfo->nproc; + PMRRR_Int n = Wstruct->n; double *restrict W = Wstruct->W; double *restrict work = Wstruct->gersch; - int ldz = Zstruct->ldz; + PMRRR_Int ldz = Zstruct->ldz; double *restrict Z = Zstruct->Z; - int *restrict Zsupp = Zstruct->Zsupp; + PMRRR_Int *restrict Zsupp = Zstruct->Zsupp; double *minW, *maxW, *minmax; - int i, p, lp, itmp[2]; + PMRRR_Int i, p, lp, itmp[2]; bool sorted; MPI_Status status; double nan_value = 0.0/0.0; - minW = (double *) malloc( nproc*sizeof(double)); + minW = (double *) malloc( (size_t)nproc*sizeof(double)); assert(minW != NULL); - maxW = (double *) malloc( nproc*sizeof(double)); + maxW = (double *) malloc( (size_t)nproc*sizeof(double)); assert(maxW != NULL); - minmax = (double *) malloc(2*nproc*sizeof(double)); + minmax = (double *) malloc(2*(size_t)nproc*sizeof(double)); assert(minmax != NULL); if (m == 0) { @@ -704,17 +708,17 @@ int sort_eigenpairs_global(proc_t *procinfo, int m, val_t *Wstruct, if ((pid == lp || pid == p) && minW[p] < maxW[lp]) { if (pid == lp) { W[m-1] = minW[p]; - MPI_Sendrecv(&Z[(m-1)*ldz], n, MPI_DOUBLE, p, lp, - work, n, MPI_DOUBLE, p, p, + MPI_Sendrecv(&Z[(m-1)*ldz], (int)n, MPI_DOUBLE, (int)p, (int)lp, + work, (int)n, MPI_DOUBLE, (int)p, (int)p, procinfo->comm, &status); - memcpy(&Z[(m-1)*ldz], work, n*sizeof(double)); + memcpy(&Z[(m-1)*ldz], work, (size_t)n*sizeof(double)); } if (pid == p) { W[0] = maxW[p-1]; - MPI_Sendrecv(&Z[0], n, MPI_DOUBLE, lp, p, - work, n, MPI_DOUBLE, lp, lp, + MPI_Sendrecv(&Z[0], (int)n, MPI_DOUBLE, (int)lp, (int)p, + work, (int)n, MPI_DOUBLE, (int)lp, (int)lp, procinfo->comm, &status); - memcpy(&Z[0], work, n*sizeof(double)); + memcpy(&Z[0], work, (size_t)n*sizeof(double)); } } @@ -722,15 +726,15 @@ int sort_eigenpairs_global(proc_t *procinfo, int m, val_t *Wstruct, * (would better be recomputed here though) */ if ((pid == lp || pid == p) && minW[p] < maxW[lp]) { if (pid == lp) { - MPI_Sendrecv(&Zsupp[2*(m-1)], 2, MPI_INT, p, lp, - itmp, 2, MPI_INT, p, p, + MPI_Sendrecv(&Zsupp[2*(m-1)], 2, PMRRR_MPI_INT_TYPE, (int)p, (int)lp, + itmp, 2, PMRRR_MPI_INT_TYPE, (int)p, (int)p, procinfo->comm, &status); Zsupp[2*(m-1)] = itmp[0]; Zsupp[2*(m-1) + 1] = itmp[1]; } if (pid == p) { - MPI_Sendrecv(&Zsupp[0], 2, MPI_INT, lp, p, - itmp, 2, MPI_INT, lp, lp, + MPI_Sendrecv(&Zsupp[0], 2, PMRRR_MPI_INT_TYPE, (int)lp, (int)p, + itmp, 2, PMRRR_MPI_INT_TYPE, (int)lp, (int)lp, procinfo->comm, &status); Zsupp[0] = itmp[0]; Zsupp[1] = itmp[1]; @@ -777,18 +781,18 @@ int sort_eigenpairs_global(proc_t *procinfo, int m, val_t *Wstruct, /* Routine to sort the eigenpairs */ static -int sort_eigenpairs(proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct) +PMRRR_Int sort_eigenpairs(proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct) { /* From inputs */ - int pid = procinfo->pid; - int n = Wstruct->n; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int n = Wstruct->n; double *restrict W = Wstruct->W; - int *restrict Windex = Wstruct->Windex; - int *restrict iproc = Wstruct->iproc; - int *restrict Zindex = Zstruct->Zindex; + PMRRR_Int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict iproc = Wstruct->iproc; + PMRRR_Int *restrict Zindex = Zstruct->Zindex; /* Others */ - int im, j; + PMRRR_Int im, j; sort_struct_t *sort_array; /* Make the first nz elements of W contains the eigenvalues @@ -803,7 +807,7 @@ int sort_eigenpairs(proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct) } } - sort_array = (sort_struct_t *) malloc(im*sizeof(sort_struct_t)); + sort_array = (sort_struct_t *) malloc((size_t)im*sizeof(sort_struct_t)); for (j=0; jpid; - bool wantZ = (jobz[0] == 'V' || jobz[0] == 'v'); - int n = Dstruct->n; - int nsplit = Dstruct->nsplit; - int *restrict isplit = Dstruct->isplit; + PMRRR_Int n = Dstruct->n; + PMRRR_Int nsplit = Dstruct->nsplit; + PMRRR_Int *restrict isplit = Dstruct->isplit; double spdiam = Dstruct->spdiam; - int nz = *nzp; double *restrict W = Wstruct->W; double *restrict Werr = Wstruct->Werr; - int *restrict Windex = Wstruct->Windex; - int *restrict iblock = Wstruct->iblock; - int *restrict iproc = Wstruct->iproc; double pivmin = tolstruct->pivmin; double tol = 4 * DBL_EPSILON; double *work; - int *iwork; - int ifirst, ilast, offset, info; - int i, j, k; - int ibegin, iend, isize, nbl; + PMRRR_Int *iwork; + PMRRR_Int ifirst, ilast, offset, info; + PMRRR_Int j; + PMRRR_Int ibegin, iend, isize, nbl; - work = (double *) malloc( 2*n * sizeof(double) ); + work = (double *) malloc( 2*(size_t)n * sizeof(double) ); assert (work != NULL); - iwork = (int *) malloc( 2*n * sizeof(int) ); + iwork = (PMRRR_Int *) malloc( 2*(size_t)n * sizeof(PMRRR_Int) ); assert (iwork != NULL); ibegin = 0; @@ -957,7 +955,7 @@ int cmp(const void *a1, const void *a2) * all computed eigenvalues (iu-il+1) in W; this routine is designed * to be called right after 'pmrrr'. */ -int PMR_comm_eigvals(MPI_Comm comm, int *nz, int *myfirstp, double *W) +PMRRR_Int PMR_comm_eigvals(MPI_Comm comm, PMRRR_Int *nz, PMRRR_Int *myfirstp, double *W) { MPI_Comm comm_dup; int nproc; @@ -967,22 +965,25 @@ int PMR_comm_eigvals(MPI_Comm comm, int *nz, int *myfirstp, double *W) MPI_Comm_dup(comm, &comm_dup); MPI_Comm_size(comm_dup, &nproc); - rcount = (int *) malloc( nproc * sizeof(int) ); + rcount = (int *) malloc( (size_t)nproc * sizeof(int) ); assert(rcount != NULL); - rdispl = (int *) malloc( nproc * sizeof(int) ); + rdispl = (int *) malloc( (size_t)nproc * sizeof(int) ); assert(rdispl != NULL); - work = (double *) malloc((*nz+1) * sizeof(double)); + work = (double *) malloc((size_t)(*nz+1) * sizeof(double)); assert(work != NULL); if (*nz > 0) { - memcpy(work, W, (*nz) * sizeof(double) ); + memcpy(work, W, (size_t)(*nz) * sizeof(double) ); } - MPI_Allgather(nz, 1, MPI_INT, rcount, 1, MPI_INT, comm_dup); + int nz32 = (int)(*nz); + int myfirstp32 = (int)(*myfirstp); + + MPI_Allgather(&nz32, 1, MPI_INT, rcount, 1, MPI_INT, comm_dup); - MPI_Allgather(myfirstp, 1, MPI_INT, rdispl, 1, MPI_INT, comm_dup); + MPI_Allgather(&myfirstp32, 1, MPI_INT, rdispl, 1, MPI_INT, comm_dup); - MPI_Allgatherv(work, *nz, MPI_DOUBLE, W, rcount, rdispl, + MPI_Allgatherv(work, (int)(*nz), MPI_DOUBLE, W, rcount, rdispl, MPI_DOUBLE, comm_dup); MPI_Comm_free(&comm_dup); @@ -997,10 +998,10 @@ int PMR_comm_eigvals(MPI_Comm comm, int *nz, int *myfirstp, double *W) /* Fortran function prototype */ -void pmrrr_(char *jobz, char *range, int *n, double *D, - double *E, double *vl, double *vu, int *il, int *iu, - int *tryracp, MPI_Fint *comm, int *nz, int *myfirst, - double *W, double *Z, int *ldz, int *Zsupp, int* info) +void pmrrr_(char *jobz, char *range, PMRRR_Int *n, double *D, + double *E, double *vl, double *vu, PMRRR_Int *il, PMRRR_Int *iu, + PMRRR_Int *tryracp, MPI_Fint *comm, PMRRR_Int *nz, PMRRR_Int *myfirst, + double *W, double *Z, PMRRR_Int *ldz, PMRRR_Int *Zsupp, PMRRR_Int* info) { MPI_Comm c_comm = MPI_Comm_f2c(*comm); @@ -1008,8 +1009,8 @@ void pmrrr_(char *jobz, char *range, int *n, double *D, c_comm, nz, myfirst, W, Z, ldz, Zsupp); } -void pmr_comm_eigvals_(MPI_Fint *comm, int *nz, int *myfirstp, - double *W, int *info) +void pmr_comm_eigvals_(MPI_Fint *comm, PMRRR_Int *nz, PMRRR_Int *myfirstp, + double *W, PMRRR_Int *info) { MPI_Comm c_comm = MPI_Comm_f2c(*comm); diff --git a/SRC/process_c_task.c b/SRC/process_c_task.c index c8963c4..a1fe110 100644 --- a/SRC/process_c_task.c +++ b/SRC/process_c_task.c @@ -61,27 +61,27 @@ static inline -rrr_t* compute_new_rrr(cluster_t *cl, int tid, proc_t *procinfo, +rrr_t* compute_new_rrr(cluster_t *cl, PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, - tol_t *tolstruct, double *work, int *iwork); + tol_t *tolstruct, double *work, PMRRR_Int *iwork); static inline -int refine_eigvals(cluster_t *cl, int rf_begin, int rf_end, - int tid, proc_t *procinfo, +PMRRR_Int refine_eigvals(cluster_t *cl, PMRRR_Int rf_begin, PMRRR_Int rf_end, + PMRRR_Int tid, proc_t *procinfo, rrr_t *RRR, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, counter_t *num_left, workQ_t *workQ, double *work, - int *iwork); + PMRRR_Int *iwork); static inline -int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, - int tid, val_t *Wstruct, rrr_t *RRR); +PMRRR_Int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, + PMRRR_Int tid, val_t *Wstruct, rrr_t *RRR); static inline -int test_comm_status(cluster_t *cl, val_t *Wstruct); +PMRRR_Int test_comm_status(cluster_t *cl, val_t *Wstruct); static inline -int create_subtasks(cluster_t *cl, int tid, proc_t *procinfo, +PMRRR_Int create_subtasks(cluster_t *cl, PMRRR_Int tid, proc_t *procinfo, rrr_t *RRR, val_t *Wstruct, vec_t *Zstruct, workQ_t *workQ, counter_t *num_left); @@ -89,22 +89,22 @@ int create_subtasks(cluster_t *cl, int tid, proc_t *procinfo, -int PMR_process_c_task(cluster_t *cl, int tid, proc_t *procinfo, +PMRRR_Int PMR_process_c_task(cluster_t *cl, PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, workQ_t *workQ, - counter_t *num_left, double *work, int *iwork) + counter_t *num_left, double *work, PMRRR_Int *iwork) { /* From inputs */ - int depth = cl->depth; - int left_pid = cl->left_pid; - int right_pid = cl->right_pid; - int pid = procinfo->pid; - int n = Wstruct->n; + PMRRR_Int depth = cl->depth; + PMRRR_Int left_pid = cl->left_pid; + PMRRR_Int right_pid = cl->right_pid; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int n = Wstruct->n; /* Others */ rrr_t *RRR; - int rf_begin, rf_end; - int status; + PMRRR_Int rf_begin, rf_end; + PMRRR_Int status; /* Protection against infinitely deep trees */ assert(depth < n); @@ -128,6 +128,8 @@ int PMR_process_c_task(cluster_t *cl, int tid, proc_t *procinfo, RRR = compute_new_rrr(cl, tid, procinfo, Wstruct, Zstruct, tolstruct, work, iwork); + rf_begin = rf_end = 0; + /* Refine eigenvalues 'rf_begin' to 'rf_end' */ if (left_pid != right_pid) { rf_begin = imax(cl->begin, cl->proc_W_begin); @@ -166,24 +168,24 @@ int PMR_process_c_task(cluster_t *cl, int tid, proc_t *procinfo, static inline -rrr_t* compute_new_rrr(cluster_t *cl, int tid, proc_t *procinfo, +rrr_t* compute_new_rrr(cluster_t *cl, PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, - tol_t *tolstruct, double *work, int *iwork) + tol_t *tolstruct, double *work, PMRRR_Int *iwork) { /* From inputs */ - int cl_begin = cl->begin; - int cl_end = cl->end; - int cl_size = cl_end - cl_begin + 1; - int depth = cl->depth; - int bl_begin = cl->bl_begin; - int bl_end = cl->bl_end; - int bl_size = bl_end - bl_begin + 1; + PMRRR_Int cl_begin = cl->begin; + PMRRR_Int cl_end = cl->end; + PMRRR_Int cl_size = cl_end - cl_begin + 1; + PMRRR_Int depth = cl->depth; + PMRRR_Int bl_begin = cl->bl_begin; + PMRRR_Int bl_end = cl->bl_end; + PMRRR_Int bl_size = bl_end - bl_begin + 1; double bl_spdiam = cl->bl_spdiam; rrr_t *RRR_parent = cl->RRR; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; - int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict Windex = Wstruct->Windex; double *restrict Wshifted = Wstruct->Wshifted; double pivmin = tolstruct->pivmin; @@ -199,21 +201,21 @@ rrr_t* compute_new_rrr(cluster_t *cl, int tid, proc_t *procinfo, double savegap; /* Others */ - int i, k, p, info; + PMRRR_Int i, k, p, info; double tmp; - int offset, IONE=1; + PMRRR_Int offset, IONE=1; /* Allocate memory for new representation for cluster */ - D = (double *) malloc(bl_size * sizeof(double)); + D = (double *) malloc((size_t)bl_size * sizeof(double)); assert(D != NULL); - L = (double *) malloc(bl_size * sizeof(double)); + L = (double *) malloc((size_t)bl_size * sizeof(double)); assert(L != NULL); - DL = (double *) malloc(bl_size * sizeof(double)); + DL = (double *) malloc((size_t)bl_size * sizeof(double)); assert(DL != NULL); - DLL = (double *) malloc(bl_size * sizeof(double)); + DLL = (double *) malloc((size_t)bl_size * sizeof(double)); assert(DLL != NULL); /* Recompute DL and DLL */ @@ -306,21 +308,21 @@ rrr_t* compute_new_rrr(cluster_t *cl, int tid, proc_t *procinfo, * Refine eigenvalues with respect to new rrr */ static inline -int refine_eigvals(cluster_t *cl, int rf_begin, int rf_end, - int tid, proc_t *procinfo, rrr_t *RRR, +PMRRR_Int refine_eigvals(cluster_t *cl, PMRRR_Int rf_begin, PMRRR_Int rf_end, + PMRRR_Int tid, proc_t *procinfo, rrr_t *RRR, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, counter_t *num_left, workQ_t *workQ, double *work, - int *iwork) + PMRRR_Int *iwork) { /* From inputs */ - int rf_size = rf_end-rf_begin+1; - int bl_begin = cl->bl_begin; - int bl_end = cl->bl_end; - int bl_size = bl_end - bl_begin + 1; + PMRRR_Int rf_size = rf_end-rf_begin+1; + PMRRR_Int bl_begin = cl->bl_begin; + PMRRR_Int bl_end = cl->bl_end; + PMRRR_Int bl_size = bl_end - bl_begin + 1; double bl_spdiam = cl->bl_spdiam; - int nthreads = procinfo->nthreads; + PMRRR_Int nthreads = procinfo->nthreads; double *restrict D = RRR->D; double *restrict L = RRR->L; @@ -329,29 +331,31 @@ int refine_eigvals(cluster_t *cl, int rf_begin, int rf_end, double *restrict W = Wstruct->W; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; - int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict Windex = Wstruct->Windex; double *restrict Wshifted = Wstruct->Wshifted; - int nz = Zstruct->nz; + PMRRR_Int nz = Zstruct->nz; double pivmin = tolstruct->pivmin; double rtol1 = tolstruct->rtol1; double rtol2 = tolstruct->rtol2; /* Others */ - int info, i, p, q, offset; + PMRRR_Int info, i, p, q, offset; double sigma, savegap; - int MIN_REFINE_CHUNK = fmax(2,nz/(4*nthreads)); - int left, own_part, others_part, num_tasks; - int ts_begin, ts_end, chunk, count; + PMRRR_Int MIN_REFINE_CHUNK = fmax(2,nz/(4*nthreads)); + PMRRR_Int left, own_part, others_part, num_tasks; + PMRRR_Int ts_begin, ts_end, chunk, count; task_t *task; sem_t sem; - int num_iter; + PMRRR_Int num_iter; /* Determine if refinement should be split into tasks */ left = PMR_get_counter_value(num_left); - own_part = (int) fmax( ceil( (double) left / nthreads ), - MIN_REFINE_CHUNK); + own_part = (PMRRR_Int) fmax( ceil( (double) left / (double) nthreads ), + (double)MIN_REFINE_CHUNK); + + savegap = 0.0; if (own_part < rf_size) { @@ -468,33 +472,33 @@ int refine_eigvals(cluster_t *cl, int rf_begin, int rf_end, static inline -int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, - int tid, val_t *Wstruct, rrr_t *RRR) +PMRRR_Int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, + PMRRR_Int tid, val_t *Wstruct, rrr_t *RRR) { /* From inputs */ - int cl_begin = cl->begin; - int cl_end = cl->end; - int bl_begin = cl->bl_begin; - int bl_end = cl->bl_end; - int proc_W_begin = cl->proc_W_begin; - int proc_W_end = cl->proc_W_end; - int left_pid = cl->left_pid; - int right_pid = cl->right_pid; - int num_messages; - // int num_messages = 4*(right_pid - left_pid); - - int pid = procinfo->pid; + PMRRR_Int cl_begin = cl->begin; + PMRRR_Int cl_end = cl->end; + PMRRR_Int bl_begin = cl->bl_begin; + PMRRR_Int bl_end = cl->bl_end; + PMRRR_Int proc_W_begin = cl->proc_W_begin; + PMRRR_Int proc_W_end = cl->proc_W_end; + PMRRR_Int left_pid = cl->left_pid; + PMRRR_Int right_pid = cl->right_pid; + PMRRR_Int num_messages; + // PMRRR_Int num_messages = 4*(right_pid - left_pid); + + PMRRR_Int pid = procinfo->pid; double *restrict W = Wstruct->W; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; double *restrict Wshifted = Wstruct->Wshifted; - int *restrict iproc = Wstruct->iproc; + PMRRR_Int *restrict iproc = Wstruct->iproc; /* Others */ - int p, i_msg, u, k, i; - int my_begin, my_end, my_size; - int other_begin, other_end, other_size; + PMRRR_Int p, i_msg, u, k, i; + PMRRR_Int my_begin, my_end, my_size; + PMRRR_Int other_begin, other_end, other_size; double sigma; int status, communication_done; MPI_Request *requests; @@ -507,6 +511,7 @@ int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, if (pid == left_pid ) my_begin = cl_begin; if (pid == right_pid) my_end = cl_end; my_size = my_end - my_begin + 1; + other_begin = 0; num_messages = 0; for (i=left_pid; i<=right_pid; i++) { @@ -518,9 +523,9 @@ int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, } } - requests = (MPI_Request *) malloc( num_messages * + requests = (MPI_Request *) malloc( (size_t)num_messages * sizeof(MPI_Request) ); - stats = (MPI_Status *) malloc( num_messages * + stats = (MPI_Status *) malloc( (size_t)num_messages * sizeof(MPI_Status) ); i_msg = 0; @@ -537,11 +542,11 @@ int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, if (p != pid && proc_involved == true) { /* send message to process p (non-blocking) */ - MPI_Isend(&Wshifted[my_begin], my_size, MPI_DOUBLE, p, - my_begin, procinfo->comm, &requests[4*i_msg]); + MPI_Isend(&Wshifted[my_begin], (int)my_size, MPI_DOUBLE, (int)p, + (int)my_begin, procinfo->comm, &requests[4*i_msg]); - MPI_Isend(&Werr[my_begin], my_size, MPI_DOUBLE, p, - my_begin, procinfo->comm, &requests[4*i_msg+1]); + MPI_Isend(&Werr[my_begin], (int)my_size, MPI_DOUBLE, (int)p, + (int)my_begin, procinfo->comm, &requests[4*i_msg+1]); /* Find eigenvalues in of process p */ other_size = 0; @@ -576,11 +581,11 @@ int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, } /* receive message from process p (non-blocking) */ - MPI_Irecv(&Wshifted[other_begin], other_size, MPI_DOUBLE, p, - other_begin, procinfo->comm, &requests[4*i_msg+2]); + MPI_Irecv(&Wshifted[other_begin], (int)other_size, MPI_DOUBLE, (int)p, + (int)other_begin, procinfo->comm, &requests[4*i_msg+2]); - MPI_Irecv(&Werr[other_begin], other_size, MPI_DOUBLE, p, - other_begin, procinfo->comm, &requests[4*i_msg+3]); + MPI_Irecv(&Werr[other_begin], (int)other_size, MPI_DOUBLE, (int)p, + (int)other_begin, procinfo->comm, &requests[4*i_msg+3]); i_msg++; } @@ -588,7 +593,7 @@ int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, } /* end for p */ num_messages = 4*i_msg; /* messages actually send */ - status = MPI_Testall(num_messages, requests, + status = MPI_Testall((int)num_messages, requests, &communication_done, stats); assert(status == MPI_SUCCESS); @@ -628,15 +633,15 @@ int communicate_refined_eigvals(cluster_t *cl, proc_t *procinfo, static inline -int test_comm_status(cluster_t *cl, val_t *Wstruct) +PMRRR_Int test_comm_status(cluster_t *cl, val_t *Wstruct) { - int cl_begin = cl->begin; - int cl_end = cl->end; - int bl_begin = cl->bl_begin; - int bl_end = cl->bl_end; + PMRRR_Int cl_begin = cl->begin; + PMRRR_Int cl_end = cl->end; + PMRRR_Int bl_begin = cl->bl_begin; + PMRRR_Int bl_end = cl->bl_end; rrr_t *RRR = cl->RRR; comm_t *comm = cl->messages; - int num_messages = comm->num_messages; + PMRRR_Int num_messages = comm->num_messages; MPI_Request *requests = comm->requests; MPI_Status *stats = comm->stats; double *restrict W = Wstruct->W; @@ -644,11 +649,12 @@ int test_comm_status(cluster_t *cl, val_t *Wstruct) double *restrict Wgap = Wstruct->Wgap; double *restrict Wshifted = Wstruct->Wshifted; - int status, k, communication_done; + int status, communication_done; + PMRRR_Int k; double sigma; /* Test if communication complete */ - status = MPI_Testall(num_messages, requests, + status = MPI_Testall((int)num_messages, requests, &communication_done, stats); assert(status == MPI_SUCCESS); @@ -680,45 +686,45 @@ int test_comm_status(cluster_t *cl, val_t *Wstruct) static inline -int create_subtasks(cluster_t *cl, int tid, proc_t *procinfo, +PMRRR_Int create_subtasks(cluster_t *cl, PMRRR_Int tid, proc_t *procinfo, rrr_t *RRR, val_t *Wstruct, vec_t *Zstruct, workQ_t *workQ, counter_t *num_left) { /* From inputs */ - int cl_begin = cl->begin; - int cl_end = cl->end; - int depth = cl->depth; - int bl_begin = cl->bl_begin; - int bl_end = cl->bl_end; - int bl_size = bl_end - bl_begin + 1; + PMRRR_Int cl_begin = cl->begin; + PMRRR_Int cl_end = cl->end; + PMRRR_Int depth = cl->depth; + PMRRR_Int bl_begin = cl->bl_begin; + PMRRR_Int bl_end = cl->bl_end; + PMRRR_Int bl_size = bl_end - bl_begin + 1; double bl_spdiam = cl->bl_spdiam; double lgap; - int pid = procinfo->pid; - int nproc = procinfo->nproc; - int nthreads = procinfo->nthreads; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int nproc = procinfo->nproc; + PMRRR_Int nthreads = procinfo->nthreads; bool proc_involved=true; double *restrict Wgap = Wstruct->Wgap; double *restrict Wshifted = Wstruct->Wshifted; - int *restrict iproc = Wstruct->iproc; + PMRRR_Int *restrict iproc = Wstruct->iproc; - int ldz = Zstruct->ldz; + PMRRR_Int ldz = Zstruct->ldz; double *restrict Z = Zstruct->Z; - int *restrict Zindex = Zstruct->Zindex; + PMRRR_Int *restrict Zindex = Zstruct->Zindex; /* others */ - int i, l, k; - int max_size; + PMRRR_Int i, l, k; + PMRRR_Int max_size; task_t *task; bool task_inserted; - int new_first, new_last, new_size, new_ftt1, new_ftt2; - int sn_first, sn_last, sn_size; + PMRRR_Int new_first, new_last, new_size, new_ftt1, new_ftt2; + PMRRR_Int sn_first, sn_last, sn_size; rrr_t *RRR_parent; - int new_lpid, new_rpid; + PMRRR_Int new_lpid, new_rpid; double *restrict D_parent; double *restrict L_parent; - int my_first, my_last; + PMRRR_Int my_first, my_last; bool copy_parent_rrr; @@ -726,6 +732,7 @@ int create_subtasks(cluster_t *cl, int tid, proc_t *procinfo, (fmin(depth+1,4)*nthreads) ); task_inserted = true; new_first = cl_begin; + sn_size = sn_first = sn_last = 0; for (i=cl_begin; i<=cl_end; i++) { if ( i == cl_end ) @@ -830,14 +837,14 @@ int create_subtasks(cluster_t *cl, int tid, proc_t *procinfo, if (copy_parent_rrr == true) { /* Copy parent RRR into alloceted arrays and mark them * for freeing later */ - D_parent = (double *) malloc(bl_size * sizeof(double)); + D_parent = (double *) malloc((size_t)bl_size * sizeof(double)); assert(D_parent != NULL); - L_parent = (double *) malloc(bl_size * sizeof(double)); + L_parent = (double *) malloc((size_t)bl_size * sizeof(double)); assert(L_parent != NULL); - memcpy(D_parent, RRR->D, bl_size*sizeof(double)); - memcpy(L_parent, RRR->L, bl_size*sizeof(double)); + memcpy(D_parent, RRR->D, (size_t)bl_size*sizeof(double)); + memcpy(L_parent, RRR->L, (size_t)bl_size*sizeof(double)); RRR_parent = PMR_create_rrr(D_parent, L_parent, NULL, NULL, bl_size, depth); @@ -846,9 +853,9 @@ int create_subtasks(cluster_t *cl, int tid, proc_t *procinfo, } else { /* copy parent RRR into Z to make cluster task independent */ memcpy(&Z[new_ftt1*ldz+bl_begin], RRR->D, - bl_size*sizeof(double)); + (size_t)bl_size*sizeof(double)); memcpy(&Z[new_ftt2*ldz+bl_begin], RRR->L, - bl_size*sizeof(double)); + (size_t)bl_size*sizeof(double)); RRR_parent = PMR_create_rrr(&Z[new_ftt1*ldz + bl_begin], &Z[new_ftt2*ldz + bl_begin], diff --git a/SRC/process_r_task.c b/SRC/process_r_task.c index 5c10edf..91dc8ae 100644 --- a/SRC/process_r_task.c +++ b/SRC/process_r_task.c @@ -54,8 +54,8 @@ #include "process_task.h" -int PMR_process_r_task(refine_t *rf, proc_t *procinfo, val_t *Wstruct, - tol_t *tolstruct, double *work, int *iwork); +PMRRR_Int PMR_process_r_task(refine_t *rf, proc_t *procinfo, val_t *Wstruct, + tol_t *tolstruct, double *work, PMRRR_Int *iwork); /* @@ -63,14 +63,14 @@ int PMR_process_r_task(refine_t *rf, proc_t *procinfo, val_t *Wstruct, * call. This routine is called to make sure that all tasks in the * queue are dequeued before continueing with other tasks. */ -void PMR_process_r_queue(int tid, proc_t *procinfo, val_t *Wstruct, +void PMR_process_r_queue(PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, workQ_t *workQ, counter_t *num_left, - double *work, int *iwork) + double *work, PMRRR_Int *iwork) { - int thread_support = procinfo->thread_support; - int t, num_tasks; - int status; + PMRRR_Int thread_support = procinfo->thread_support; + PMRRR_Int t, num_tasks; + PMRRR_Int status; task_t *task; num_tasks = PMR_get_num_tasks(workQ->r_queue); @@ -118,23 +118,23 @@ void PMR_process_r_queue(int tid, proc_t *procinfo, val_t *Wstruct, /* * Process the task of refining a subset of eigenvalues. */ -int PMR_process_r_task(refine_t *rf, proc_t *procinfo, +PMRRR_Int PMR_process_r_task(refine_t *rf, proc_t *procinfo, val_t *Wstruct, tol_t *tolstruct, - double *work, int *iwork) + double *work, PMRRR_Int *iwork) { /* From inputs */ - int ts_begin = rf->begin; + PMRRR_Int ts_begin = rf->begin; double *restrict D = rf->D; double *restrict DLL = rf->DLL; - int p = rf->p; - int q = rf->q; - int bl_size = rf->bl_size; + PMRRR_Int p = rf->p; + PMRRR_Int q = rf->q; + PMRRR_Int bl_size = rf->bl_size; double bl_spdiam = rf->bl_spdiam; sem_t *sem = rf->sem; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; - int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict Windex = Wstruct->Windex; double *restrict Wshifted = Wstruct->Wshifted; double rtol1 = tolstruct->rtol1; @@ -142,8 +142,8 @@ int PMR_process_r_task(refine_t *rf, proc_t *procinfo, double pivmin = tolstruct->pivmin; /* Others */ - int info, offset; - double savegap; + PMRRR_Int info, offset; + double savegap = 0.0; offset = Windex[ts_begin] - 1; diff --git a/SRC/process_s_task.c b/SRC/process_s_task.c index 57cff47..a540122 100644 --- a/SRC/process_s_task.c +++ b/SRC/process_s_task.c @@ -53,17 +53,17 @@ #include "process_task.h" -int PMR_process_s_task(singleton_t *sng, int tid, proc_t *procinfo, +PMRRR_Int PMR_process_s_task(singleton_t *sng, PMRRR_Int tid, proc_t *procinfo, val_t *Wstruct, vec_t *Zstruct, tol_t *tolstruct, counter_t *num_left, - double *work, int *iwork) + double *work, PMRRR_Int *iwork) { /* Inputs */ - int begin = sng->begin; - int end = sng->end; - int bl_begin = sng->bl_begin; - int bl_end = sng->bl_end; - int bl_size = bl_end - bl_begin + 1; + PMRRR_Int begin = sng->begin; + PMRRR_Int end = sng->end; + PMRRR_Int bl_begin = sng->bl_begin; + PMRRR_Int bl_end = sng->bl_end; + PMRRR_Int bl_size = bl_end - bl_begin + 1; double bl_spdiam = sng->bl_spdiam; rrr_t *RRR = sng->RRR; double *restrict D = RRR->D; @@ -71,43 +71,45 @@ int PMR_process_s_task(singleton_t *sng, int tid, proc_t *procinfo, double *restrict DL = RRR->DL; double *restrict DLL = RRR->DLL; - int pid = procinfo->pid; - int n = Wstruct->n; + PMRRR_Int pid = procinfo->pid; + PMRRR_Int n = Wstruct->n; double *restrict W = Wstruct->W; double *restrict Werr = Wstruct->Werr; double *restrict Wgap = Wstruct->Wgap; - int *restrict Windex = Wstruct->Windex; - int *restrict iproc = Wstruct->iproc; + PMRRR_Int *restrict Windex = Wstruct->Windex; + PMRRR_Int *restrict iproc = Wstruct->iproc; double *restrict Wshifted = Wstruct->Wshifted; - int ldz = Zstruct->ldz; + PMRRR_Int ldz = Zstruct->ldz; double *restrict Z = Zstruct->Z; - int *restrict isuppZ = Zstruct->Zsupp;; - int *restrict Zindex = Zstruct->Zindex; + PMRRR_Int *restrict isuppZ = Zstruct->Zsupp;; + PMRRR_Int *restrict Zindex = Zstruct->Zindex; double pivmin = tolstruct->pivmin; /* others */ - int info, i, k, itmp, num_decrement=0; - int IONE = 1; + PMRRR_Int info, i, k, itmp, num_decrement=0; + PMRRR_Int IONE = 1; double DZERO = 0.0; double tol, lambda, left, right; - int i_local, zind; + PMRRR_Int i_local, zind; double gap, lgap, rgap, gaptol, savedgap, tmp; bool usedBS, usedRQ, needBS, wantNC, step2II; - int r, offset; + PMRRR_Int r, offset; double twoeps = 2*DBL_EPSILON, RQtol = 2*DBL_EPSILON; double residual, bstres, bstw; - int i_supmn, i_supmx; + PMRRR_Int i_supmn, i_supmx; double RQcorr; - int negcount; - int sgndef, suppsize; + PMRRR_Int negcount; + PMRRR_Int sgndef, suppsize; double sigma; - int i_Zfrom, i_Zto; + PMRRR_Int i_Zfrom, i_Zto; double ztz, norminv, mingma; /* set tolerance parameter */ tol = 4.0 * log( (double) bl_size ) * DBL_EPSILON; + bstres = bstw = 0.0; + /* loop over all singletons in the task */ for (i=begin; i<=end; i++) { @@ -120,7 +122,7 @@ int PMR_process_s_task(singleton_t *sng, int tid, proc_t *procinfo, if (bl_size == 1) { /* set eigenvector to column of identity matrix */ zind = Zindex[i]; - memset(&Z[zind*ldz], 0.0, n*sizeof(double) ); + memset(&Z[zind*ldz], 0.0, (size_t)n*sizeof(double) ); Z[zind*ldz + bl_begin] = 1.0; isuppZ[2*zind ] = bl_begin + 1; isuppZ[2*zind + 1] = bl_begin + 1; @@ -174,7 +176,7 @@ int PMR_process_s_task(singleton_t *sng, int tid, proc_t *procinfo, /* IEEE floating point is assumed, so that all 0 bits are 0.0 */ zind = Zindex[i]; - memset(&Z[zind*ldz], 0.0, n*sizeof(double)); + memset(&Z[zind*ldz], 0.0, (size_t)n*sizeof(double)); /* inverse iteration with twisted factorization */ for (k=1; k<=MAXITER; k++) { @@ -228,7 +230,7 @@ int PMR_process_s_task(singleton_t *sng, int tid, proc_t *procinfo, sgndef = 1; /* wanted eigenvalue lies to the right */ } - if ( RQcorr*sgndef >= 0.0 + if ( RQcorr*(double)sgndef >= 0.0 && lambda+RQcorr <= right && lambda+RQcorr >= left ) { usedRQ = true; diff --git a/SRC/queue.c b/SRC/queue.c index 68c2b76..ac44140 100644 --- a/SRC/queue.c +++ b/SRC/queue.c @@ -49,7 +49,7 @@ queue_t *PMR_create_empty_queue(void) { - int info; + PMRRR_Int info; queue_t *queue; queue = (queue_t *) malloc(sizeof(queue_t)); @@ -83,9 +83,9 @@ void PMR_destroy_queue(queue_t *queue) -int PMR_insert_task_at_front(queue_t *queue, task_t *task) +PMRRR_Int PMR_insert_task_at_front(queue_t *queue, task_t *task) { - int info; + PMRRR_Int info; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&queue->lock); @@ -115,9 +115,9 @@ int PMR_insert_task_at_front(queue_t *queue, task_t *task) -int PMR_insert_task_at_back(queue_t *queue, task_t *task) +PMRRR_Int PMR_insert_task_at_back(queue_t *queue, task_t *task) { - int info; + PMRRR_Int info; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&queue->lock); @@ -150,7 +150,7 @@ int PMR_insert_task_at_back(queue_t *queue, task_t *task) task_t *PMR_remove_task_at_front(queue_t *queue) { - int info; + PMRRR_Int info; task_t *task; #ifdef NOSPINLOCKS @@ -191,7 +191,7 @@ task_t *PMR_remove_task_at_front(queue_t *queue) task_t *PMR_remove_task_at_back (queue_t *queue) { - int info; + PMRRR_Int info; task_t *task; #ifdef NOSPINLOCKS @@ -230,9 +230,9 @@ task_t *PMR_remove_task_at_back (queue_t *queue) -int PMR_get_num_tasks(queue_t *queue) +PMRRR_Int PMR_get_num_tasks(queue_t *queue) { - int info, num_tasks; + PMRRR_Int info, num_tasks; #ifdef NOSPINLOCKS info = pthread_mutex_lock(&queue->lock); diff --git a/SRC/rrr.c b/SRC/rrr.c index 8abfa87..d093eb3 100644 --- a/SRC/rrr.c +++ b/SRC/rrr.c @@ -48,9 +48,9 @@ rrr_t *PMR_create_rrr(double *restrict D, double *restrict L, double *restrict DL, double *restrict DLL, - int size, int depth) + PMRRR_Int size, PMRRR_Int depth) { - int info; + PMRRR_Int info; rrr_t *RRR; RRR = (rrr_t *) malloc( sizeof(rrr_t) ); @@ -76,7 +76,7 @@ rrr_t *PMR_create_rrr(double *restrict D, double *restrict L, rrr_t *PMR_reset_rrr(rrr_t *RRR, double *restrict D, double *restrict L, double *restrict DL, - double *restrict DLL, int size, int depth) + double *restrict DLL, PMRRR_Int size, PMRRR_Int depth) { RRR->D = D; RRR->L = L; @@ -91,10 +91,10 @@ rrr_t *PMR_reset_rrr(rrr_t *RRR, double *restrict D, -int PMR_increment_rrr_dependencies(rrr_t *RRR) +PMRRR_Int PMR_increment_rrr_dependencies(rrr_t *RRR) { /* returns number of dependencies */ - int i, info; + PMRRR_Int i, info; info = pthread_mutex_lock(&RRR->mutex); assert(info == 0); @@ -110,9 +110,9 @@ int PMR_increment_rrr_dependencies(rrr_t *RRR) -int PMR_set_parent_processed_flag(rrr_t *RRR) +PMRRR_Int PMR_set_parent_processed_flag(rrr_t *RRR) { - int info; + PMRRR_Int info; info = pthread_mutex_lock(&RRR->mutex); assert(info == 0); @@ -127,9 +127,9 @@ int PMR_set_parent_processed_flag(rrr_t *RRR) -int PMR_set_copied_parent_rrr_flag(rrr_t *RRR, bool val) +PMRRR_Int PMR_set_copied_parent_rrr_flag(rrr_t *RRR, bool val) { - int info; + PMRRR_Int info; info = pthread_mutex_lock(&RRR->mutex); assert(info == 0); @@ -144,11 +144,11 @@ int PMR_set_copied_parent_rrr_flag(rrr_t *RRR, bool val) -int PMR_try_destroy_rrr(rrr_t *RRR) +PMRRR_Int PMR_try_destroy_rrr(rrr_t *RRR) { /* return 0 on success, otherwise 1 */ - int info, tmp=0; + PMRRR_Int info, tmp=0; info = pthread_mutex_lock(&RRR->mutex); assert(info == 0); diff --git a/SRC/tasks.c b/SRC/tasks.c index 1d58905..dce87ec 100644 --- a/SRC/tasks.c +++ b/SRC/tasks.c @@ -46,8 +46,8 @@ #include "rrr.h" -task_t *PMR_create_s_task(int first, int last, int depth, - int bl_begin, int bl_end, +task_t *PMR_create_s_task(PMRRR_Int first, PMRRR_Int last, PMRRR_Int depth, + PMRRR_Int bl_begin, PMRRR_Int bl_end, double spdiam, double lgap, rrr_t *RRR) { task_t *t; @@ -77,10 +77,10 @@ task_t *PMR_create_s_task(int first, int last, int depth, } -task_t *PMR_create_c_task(int first, int last, int depth, - int bl_begin, int bl_end, double spdiam, - double lgap, int proc_W_begin, - int proc_W_end, int lpid, int rpid, +task_t *PMR_create_c_task(PMRRR_Int first, PMRRR_Int last, PMRRR_Int depth, + PMRRR_Int bl_begin, PMRRR_Int bl_end, double spdiam, + double lgap, PMRRR_Int proc_W_begin, + PMRRR_Int proc_W_end, PMRRR_Int lpid, PMRRR_Int rpid, rrr_t *RRR) { task_t *t; @@ -116,9 +116,9 @@ task_t *PMR_create_c_task(int first, int last, int depth, -task_t *PMR_create_r_task(int begin, int end, double *D, - double *DLL, int p, int q, int bl_size, - double bl_spdiam, int tid, sem_t *sem) +task_t *PMR_create_r_task(PMRRR_Int begin, PMRRR_Int end, double *D, + double *DLL, PMRRR_Int p, PMRRR_Int q, PMRRR_Int bl_size, + double bl_spdiam, PMRRR_Int tid, sem_t *sem) { task_t *t; refine_t *r; diff --git a/make.inc b/make.inc index 2b35c2f..1766b2f 100644 --- a/make.inc +++ b/make.inc @@ -1,8 +1,9 @@ # Compiler for C and Fortran +MPICC = mpicc CC = $(MPICC) # Compiler flags -CFLAGS = -pthread -O3 -restrict +CFLAGS = -pthread -O3 # Archiver and flags used when building the archive AR = /usr/bin/ar