From 91a45cad9082ed9aa2725f12a657689a5f0465f7 Mon Sep 17 00:00:00 2001 From: Andreas Suter Date: Tue, 27 May 2025 15:59:58 +0200 Subject: [PATCH] added ZF, LF geometry, dynamic approximation width/hopp << 1. --- src/external/LF_GL/main.cpp | 122 ++++++++++++++++++++++++++++++------ 1 file changed, 102 insertions(+), 20 deletions(-) diff --git a/src/external/LF_GL/main.cpp b/src/external/LF_GL/main.cpp index 51a0e068..52fe7c62 100644 --- a/src/external/LF_GL/main.cpp +++ b/src/external/LF_GL/main.cpp @@ -13,10 +13,12 @@ void lf_gl_syntax() { std::cout << std::endl; - std::cout << "usage lg_gl [-p field width hopp [-t tmax] -g [1|0] -o flnOut] | [-h]" << std::endl; + std::cout << "usage lg_gl [[-p field width hopp | -a width hopp] [-t tmax] -g [1|0] -o flnOut] | [-h]" << std::endl; std::cout << " -p: field in (G)" << std::endl; std::cout << " width in (1/us)" << std::endl; std::cout << " hopp in (1/us)" << std::endl; + std::cout << " -a: width in (1/us)" << std::endl; + std::cout << " hopp in (1/us)" << std::endl; std::cout << " -t: tmax in (us). Default: tmax=10 (us)." << std::endl; std::cout << " -g: 1=Gaussian LF data, 0=Gaussian averaged -> Lorentz LF data" << std::endl; std::cout << " -o: flnOut = output file name." << std::endl; @@ -25,11 +27,45 @@ void lf_gl_syntax() } // ---------------------------------------------------------------------------- -int lf_gl_write(const std::string fln, const bool gaussian_only, const std::vector ¶m, +/** + *

ZF Approximation for dynamic (LF geometry). + * Gaussian approximation: Abragam function, A. Yaouanc and P. Dalmas, p. 123, Eq. (6.62) + * Lorentzian approximation: A. Yaouanc and P. Dalmas, p. 130, Eq. (6.89) + * + * @param gaussian if true, Gaussian approximation, otherwise Lorentzian + * @param param [0]: field = 0, [1]: width (1/us), [2]: hopp (1/us) + * @param tmax calculate upto tmax + * @param tt time vector + * @param pol polarization vector + */ +void lf_zf_approx(const bool gaussian, const std::vector ¶m, + const double tmax, std::vector &tt, std::vector &pol) +{ + tt.clear(); + pol.clear(); + double t=0, dval; + do { + tt.push_back(t); + t += 0.001; + if (gaussian) { + dval = exp(-2.0*pow(param[1]/param[2], 2.0)*(exp(-param[2]*t)-1.0+param[2]*t)); + } else { + dval = exp(-sqrt(4.0*pow(param[1]/param[2], 2.0)*(exp(-param[2]*t)-1.0+param[2]*t))); + } + pol.push_back(dval); + } while (t <= tmax); +} + +// ---------------------------------------------------------------------------- +int lf_gl_write(const std::string fln, const bool zf_approx, + const bool gaussian_only, + const std::vector ¶m, const std::vector &tt, const std::vector &pol) { std::ofstream fout(fln, std::ofstream::out); + if (zf_approx) + fout << "# ZF Approximation." << std::endl; if (gaussian_only) fout << "# Gaussian only" << std::endl; else @@ -50,9 +86,10 @@ int main(int argc, char *argv[]) std::vector param{0.0, 0.0, 0.0}; double tmax{10.0}; bool gaussian_only = false; + bool zf_approx = false; std::string flnOut=""; - if (argc < 9) { + if (argc < 8) { lf_gl_syntax(); return 0; } @@ -62,7 +99,7 @@ int main(int argc, char *argv[]) if (!strcmp(argv[i], "-p")) { if (i+3 >= argc) { std::cout << std::endl; - std::cout << "**ERROR** in handling parameters. Not enough input present." << std::endl; + std::cout << "**ERROR** in handling parameters '-p'. Not enough input present." << std::endl; lf_gl_syntax(); return 1; } @@ -112,6 +149,46 @@ int main(int argc, char *argv[]) } param[2] = dval; i += 3; + } else if (!strcmp(argv[i], "-a")) { + if (i+2 >= argc) { + std::cout << std::endl; + std::cout << "**ERROR** in handling parameters for '-a'. Not enough input present." << std::endl; + lf_gl_syntax(); + return 1; + } + param[0] = 0; // ZF approximation + try { + dval = std::stod(argv[i+1], &pos); + } catch (...) { + std::cout << std::endl; + std::cout << "**ERROR** in handling of width." << std::endl; + lf_gl_syntax(); + return 1; + } + if (pos != strlen(argv[i+1])) { + std::cout << std::endl; + std::cout << "**ERROR** in handling of field: '" << argv[i+1] << "'" << std::endl; + lf_gl_syntax(); + return 1; + } + param[1] = dval; + try { + dval = std::stod(argv[i+2], &pos); + } catch (...) { + std::cout << std::endl; + std::cout << "**ERROR** in handling of hopp." << std::endl; + lf_gl_syntax(); + return 1; + } + if (pos != strlen(argv[i+2])) { + std::cout << std::endl; + std::cout << "**ERROR** in handling of field: '" << argv[i+1] << "'" << std::endl; + lf_gl_syntax(); + return 1; + } + param[2] = dval; + i += 2; + zf_approx = true; } else if (!strcmp(argv[i], "-t")) { if (i+1 >= argc) { std::cout << std::endl; @@ -176,29 +253,34 @@ int main(int argc, char *argv[]) std::cout << "Gaussian LF only" << std::endl; else std::cout << "Gaussian averaged -> Lorentz LF" << std::endl; + std::cout << "ZF Approx: " << zf_approx << std::endl; std::cout << "flnOut: " << flnOut << std::endl; std::vector tt; std::vector pol; - if (gaussian_only) { - PGKT_LF gkt_lf(param, tmax); - if (!gkt_lf.IsValid()) { - std::cout << "**ERROR** in Gaussian LF calculation" << std::endl; - return 2; + if (zf_approx) { + lf_zf_approx(gaussian_only, param, tmax, tt, pol); + } else { // full calculation + if (gaussian_only) { + PGKT_LF gkt_lf(param, tmax); + if (!gkt_lf.IsValid()) { + std::cout << "**ERROR** in Gaussian LF calculation" << std::endl; + return 2; + } + tt = gkt_lf.GetTime(); + pol = gkt_lf.GetPol(); + } else { + PLGKT_LF lgkt_lf(param, tmax); + if (!lgkt_lf.IsValid()) { + std::cout << "**ERROR** in Gaussian/Lorentzian LF calculation" << std::endl; + return 2; + } + tt = lgkt_lf.GetTime(); + pol = lgkt_lf.GetPol(); } - tt = gkt_lf.GetTime(); - pol = gkt_lf.GetPol(); - } else { - PLGKT_LF lgkt_lf(param, tmax); - if (!lgkt_lf.IsValid()) { - std::cout << "**ERROR** in Gaussian/Lorentzian LF calculation" << std::endl; - return 2; - } - tt = lgkt_lf.GetTime(); - pol = lgkt_lf.GetPol(); } - lf_gl_write(flnOut, gaussian_only, param, tt, pol); + lf_gl_write(flnOut, zf_approx, gaussian_only, param, tt, pol); return 0; }