Skip to content

Commit

Permalink
fix issue 110
Browse files Browse the repository at this point in the history 
problem was a wrong datatype MPI_Allreduce in scalar_prod_bi.c

also removed CRAY ifdefs, as the machine doesn't exist anymore
  • Loading branch information
Carsten Urbach committed May 27, 2012
1 parent 2a4ae71 commit 7a853c9
Show file tree
Hide file tree
Showing 10 changed files with 8 additions and 1,246 deletions.
16 changes: 4 additions & 12 deletions eigenvalues_bi.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,6 @@
*
* You should have received a copy of the GNU General Public License
* along with tmLQCD. If not, see <http://www.gnu.org/licenses/>.
***********************************************************************/
/*******************************************************************************
*
* Here we compute the nr_of_eigenvalues lowest eigenvalues
* of (gamma5*D)^2. Therefore we use the arnoldi routines.
Expand Down Expand Up @@ -71,11 +69,9 @@ double eigenvalues_bi(int * nr_of_eigenvalues,
**********************/
int verbosity = g_debug_level, converged = 0, blocksize = 1, blockwise = 0;
int solver_it_max = 200, j_max, j_min;
/*int it_max = 10000;*/
double decay_min = 1.7, decay_max = 1.5, prec,
threshold_min = 1.e-3, threshold_max = 5.e-2,
startvalue, threshold, decay, returnvalue;
/* static int v0dim = 0; */
int v0dim = 0;

/**********************
Expand Down Expand Up @@ -139,7 +135,10 @@ double eigenvalues_bi(int * nr_of_eigenvalues,
if((g_proc_id==0) && (g_debug_level > 4)) {
printf(" Values of mu = %e mubar = %e eps = %e precision = %e \n \n", g_mu, g_mubar, g_epsbar, precision);
}


/* here n and lda are equal, because Q_Qdagger_ND_BI does an internal */
/* conversion to non _bi fields which are subject to xchange_fields */
/* so _bi fields do not need boundary */
jdher_bi((VOLUME)/2*sizeof(bispinor)/sizeof(_Complex double), (VOLUME)/2*sizeof(bispinor)/sizeof(_Complex double),
startvalue, prec,
(*nr_of_eigenvalues), j_max, j_min,
Expand All @@ -150,13 +149,6 @@ double eigenvalues_bi(int * nr_of_eigenvalues,
&returncode, maxmin, 1,
&Q_Qdagger_ND_BI);

/* IN THE LAST LINE, INSERT:
Q_Qdagger_ND_BI; Non-degenerate case - on 1 bispinor
Q_Qdagger_ND; Non-degenerate case - on 2 spinors
Qtm_pm_psi; Degenerate case - on 1 spinor
*/


*nr_of_eigenvalues = converged;

returnvalue = eigenvls_bi[0];
Expand Down
2 changes: 1 addition & 1 deletion linalg/scalar_prod_bi.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -75,7 +75,7 @@ _Complex double scalar_prod_bi(bispinor * const S, bispinor * const R, const int
kc = ks + kc;

#if defined MPI
MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, MPI_COMM_WORLD);
kc = ks;
#endif

Expand Down
66 changes: 0 additions & 66 deletions solver/gram-schmidt.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -42,45 +42,6 @@ static int ONE = 1;
*
*/

void IteratedClassicalGS_old(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[]) {
const double alpha = 0.5;

double vnrm_old;
int i, n2, isorth = 0;
_Complex double CMONE, CONE, CZERO;
#ifdef CRAY
char * cupl_c = "C", *cupl_n = "N";
_fcd fupl_c, fupl_n;
fupl_c = _cptofcd(cupl_c, strlen(cupl_c));
fupl_n = _cptofcd(cupl_n, strlen(cupl_n));
#else
char * fupl_c = "C", *fupl_n = "N";
#endif

n2 = 2*n;
CMONE = -1.;
CONE = 1.;
CZERO = 0.;

#ifdef HAVE_LAPACK
vnrm_old = _FT(dnrm2)(&n2, (double*) v, &ONE);
for (i = 0; !isorth && i < max_cgs_it; i ++) {
_FT(zgemv)(fupl_c, &n, &m, &CONE, A, &n, v, &ONE, &CZERO, work1, &ONE, 1);
_FT(zgemv)(fupl_n, &n, &m, &CMONE, A, &n, work1, &ONE, &CONE, v, &ONE, 1);

(*vnrm) = _FT(dnrm2)(&n2, (double*) v, &ONE);

isorth=((*vnrm) > alpha*vnrm_old);
vnrm_old = *vnrm;
}
#endif
if (i >= max_cgs_it) {
/* errorhandler(400,""); */
}
}


void IteratedClassicalGS(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[], int lda) {
const double alpha = 0.5;
Expand All @@ -89,13 +50,7 @@ void IteratedClassicalGS(_Complex double v[], double *vnrm, int n, int m, _Compl
int i, isorth = 0;
int j;
_Complex double CMONE, CONE;
#ifdef CRAY
char *cupl_n = "N";
_fcd fupl_n;
fupl_n = _cptofcd(cupl_n, strlen(cupl_n));
#else
char *fupl_n = "N";
#endif

CMONE = -1.;
CONE = 1.;
Expand Down Expand Up @@ -131,13 +86,7 @@ void IteratedClassicalGS_su3vect(_Complex double v[], double *vnrm, int n, int m
int j;
_Complex double CMONE, CONE;

#ifdef CRAY
char *cupl_n = "N";
_fcd fupl_n;
fupl_n = _cptofcd(cupl_n, strlen(cupl_n));
#else
char *fupl_n = "N";
#endif

CMONE = -1.;
CONE = 1.;
Expand Down Expand Up @@ -170,21 +119,6 @@ void IteratedClassicalGS_su3vect(_Complex double v[], double *vnrm, int n, int m
* Orthogonlaizes v with respect to span{A[:,1:m]}
*/

void ModifiedGS_old(_Complex double v[], int n, int m, _Complex double A[]){

int i;
_Complex double s;

for (i = 0; i < m; i ++) {
s = scalar_prod((spinor*)(A+i*n), (spinor*) v, n*sizeof(_Complex double)/sizeof(spinor), 1);
s = -s;
#ifdef HAVE_LAPACK
_FT(zaxpy)(&n, &s, A+i*n, &ONE, v, &ONE);
#endif
}
}


void ModifiedGS(_Complex double v[], int n, int m, _Complex double A[], int lda) {

int i;
Expand Down
2 changes: 0 additions & 2 deletions solver/gram-schmidt.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,13 +20,11 @@
#define _GRAM_SCHMIDT_H
#include <complex.h>

void IteratedClassicalGS_old(_Complex double v[], double *vnrm, int n, int m, _Complex double A[], _Complex double work1[]);
void IteratedClassicalGS(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[], int lda) ;
void IteratedClassicalGS_su3vect(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[], int lda);

void ModifiedGS_old(_Complex double v[], int n, int m, _Complex double A[]);
void ModifiedGS(_Complex double v[], int n, int m, _Complex double A[], int lda);
void ModifiedGS_su3vect(_Complex double v[], int n, int m, _Complex double A[], int lda);

Expand Down
63 changes: 1 addition & 62 deletions solver/gram-schmidt_bi.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,9 +38,6 @@
#include <complex.h>
#include "linalg_eo.h"
#include "linalg/blas.h"
#ifdef CRAY
#include <fortran.h>
#endif
#include "gram-schmidt_bi.h"

const int max_cgs_it_bi=5;
Expand All @@ -54,44 +51,6 @@ static int ONE = 1;
*
*/

void IteratedClassicalGS_bi_old(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[]) {
const double alpha = 0.5;

double vnrm_old;
int i, n2, isorth = 0;
_Complex double CMONE, CONE, CZERO;
#ifdef CRAY
char * cupl_c = "C", *cupl_n = "N";
_fcd fupl_c, fupl_n;
fupl_c = _cptofcd(cupl_c, strlen(cupl_c));
fupl_n = _cptofcd(cupl_n, strlen(cupl_n));
#else
char * fupl_c = "C", *fupl_n = "N";
#endif

n2 = 2*n;
CMONE = -1.;
CONE = 1.;
CZERO = 0.;

vnrm_old = _FT(dnrm2)(&n2, (double*) v, &ONE);

for (i = 0; !isorth && i < max_cgs_it_bi; i ++) {
_FT(zgemv)(fupl_c, &n, &m, &CONE, A, &n, v, &ONE, &CZERO, work1, &ONE, 1);
_FT(zgemv)(fupl_n, &n, &m, &CMONE, A, &n, work1, &ONE, &CONE, v, &ONE, 1);

(*vnrm) = _FT(dnrm2)(&n2, (double*) v, &ONE);

isorth=((*vnrm) > alpha*vnrm_old);
vnrm_old = *vnrm;
}
if (i >= max_cgs_it_bi) {
/* errorhandler(400,""); */
}
}


void IteratedClassicalGS_bi(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[], int lda) {
const double alpha = 0.5;
Expand All @@ -100,13 +59,7 @@ void IteratedClassicalGS_bi(_Complex double v[], double *vnrm, int n, int m, _Co
int i, isorth = 0;
int j;
_Complex double CMONE, CONE;
#ifdef CRAY
char *cupl_n = "N";
_fcd fupl_n;
fupl_n = _cptofcd(cupl_n, strlen(cupl_n));
#else
char *fupl_n = "N";
#endif

CMONE = -1.;
CONE = 1.;
Expand Down Expand Up @@ -136,27 +89,13 @@ void IteratedClassicalGS_bi(_Complex double v[], double *vnrm, int n, int m, _Co
* Orthogonlaizes v with respect to span{A[:,1:m]}
*/

void ModifiedGS_bi_old(_Complex double v[], int n, int m, _Complex double A[]){

int i;
_Complex double s;

for (i = 0; i < m; i ++) {
s = scalar_prod_bi((bispinor*)(A+i*n), (bispinor*) v, n*sizeof(_Complex double)/sizeof(bispinor));
s = -s;
_FT(zaxpy)(&n, &s, A+i*n, &ONE, v, &ONE);
}
}


void ModifiedGS_bi(_Complex double v[], int n, int m, _Complex double A[], int lda){

int i;
_Complex double s;

for (i = 0; i < m; i ++) {
s = scalar_prod_bi((bispinor*)(A+i*lda), (bispinor*) v, n*sizeof(_Complex double)/sizeof(bispinor));
s = -s;
s = -scalar_prod_bi((bispinor*)(A+i*lda), (bispinor*) v, n*sizeof(_Complex double)/sizeof(bispinor));
_FT(zaxpy)(&n, &s, A+i*lda, &ONE, v, &ONE);
}
}
Expand Down
2 changes: 0 additions & 2 deletions solver/gram-schmidt_bi.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,11 +20,9 @@
#define _GRAM_SCHMIDT_BI_H
#include <complex.h>

void IteratedClassicalGS_bi_old(_Complex double v[], double *vnrm, int n, int m, _Complex double A[], _Complex double work1[]);
void IteratedClassicalGS_bi(_Complex double v[], double *vnrm, int n, int m, _Complex double A[],
_Complex double work1[], int lda) ;

void ModifiedGS_bi_old(_Complex double v[], int n, int m, _Complex double A[]);
void ModifiedGS_bi(_Complex double v[], int n, int m, _Complex double A[], int lda);

#endif
34 changes: 1 addition & 33 deletions solver/jdher.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,9 +62,6 @@
#include "solver/gram-schmidt.h"
#include "solver/quicksort.h"
#include "jdher.h"
#ifdef CRAY
#include <fortran.h>
#endif

#define min(a,b)((a)<(b) ? (a) : (b))
#define max(a,b)((a)<(b) ? (b) : (a))
Expand Down Expand Up @@ -110,19 +107,7 @@ _Complex double * p_Q;
_Complex double * p_work;
matrix_mult p_A_psi;

/****************************************************************************
* *
* Of course on the CRAY everything is different :( !! *
* that's why we need something more.
* *
****************************************************************************/

#ifdef CRAY
char * cupl_u = "U", * cupl_c = "C", *cupl_n = "N", * cupl_a = "A", *cupl_v = "V", *cilaenv = "zhetrd", *cvu = "VU";
_fcd fupl_u, fupl_c, fupl_a, fupl_n, fupl_v, filaenv, fvu;
#else
static char * fupl_u = "U", *fupl_n = "N", * fupl_a = "A", *fupl_v = "V", *filaenv = "zhetrd", *fvu = "VU";
#endif
static char * fupl_u = "U", *fupl_n = "N", * fupl_a = "A", *fupl_v = "V", *filaenv = "zhetrd", *fvu = "VU";

/****************************************************************************
* *
Expand Down Expand Up @@ -177,23 +162,6 @@ void jdher(int n, int lda, double tau, double tol,
int ISEED[4] = {2, 3, 5, 7};
ISEED[0] = g_proc_id+2;

/****************************************************************************
* *
* Of course on the CRAY everything is different :( !! *
* that's why we need something more.
* *
****************************************************************************/

#ifdef CRAY
fupl_u = _cptofcd(cupl_u, strlen(cupl_u));
fupl_c = _cptofcd(cupl_c, strlen(cupl_c));
fupl_n = _cptofcd(cupl_n, strlen(cupl_n));
fupl_a = _cptofcd(cupl_a, strlen(cupl_a));
fupl_v = _cptofcd(cupl_v, strlen(cupl_v));
filaenv = _cptofcd(cilaenv, strlen(cilaenv));
fvu = _cptofcd(cvu, strlen(cvu));
#endif

/****************************************************************************
* *
* Execution starts here... *
Expand Down
Loading

0 comments on commit 7a853c9

Please sign in to comment.