From 0809bc060f289b23f21683b0a35e43b9f067a24e Mon Sep 17 00:00:00 2001 From: Bartosz Kostrzewa Date: Fri, 15 Mar 2024 16:42:21 +0100 Subject: [PATCH 1/2] WIP: eigenvalue online measurement --- meas/eigenvalues.c | 101 ++++++++++++++++++++++++++++++++++++++++++++ meas/eigenvalues.h | 32 ++++++++++++++ meas/measurements.h | 20 ++++++++- 3 files changed, 152 insertions(+), 1 deletion(-) create mode 100644 meas/eigenvalues.c create mode 100644 meas/eigenvalues.h diff --git a/meas/eigenvalues.c b/meas/eigenvalues.c new file mode 100644 index 000000000..41dfce6e9 --- /dev/null +++ b/meas/eigenvalues.c @@ -0,0 +1,101 @@ +/*********************************************************************** + * + * Copyright (C) 2024 Bartosz Kostrzewa + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + * + ************************************************************************/ + +#ifdef HAVE_CONFIG_H +# include +#endif +#ifdef TM_USE_OMP +# include +#endif + +#include +#include + +#include "global.h" +#include "geometry_eo.h" +#include "fatal_error.h" +#include "measurements.h" +#include "operator.h" +#include "gettime.h" +#include "tm_debug_printf.h" +#include "misc_types.h" + +#ifdef TM_USE_QUDA +#include "quda_interface.h" +#endif + +void eigenvalues_measurement(const int traj, const int id, const int ieo) { + tm_stopwatch_push(&g_timers, __func__, ""); + init_operators(); + + if(no_operators < 1){ + tm_debug_printf(0, 0, "Error: no operators defined in input file, cannot perform eigenvalues online measurement!\n"); +#ifdef TM_USE_MPI + MPI_Finallize(); +#endif + exit(1); + } + + eig_param_t eig = measurement_list[id].eig; + + double * evals = calloc(eig.n_evals, sizeof(double)); + + for( int op_id = 0; op_id < no_operators; op_id++ ){ + operator * optr = &operator_list[op_id]; + + if( optr->type != TMWILSON || optr->type != WILSON || + optr->type != CLOVER || optr->type != DBTMWILSON || + optr->type != DBCLOVER ){ + tm_debug_printf(0, 0, "Error: only operator types WILSON, TMWILSON, CLOVER, DBTMWILSON and DBCLOVER are supported.\n"); +#ifdef TM_USE_MPI + MPI_Finalize(); +#endif + exit(1); + } + + op_backup_restore_globals(TM_BACKUP_GLOBALS); + op_set_globals(op_id); + + if( measurement_list[id].external_library == QUDA_LIB ){ +#ifdef TM_USE_QUDA + eigsolveQuda(evals, eig.n_evals, eig.tol, 1, 0, eig.max_iterations, eig.maxmin, + optr->eps_sq, optr->maxiter, eig.polydeg, eig.amin, eig.amax, eig.n_kr, + optr->solver, optr->solver, optr->rel_prec, ieo, + optr->sloppy_precision, optr->sloppy_precision, optr->compression_type, + optr->type == CLOVER || optr->type == TMWILSON || optr->type == WILSON); +#else + tm_debug_printf(0, 0, "Error: Attempted to use QUDA eigensolver but this build was not configured for QUDA usage.\n"); +#ifdef TM_USE_MPI + MPI_Finalize(); +#endif + exit(1); +#endif + } + + op_backup_restore_globals(TM_RESTORE_GLOBALS); + + } // loop over operators + + free(evals); + tm_stopwatch_pop(&g_timers, 0, 1, ""); + return; +} + diff --git a/meas/eigenvalues.h b/meas/eigenvalues.h new file mode 100644 index 000000000..fbf24284c --- /dev/null +++ b/meas/eigenvalues.h @@ -0,0 +1,32 @@ +/*********************************************************************** + * + * Copyright (C) 2024 Bartosz Kostrzewa + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + ***********************************************************************/ + +#ifndef _MEASURE_EIGENVALUES_H +#define _MEASURE_EIGENVALUES_H + +#include "su3.h" + +/* implements the eigenvalue online measurement + * writes into eigenvalues.mXY.IJKLMN + * where XY corresponds to the operator index and IJKLMN refers to the trajectory counter */ + +void eigenvalues_measurement(const int traj, const int id, const int ieo); + +#endif diff --git a/meas/measurements.h b/meas/measurements.h index 0f5d6d6ab..268960e37 100644 --- a/meas/measurements.h +++ b/meas/measurements.h @@ -22,6 +22,8 @@ #ifndef _MEASUREMENTS_H #define _MEASUREMENTS_H +#include "misc_types.h" + #define max_no_measurements 20 /* Give the measurement types an unambiguous ID*/ @@ -30,9 +32,23 @@ enum MEAS_TYPE { PIONNORM, POLYAKOV, ORIENTED_PLAQUETTES, - GRADIENT_FLOW + GRADIENT_FLOW, + EIGENVALUES }; +typedef struct { + int n_evals; + int blksize; + int blkwise; + int max_iterations; + double tol; + int maxmin; + int polydeg; + double amin; + double amax; + int n_kr; +} eig_param_t; + typedef struct { enum MEAS_TYPE type; int initialised; @@ -70,6 +86,8 @@ typedef struct { void (*measurefunc) (const int traj, const int id, const int ieo); ExternalLibrary external_library; + + eig_param_t eig; } measurement; From ab5005253a7893e756ebf40ce73f1546168a7eab Mon Sep 17 00:00:00 2001 From: Bartosz Kostrzewa Date: Fri, 15 Mar 2024 17:32:35 +0100 Subject: [PATCH 2/2] basics for eigenvalues online / offline measurement wired up --- meas/Makefile.in | 3 +- meas/eigenvalues.c | 10 +-- meas/measurements.c | 5 ++ meas/measurements.h | 2 +- read_input.l | 69 ++++++++++++++++++- .../sample-eigenvalues-measurement.input | 50 ++++++++++++++ 6 files changed, 130 insertions(+), 9 deletions(-) create mode 100644 sample-input/sample-eigenvalues-measurement.input diff --git a/meas/Makefile.in b/meas/Makefile.in index 1d04ce4e8..fa7333375 100644 --- a/meas/Makefile.in +++ b/meas/Makefile.in @@ -37,7 +37,8 @@ libmeas_TARGETS = measurements \ pion_norm \ polyakov_loop \ measure_clover_field_strength_observables \ - gradient_flow + gradient_flow \ + eigenvalues libmeas_OBJECTS = $(addsuffix .o, ${libmeas_TARGETS}) diff --git a/meas/eigenvalues.c b/meas/eigenvalues.c index 41dfce6e9..34d8a204c 100644 --- a/meas/eigenvalues.c +++ b/meas/eigenvalues.c @@ -49,7 +49,7 @@ void eigenvalues_measurement(const int traj, const int id, const int ieo) { if(no_operators < 1){ tm_debug_printf(0, 0, "Error: no operators defined in input file, cannot perform eigenvalues online measurement!\n"); #ifdef TM_USE_MPI - MPI_Finallize(); + MPI_Finalize(); #endif exit(1); } @@ -61,9 +61,9 @@ void eigenvalues_measurement(const int traj, const int id, const int ieo) { for( int op_id = 0; op_id < no_operators; op_id++ ){ operator * optr = &operator_list[op_id]; - if( optr->type != TMWILSON || optr->type != WILSON || - optr->type != CLOVER || optr->type != DBTMWILSON || - optr->type != DBCLOVER ){ + if( !(optr->type == TMWILSON || optr->type == WILSON || + optr->type == CLOVER || optr->type == DBTMWILSON || + optr->type == DBCLOVER) ){ tm_debug_printf(0, 0, "Error: only operator types WILSON, TMWILSON, CLOVER, DBTMWILSON and DBCLOVER are supported.\n"); #ifdef TM_USE_MPI MPI_Finalize(); @@ -76,7 +76,7 @@ void eigenvalues_measurement(const int traj, const int id, const int ieo) { if( measurement_list[id].external_library == QUDA_LIB ){ #ifdef TM_USE_QUDA - eigsolveQuda(evals, eig.n_evals, eig.tol, 1, 0, eig.max_iterations, eig.maxmin, + eigsolveQuda(evals, eig.n_evals, eig.tol, 1, 0, eig.max_iter, eig.maxmin, optr->eps_sq, optr->maxiter, eig.polydeg, eig.amin, eig.amax, eig.n_kr, optr->solver, optr->solver, optr->rel_prec, ieo, optr->sloppy_precision, optr->sloppy_precision, optr->compression_type, diff --git a/meas/measurements.c b/meas/measurements.c index b18d263f9..70862cfac 100644 --- a/meas/measurements.c +++ b/meas/measurements.c @@ -30,6 +30,7 @@ #include "default_input_values.h" #include "read_input.h" +#include "eigenvalues.h" #include "pion_norm.h" #include "correlators.h" #include "polyakov_loop.h" @@ -80,6 +81,10 @@ int init_measurements(){ measurement_list[i].measurefunc = &gradient_flow_measurement; } + if(measurement_list[i].type == EIGENVALUES) { + measurement_list[i].measurefunc = &eigenvalues_measurement; + } + measurement_list[i].id = i; } return(0); diff --git a/meas/measurements.h b/meas/measurements.h index 268960e37..810c38d2a 100644 --- a/meas/measurements.h +++ b/meas/measurements.h @@ -40,7 +40,7 @@ typedef struct { int n_evals; int blksize; int blkwise; - int max_iterations; + int max_iter; double tol; int maxmin; int polydeg; diff --git a/read_input.l b/read_input.l index 3b39524dc..24f8cbefe 100644 --- a/read_input.l +++ b/read_input.l @@ -553,6 +553,7 @@ static inline double fltlist_next_token(int * const list_end){ %x PLOOP %x ORIENTEDPLAQUETTESMEAS %x GRADIENTFLOWMEAS +%x EIGENVALUESMEAS %x REWEIGH %x REWSAMPLES @@ -3058,6 +3059,10 @@ static inline double fltlist_next_token(int * const list_end){ meas->type = GRADIENT_FLOW; strcpy(meas->name, "GRADIENTFLOW"); } + else if(strcmp(yytext, "EIGENVALUES")==0) { + meas->type = EIGENVALUES; + strcpy(meas->name, "EIGENVALUES"); + } else { fprintf(stderr, "Unknown measurement type %s in line %d\n", yytext, line_of_file); exit(1); @@ -3086,9 +3091,10 @@ static inline double fltlist_next_token(int * const list_end){ else if(meas->type == POLYAKOV) BEGIN(PLOOP); else if(meas->type == ORIENTED_PLAQUETTES) BEGIN(ORIENTEDPLAQUETTESMEAS); else if(meas->type == GRADIENT_FLOW) BEGIN(GRADIENTFLOWMEAS); + else if(meas->type == EIGENVALUES) BEGIN(EIGENVALUESMEAS); } -{ +{ ^EndMeasurement{SPC}* { if(myverbose) printf("Measurement with id %d parsed in line %d\n\n", meas->id, line_of_file); BEGIN(0); @@ -3142,6 +3148,65 @@ static inline double fltlist_next_token(int * const list_end){ } } +{ + {SPC}*Tolerance{EQL}{FLT} { + sscanf(yytext, " %[2a-zA-Z] = %lf", name, &c); + meas->eig.tol = c; + if(myverbose) printf(" Eigenvalue measurement tolerance set to %lf line %d, measurement id=%d\n", meas->eig.tol, line_of_file, meas->id); + } + {SPC}*NoEigenvalues{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + meas->eig.n_evals = a; + if(myverbose) printf(" Number of eigenvalues to be determined set to %d, line %d measurement id=%d\n", a, line_of_file, meas->id); + } + {SPC}*BlockSize{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + meas->eig.blksize = a; + if(myverbose) printf(" Block size for eigensolver set to %d, line %d measurement id=%d\n", a, line_of_file, meas->id); + } + {SPC}*Blockwise{EQL}yes { + meas->eig.blkwise = 1; + if(myverbose) printf(" Blockwise eigenvalue determination, line %d, measurement id=%d\n", line_of_file, meas->id); + } + {SPC}*Blockwise{EQL}no { + meas->eig.blkwise = 0; + if(myverbose) printf(" NO blockwise eigenvalue deterination, line %d, measurement id=%d\n", line_of_file, meas->id); + } + {SPC}*LargestEigenvalues{EQL}yes { + meas->eig.maxmin = 1; + if(myverbose) printf(" Largest eigenvalues to be determined, line %d, measurement id=%d\n", line_of_file, meas->id); + } + {SPC}*LargestEigenvalues{EQL}no { + meas->eig.maxmin = 0; + if(myverbose) printf(" Smallest eigenvalues to be determined, line %d, measurement id=%d\n", line_of_file, meas->id); + } + {SPC}*MaxIterations{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + meas->eig.max_iter = a; + if(myverbose) printf(" Maximum number of iterations for eigensolver set to %d, line %d measurement id=%d\n", a, line_of_file, meas->id); + } + {SPC}*PolynomialDegree{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + meas->eig.polydeg = a; + if(myverbose) printf(" Degree of polynomial for eigensolver set to %d, line %d measurement id=%d\n", a, line_of_file, meas->id); + } + {SPC}*PolyMin{EQL}{FLT} { + sscanf(yytext, " %[2a-zA-Z] = %lf", name, &c); + meas->eig.amin = c; + if(myverbose) printf(" Minimum eigenvalue to exclude using polynomial acceleration in eigensolver set to %e line %d, measurement id=%d\n", c, line_of_file, meas->id); + } + {SPC}*PolyMax{EQL}{FLT} { + sscanf(yytext, " %[2a-zA-Z] = %lf", name, &c); + meas->eig.amax = c; + if(myverbose) printf(" Maximum eigenvalue to exclude using polynomial acceleration in eigensolver set to %e line %d, measurement id=%d\n", c, line_of_file, meas->id); + } + {SPC}*KrylovSubspaceSize{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + meas->eig.n_kr = a; + if(myverbose) printf(" Krylov subspace size for polynomial acceleration in eigensolver set to %d, line %d measurement id=%d\n", a, line_of_file, meas->id); + } +} + { {SPC}*Direction{EQL}[03] { sscanf(yytext, " %[a-zA-Z] = %d", name, &a); @@ -3691,7 +3756,7 @@ static inline double fltlist_next_token(int * const list_end){ BEGIN(comment_caller); } -{SPC}*\n { +{SPC}*\n { line_of_file++; } <*>{SPC}*\n { diff --git a/sample-input/sample-eigenvalues-measurement.input b/sample-input/sample-eigenvalues-measurement.input new file mode 100644 index 000000000..80965e097 --- /dev/null +++ b/sample-input/sample-eigenvalues-measurement.input @@ -0,0 +1,50 @@ +# example input file for eigenvalue measurements using "offline_measurement" +# requires 2 8^4 gauge configuration conf.0000 and conf.0002 + +L=8 +T=8 + +DebugLevel = 5 +ompnumthreads=4 + +InitialStoreCounter = 0 +Measurements = 2 +# measurements will be carried out in nsave steps +# e.g. for conf.0000 and conf.0002 in this case +nsave=2 +2kappamu = 0.05 +kappa = 0.177 +BCAngleT = 1 +GaugeConfigInputFile = conf +UseEvenOdd = yes + +# the eigenvalues measurement requires at least ONE operator to be defined +# if multiple operators are defined, the measurement will loop over them +# and produce numbered output files +BeginMeasurement EIGENVALUES + Frequency = 1 + NoEigenvalues = 10 + LargestEigenvalues = yes + BlockSize = 1 + Blockwise = no + MaxIterations = 10000 + Tolerance = 1e-14 + PolynomialDegree = 0 + PolyMin = 0.0002 + PolyMax = 0.003 + KrylovSubspaceSize = 100 + UseExternalLibrary = quda +EndMeasurement + +# note: setting the solver to CGMMS will result in the CGMMS inversion taking place +# because the solver is not properly decoupled form the rest of the code +BeginOperator TMWILSON + 2kappaMu = 0.05 + kappa = 0.177 + UseEvenOdd = yes + Solver = CG + SolverPrecision = 1e-14 + MaxSolverIterations = 1000 + AddDownPropagator = no +EndOperator +