Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improved eigvals and modenumber #94

Merged
merged 2 commits into from
May 21, 2024
Merged

Improved eigvals and modenumber #94

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
20dd4df
Improved eigvals and modenumber
sofiemartins May 18, 2024
0f2044b
Formatting and fiixed factor for max eva
sofiemartins May 18, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
19 changes: 0 additions & 19 deletions ModeNumber/input_file
View file
Open in desktop

This file was deleted.

251 changes: 34 additions & 217 deletions ModeNumber/mk_eigvals.c
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
/***************************************************************************\
* Copyright (c) 2009, Claudio Pica *
* Copyright (c) 2009-2024, Claudio Pica, Sofie Martins *
* All rights reserved. *
\***************************************************************************/

Expand All @@ -23,6 +23,8 @@
#error This code does not work with the fermion twisting !!!
#endif

static double hevamass;

/* Mesons parameters */
typedef struct input_eigval {
/* EVA parameters */
Expand All @@ -33,9 +35,10 @@ typedef struct input_eigval {
double omega1; /* absolute precision */
double omega2; /* relative precision */
double mass; /* quenched mass */
char configlist[256];

/* for the reading function */
input_record_t read[8];
input_record_t read[9];

} input_eigval;

Expand All @@ -49,202 +52,28 @@ typedef struct input_eigval {
{ "absolute precision", "eva:omega1 = %lf", DOUBLE_T, &(varname).omega1 }, \
{ "relative precision", "eva:omega2 = %lf", DOUBLE_T, &(varname).omega2 }, \
{ "quark quenched mass", "eva:mass = %lf", DOUBLE_T, &(varname).mass }, \
{ "Configuration list", "eva:configlist = %s", STRING_T, &(varname).configlist }, \
{ NULL, NULL, INT_T, NULL } \
} \
}

char cnfg_filename[256] = "";
char list_filename[256] = "";
char input_filename[255] = "input_file";
char output_filename[255] = "eigval.out";
enum { UNKNOWN_CNFG, DYNAMICAL_CNFG, QUENCHED_CNFG };

input_eigval eig_var = init_input_eigval(eig_var);

typedef struct {
char string[256];
int t, x, y, z;
int nc, nf;
double b, m;
int n;
int type;
} filename_t;

int parse_cnfg_filename(char *filename, filename_t *fn) {
int hm;
char *tmp = NULL;
char *basename;

basename = filename;
while ((tmp = strchr(basename, '/')) != NULL) {
basename = tmp + 1;
}

#ifdef REPR_FUNDAMENTAL
#define repr_name "FUN"
#elif defined REPR_SYMMETRIC
#define repr_name "SYM"
#elif defined REPR_ANTISYMMETRIC
#define repr_name "ASY"
#elif defined REPR_ADJOINT
#define repr_name "ADJ"
#endif
hm = sscanf(basename, "%*[^_]_%dx%dx%dx%d%*[Nn]c%dr" repr_name "%*[Nn]f%db%lfm%lfn%d", &(fn->t), &(fn->x), &(fn->y),
&(fn->z), &(fn->nc), &(fn->nf), &(fn->b), &(fn->m), &(fn->n));
if (hm == 9) {
fn->m = -fn->m; /* invert sign of mass */
fn->type = DYNAMICAL_CNFG;
return DYNAMICAL_CNFG;
}
#undef repr_name

double kappa;
hm = sscanf(basename, "%dx%dx%dx%d%*[Nn]c%d%*[Nn]f%db%lfk%lfn%d", &(fn->t), &(fn->x), &(fn->y), &(fn->z), &(fn->nc),
&(fn->nf), &(fn->b), &kappa, &(fn->n));
if (hm == 9) {
fn->m = .5 / kappa - 4.;
fn->type = DYNAMICAL_CNFG;
return DYNAMICAL_CNFG;
}

hm = sscanf(basename, "%*[^_]_%dx%dx%dx%d%*[Nn]c%db%lfn%d", &(fn->t), &(fn->x), &(fn->y), &(fn->z), &(fn->nc), &(fn->b),
&(fn->n));
if (hm == 7) {
fn->type = QUENCHED_CNFG;
return QUENCHED_CNFG;
}

hm = sscanf(basename, "%dx%dx%dx%d%*[Nn]c%db%lfn%d", &(fn->t), &(fn->x), &(fn->y), &(fn->z), &(fn->nc), &(fn->b), &(fn->n));
if (hm == 7) {
fn->type = QUENCHED_CNFG;
return QUENCHED_CNFG;
}

fn->type = UNKNOWN_CNFG;
return UNKNOWN_CNFG;
}

void read_cmdline(int argc, char *argv[]) {
int i, ai = 0, ao = 0, ac = 0, al = 0, am = 0;
FILE *list = NULL;

for (i = 1; i < argc; i++) {
if (strcmp(argv[i], "-i") == 0) {
ai = i + 1;
} else if (strcmp(argv[i], "-o") == 0) {
ao = i + 1;
} else if (strcmp(argv[i], "-c") == 0) {
ac = i + 1;
} else if (strcmp(argv[i], "-l") == 0) {
al = i + 1;
} else if (strcmp(argv[i], "-m") == 0) {
am = i;
}
}

if (am != 0) {
print_compiling_info();
exit(0);
}

if (ao != 0) { strcpy(output_filename, argv[ao]); }
if (ai != 0) { strcpy(input_filename, argv[ai]); }

error((ac == 0 && al == 0) || (ac != 0 && al != 0), 1, "parse_cmdline [eigval.c]",
"Syntax: eigval { -c <config file> | -l <list file> } [-i <input file>] [-o <output file>] [-m]");

if (ac != 0) {
strcpy(cnfg_filename, argv[ac]);
strcpy(list_filename, "");
} else if (al != 0) {
strcpy(list_filename, argv[al]);
error((list = fopen(list_filename, "r")) == NULL, 1, "parse_cmdline [eigval.c]", "Failed to open list file\n");
error(fscanf(list, "%s", cnfg_filename) == 0, 1, "parse_cmdline [eigval.c]", "Empty list file\n");
fclose(list);
}
}

double hevamass = 0.;
void H2EVA(spinor_field *out, spinor_field *in) {
g5Dphi_sq(hevamass, out, in);
}
void HEVA(spinor_field *out, spinor_field *in) {
g5Dphi(hevamass, out, in);
}

int main(int argc, char *argv[]) {
int i, n;
char tmp[256];
FILE *list;
filename_t fpars;
FILE *list = NULL;

/* setup process id and communications */
read_cmdline(argc, argv);
setup_process(&argc, &argv);

read_input(glb_var.read, input_filename);

/* logger setup */
/* disable logger for MPI processes != 0 */
logger_setlevel(0, 50);
if (PID != 0) { logger_disable(); }
if (PID == 0) {
sprintf(tmp, ">%s", output_filename);
logger_stdout(tmp);
sprintf(tmp, "err_%d", PID);
freopen(tmp, "w", stderr);
}

print_compiling_info_short();
lprintf("MAIN", 0, "PId = %d [world_size: %d]\n\n", PID, WORLD_SIZE);
lprintf("MAIN", 0, "input file [%s]\n", input_filename);
lprintf("MAIN", 0, "output file [%s]\n", output_filename);
if (strcmp(list_filename, "") != 0) {
lprintf("MAIN", 0, "list file [%s]\n", list_filename);
} else {
lprintf("MAIN", 0, "cnfg file [%s]\n", cnfg_filename);
}

/* read & broadcast parameters */
parse_cnfg_filename(cnfg_filename, &fpars);

eig_var.mass = -20.;
read_input(eig_var.read, input_filename);
GLB_T = fpars.t;
GLB_X = fpars.x;
GLB_Y = fpars.y;
GLB_Z = fpars.z;
error(fpars.type == UNKNOWN_CNFG, 1, "eigval.c", "Bad name for a configuration file");
error(fpars.nc != NG, 1, "eigval.c", "Bad NG");

lprintf("MAIN", 0, "Gauge group: SU(%d)\n", NG);
lprintf("MAIN", 0, "Fermion representation: " REPR_NAME " [dim=%d]\n", NF);

if (fpars.type == DYNAMICAL_CNFG) { hevamass = fpars.m; }
if (fpars.type == QUENCHED_CNFG || eig_var.mass > -10.) { hevamass = eig_var.mass; }

/* setup communication geometry */
if (geometry_init() == 1) {
finalize_process();
return 0;
setup_gauge_fields();
read_input(glb_var.read, get_input_filename());
read_input(eig_var.read, get_input_filename());
lprintf("MAIN", 0, "list file: [%s]\n", eig_var.configlist);
if (strcmp(eig_var.configlist, "") != 0) {
error((list = fopen(eig_var.configlist, "r")) == NULL, 1, "main [mk_eigvals.c]", "Failed to open list file\n");
}

/* setup lattice geometry */
geometry_mpi_eo();
/* test_geometry_mpi_eo(); */

/* setup random numbers */
read_input(rlx_var.read, input_filename);
lprintf("MAIN", 0, "RLXD [%d,%d]\n", rlx_var.rlxd_level, rlx_var.rlxd_seed + MPI_PID);
rlxd_init(rlx_var.rlxd_level, rlx_var.rlxd_seed);

init_BCs(NULL);

/* alloc global gauge fields */
u_gauge = alloc_suNg_field(&glattice);
#ifdef ALLOCATE_REPR_GAUGE_FIELD
u_gauge_f = alloc_suNf_field(&glattice);
#endif
hevamass = eig_var.mass;

lprintf("MAIN", 0, "EVA Parameters:\n");
lprintf("MAIN", 0, "search space dimension (eva:nevt) = %d\n", eig_var.nevt);
Expand All @@ -255,80 +84,68 @@ int main(int argc, char *argv[]) {
lprintf("MAIN", 0, "relative precision (eva:omega2) = %e\n", eig_var.omega2);
lprintf("MAIN", 0, "mass = %f\n", hevamass);

list = NULL;
if (strcmp(list_filename, "") != 0) {
error((list = fopen(list_filename, "r")) == NULL, 1, "main [eigval.c]", "Failed to open list file\n");
}

/* EVA parameters */
double max, mupp;
double *eva_val;
int status, ie;
spinor_field *eva_vec, *eva_ws;

eva_val = malloc(sizeof(double) * eig_var.nevt);
#ifndef UPDATE_EO
eva_vec = alloc_spinor_field(eig_var.nevt + 1, &glattice);
#else
eva_vec = alloc_spinor_field(eig_var.nevt + 1, &glat_even);
#endif
eva_ws = eva_vec + eig_var.nevt;

mupp = fabs(hevamass + 4) + 4;
mupp *= mupp;
/* END of EVA parameters */

char cnfg_filename[256];
Timer clock;
timer_set(&clock);

i = 0;
while (1) {
if (list != NULL) {
if (fscanf(list, "%s", cnfg_filename) == 0 || feof(list)) { break; }
}

i++;

lprintf("MAIN", 0, "Configuration from %s\n", cnfg_filename);
/* NESSUN CHECK SULLA CONSISTENZA CON I PARAMETRI DEFINITI !!! */
read_gauge_field(cnfg_filename);
represent_gauge_field();

lprintf("TEST", 0, "<p> %1.6f\n", avr_plaquette());

int MVM = 0; /* counter for matrix-vector multiplications */

max_eigval(&H2EVA, &glattice, &max);
timer_lap(&clock);
#ifdef UPDATE_EO
max_eigval(&H2, &glat_even, &max);
#else
max_eigval(&H2, &glattice, &max);
#endif
lprintf("MAIN", 0, "MAXCHECK: cnfg=%e uppbound=%e diff=%e %s\n", max, mupp, mupp - max,
(mupp - max) < 0 ? "[FAILED]" : "[OK]");
max *= 1.1;

ie = eva(eig_var.nev, eig_var.nevt, 0, eig_var.kmax, eig_var.maxiter, max, eig_var.omega1, eig_var.omega2, &H2EVA,
eva_vec, eva_val, &status);
MVM += status;
ie = eva(eig_var.nev, eig_var.nevt, 0, eig_var.kmax, eig_var.maxiter, max, eig_var.omega1, eig_var.omega2, &H2, eva_vec,
eva_val, &status);
while (ie != 0) { /* if failed restart EVA */
lprintf("MAIN", 0, "Restarting EVA!\n");
ie = eva(eig_var.nev, eig_var.nevt, 2, eig_var.kmax, eig_var.maxiter, max, eig_var.omega1, eig_var.omega2, &H2EVA,
ie = eva(eig_var.nev, eig_var.nevt, 2, eig_var.kmax, eig_var.maxiter, max, eig_var.omega1, eig_var.omega2, &H2,
eva_vec, eva_val, &status);
MVM += status;
}

lprintf("MAIN", 0, "EVA MVM = %d\n", MVM);
for (n = 0; n < eig_var.nev; ++n) {
HEVA(&eva_ws[0], &eva_vec[n]);
H(&eva_ws[0], &eva_vec[n]);
lprintf("RESULT", 0, "Eig %d = %.15e %.15e\n", n, eva_val[n],
prod_re_spinor_field(&eva_ws[0], &eva_vec[n]) / sqnorm_spinor_field(&eva_vec[n]));
}

double elapsed = timer_lap(&clock) * 1e-6;
lprintf("TIMING", 0, "Eigenvalue determination for configuration [%s] done [%lf sec]", cnfg_filename, elapsed);
if (list == NULL) { break; }
}

if (list != NULL) { fclose(list); }

free_BCs();

free(eva_val);
free_spinor_field(eva_vec);

free_suNg_field(u_gauge);
#ifdef ALLOCATE_REPR_GAUGE_FIELD
free_suNf_field(u_gauge_f);
#endif

finalize_process();

return 0;
}
Loading
Loading