From 3d67264370cce89aca9872bcb72a2337cd2c825f Mon Sep 17 00:00:00 2001 From: Aniket Sen Date: Sat, 17 Jun 2023 10:51:44 +0200 Subject: [PATCH] user input for poly acc settings --- default_input_values.h | 4 ++++ monomial/monomial.c | 5 +++++ monomial/monomial.h | 2 ++ phmc.c | 6 ++++-- quda_interface.c | 14 ++++++++------ quda_interface.h | 3 ++- read_input.l | 20 ++++++++++++++++++++ 7 files changed, 45 insertions(+), 9 deletions(-) diff --git a/default_input_values.h b/default_input_values.h index 64b1a3531..820f3d98e 100644 --- a/default_input_values.h +++ b/default_input_values.h @@ -147,6 +147,10 @@ #define _default_phmc_pure_phmc 0 #define _default_stilde_max 3. #define _default_stilde_min 0.01 +#define _default_eig_polydeg 128 +#define _default_eig_amin 0.001 +#define _default_eig_amax 4 +#define _default_eig_n_kr 96 #define _default_degree_of_p 48 #define _default_propagator_splitted 1 #define _default_source_splitted 1 diff --git a/monomial/monomial.c b/monomial/monomial.c index 4214fef5d..b7e8a55c5 100644 --- a/monomial/monomial.c +++ b/monomial/monomial.c @@ -144,6 +144,11 @@ int add_monomial(const int type) { monomial_list[no_monomials].PrecisionHfinal = _default_g_acc_Hfin; monomial_list[no_monomials].PrecisionPtilde = _default_g_acc_Ptilde; + monomial_list[no_monomials].eig_polydeg = _default_eig_polydeg; + monomial_list[no_monomials].eig_amin = _default_eig_amin; + monomial_list[no_monomials].eig_amax = _default_eig_amax; + monomial_list[no_monomials].eig_n_kr = _default_eig_n_kr; + monomial_list[no_monomials].rat.order = 12; monomial_list[no_monomials].rat.range[0] = _default_stilde_min; monomial_list[no_monomials].rat.range[1] = _default_stilde_max; diff --git a/monomial/monomial.h b/monomial/monomial.h index c8f7e8b40..bbdb30dd2 100644 --- a/monomial/monomial.h +++ b/monomial/monomial.h @@ -114,6 +114,8 @@ typedef struct { double PrecisionHfinal; double StildeMin, StildeMax; double EVMin, EVMax, EVMaxInv; + int eig_polydeg, eig_n_kr; + double eig_amin, eig_amax; ExternalLibrary external_library; ExternalEigSolver external_eigsolver; double * MDPolyCoefs, * PtildeCoefs; diff --git a/phmc.c b/phmc.c index 32900382e..cf4e73b14 100644 --- a/phmc.c +++ b/phmc.c @@ -231,7 +231,8 @@ void phmc_compute_ev(const int trajectory_counter, if(mnl->external_eigsolver == QUDA_EIGSOLVER) { #ifdef TM_USE_QUDA temp = eigsolveQuda(no_eigenvalues, eigenvalue_precision, 1, 0, max_iter_ev, 0, - mnl->accprec, mnl->maxiter, mnl->solver, g_relative_precision_flag, + mnl->accprec, mnl->maxiter, mnl->eig_polydeg, mnl->eig_amin, + mnl->eig_amax, mnl->eig_n_kr, mnl->solver, g_relative_precision_flag, 1, // we only support even-odd here mnl->solver_params.refinement_precision, mnl->solver_params.sloppy_precision, @@ -257,7 +258,8 @@ void phmc_compute_ev(const int trajectory_counter, if(mnl->external_eigsolver == QUDA_EIGSOLVER) { #ifdef TM_USE_QUDA temp2 = eigsolveQuda(no_eigenvalues, eigenvalue_precision, 1, 0, max_iter_ev, 1, - mnl->accprec, mnl->maxiter, mnl->solver, g_relative_precision_flag, + mnl->accprec, mnl->maxiter, mnl->eig_polydeg, mnl->eig_amin, + mnl->eig_amax, mnl->eig_n_kr, mnl->solver, g_relative_precision_flag, 1, // we only support even-odd here mnl->solver_params.refinement_precision, mnl->solver_params.sloppy_precision, diff --git a/quda_interface.c b/quda_interface.c index 4a985641c..c17b2015f 100644 --- a/quda_interface.c +++ b/quda_interface.c @@ -2598,7 +2598,8 @@ Interface function for Eigensolver on Quda double eigsolveQuda(int n, double tol, int blksize, int blkwise, int max_iterations, int maxmin, - const double precision, const int max_iter, const int solver_flag, const int rel_prec, + const double precision, const int max_iter, const int polydeg, const double amin, + const double amax, const int n_kr, const int solver_flag, const int rel_prec, const int even_odd_flag, const SloppyPrecision refinement_precision, SloppyPrecision sloppy_precision, CompressionType compression) { @@ -2646,6 +2647,7 @@ double eigsolveQuda(int n, double tol, int blksize, int blkwise, int max_iterati QudaInvertParam eig_invert_param = newQudaInvertParam(); eig_invert_param = inv_param; eig_param.invert_param = &eig_invert_param; + eig_param.invert_param->verbosity = QUDA_VERBOSE; /* AS The following two are set to cuda_prec, otherwise * it gives an error. Such high precision might not be * necessary. But have not found a way to consistently set @@ -2689,10 +2691,10 @@ double eigsolveQuda(int n, double tol, int blksize, int blkwise, int max_iterati // at a relatively coarse lattice spacing for the eigenvalues // of the twisted-clover ND operator with values of musigma / mudelta // reproducing physical sea strange and charm quark masses - eig_param.use_poly_acc = maxmin == 1 ? QUDA_BOOLEAN_FALSE : QUDA_BOOLEAN_TRUE; - eig_param.poly_deg = 128; - eig_param.a_min = 1e-3; - eig_param.a_max = 4; + eig_param.use_poly_acc = (maxmin == 1) && (polydeg != 0) ? QUDA_BOOLEAN_FALSE : QUDA_BOOLEAN_TRUE; + eig_param.poly_deg = polydeg; + eig_param.a_min = amin; + eig_param.a_max = amax; /* Daggers the operator. Not necessary for * most cases. */ @@ -2730,7 +2732,7 @@ double eigsolveQuda(int n, double tol, int blksize, int blkwise, int max_iterati * From my understanding, QUDA automatically scales * this search space, however more testing on this * might be necessary */ - eig_param.n_kr = 96; + eig_param.n_kr = n_kr; eig_param.max_restarts = max_iterations; diff --git a/quda_interface.h b/quda_interface.h index 5320fd5ea..13fb578c9 100644 --- a/quda_interface.h +++ b/quda_interface.h @@ -176,7 +176,8 @@ void compute_WFlow_quda(const double eps ,const double tmax, const int traj, FIL double eigsolveQuda(int n, double tol, int blksize, int blkwise, int max_iterations, int maxmin, - const double precision, const int max_iter, const int solver_flag, const int rel_prec, + const double precision, const int max_iter, const int polydeg, const double amin, + const double amax, const int n_kr, const int solver_flag, const int rel_prec, const int even_odd_flag, const SloppyPrecision refinement_precision, SloppyPrecision sloppy_precision, CompressionType compression); diff --git a/read_input.l b/read_input.l index 4e7b0db83..3fde529d0 100644 --- a/read_input.l +++ b/read_input.l @@ -2607,6 +2607,26 @@ static inline double fltlist_next_token(int * const list_end){ if(myverbose) printf(" Do not use external eigensolver line %d monomial %d\n", line_of_file, current_monomial); mnl->external_eigsolver = NO_EXT_EIGSOLVER; } + {SPC}*EigAmin{EQL}{FLT} { + sscanf(yytext, " %[a-zA-Z] = %lf", name, &c); + mnl->eig_amin = c; + if(myverbose!=0) printf(" eig_amin set to %e line %d monomial %d\n", c, line_of_file, current_monomial); + } + {SPC}*EigAmax{EQL}{FLT} { + sscanf(yytext, " %[a-zA-Z] = %lf", name, &c); + mnl->eig_amax = c; + if(myverbose!=0) printf(" eig_amax set to %e line %d monomial %d\n", c, line_of_file, current_monomial); + } + {SPC}*EigPolyDeg{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + mnl->eig_polydeg = a; + if(myverbose!=0) printf(" eig_polydeg set to %d line %d monomial %d\n", a, line_of_file, current_monomial); + } + {SPC}*EigNkr{EQL}{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] = %d", name, &a); + mnl->eig_n_kr = a; + if(myverbose!=0) printf(" eig_n_kr set to %d line %d monomial %d\n", a, line_of_file, current_monomial); + } } { {SPC}*MaxPtildeDegree{EQL}{DIGIT}+ {