Skip to content

Commit eada515

Browse files
committed
add adapt flat hybridization input
1 parent 89a55ad commit eada515

7 files changed

Lines changed: 141 additions & 18 deletions

File tree

test/tools/adapt/CMakeLists.txt

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -7,7 +7,9 @@ file(STRINGS simple_tests simple_tests)
77
foreach(test ${simple_tests})
88
configure_file(${test}/run ${test}/run COPYONLY) # also ensures the test dir is created
99
configure_file(${test}/param ${test}/param COPYONLY)
10-
configure_file(${test}/Delta.dat ${test}/Delta.dat COPYONLY)
10+
if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${test}/Delta.dat)
11+
configure_file(${test}/Delta.dat ${test}/Delta.dat COPYONLY)
12+
endif()
1113
add_test(NAME ${test} COMMAND sh ${CMAKE_CURRENT_SOURCE_DIR}/${test}/run ${PROJECT_BINARY_DIR}/tools ${PROJECT_SOURCE_DIR} ${CMAKE_CURRENT_SOURCE_DIR}/${test}
1214
WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/${test})
1315
set_tests_properties(${test} PROPERTIES ENVIRONMENT "${TEST_BUILD_ENV}")
Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,11 @@
1+
[param]
2+
dos=Delta.dat
3+
Lambda=2
4+
adapt=false
5+
xmax=2.5
6+
outputstep=0.25
7+
xfine=1.5
8+
dx_fine=1e-4
9+
dx_fast=1e-3
10+
allowed_error=1e-8
11+
max_subdiv=20

test/tools/adapt/adapt14_flat/run

Lines changed: 27 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,27 @@
1+
#!/bin/sh
2+
set -eu
3+
export PATH="$1:$PATH"
4+
5+
rm -rf from_file from_flat
6+
mkdir from_file from_flat
7+
cp param from_file/param
8+
cat >from_file/Delta.dat <<'EOF'
9+
-1.0 0.01
10+
-1e-10 0.01
11+
1e-10 0.01
12+
1.0 0.01
13+
EOF
14+
15+
cp param from_flat/param
16+
17+
(cd from_file && adapt P && adapt N)
18+
(cd from_flat && adapt --flat 0.01 P && adapt --flat 0.01 N)
19+
20+
test ! -e from_flat/Delta.dat
21+
"$2/test/mycomp.pl" from_file/FSOL.dat from_flat/FSOL.dat
22+
"$2/test/mycomp.pl" from_file/FSOLNEG.dat from_flat/FSOLNEG.dat
23+
24+
if adapt --flat not-a-number P >invalid_stdout.txt 2>invalid_stderr.txt; then
25+
exit 1
26+
fi
27+
grep -q "positive finite number" invalid_stderr.txt

test/tools/adapt/simple_tests

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,7 @@ adapt10_dxfast
33
adapt11_allowederror
44
adapt12_bandrescale
55
adapt13_invalid_data
6+
adapt14_flat
67
adapt2_xrange
78
adapt3_mesh
89
adapt4_adapt=false

tools/adapt/adapt.cc

Lines changed: 81 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,9 @@
1212
#include <cstdlib>
1313
#include <ctime>
1414
#include <cctype>
15+
#include <cmath>
16+
#include <optional>
17+
#include <stdexcept>
1518

1619
using namespace std::string_literals;
1720

@@ -32,31 +35,95 @@ void help(int argc, char **argv, const std::string &help_message)
3235
}
3336
}
3437

35-
const auto usage = "Usage: adapt [-h] [P|N] [param_filename]"s;
38+
const auto usage = "Usage: adapt [-h|--help] [--flat GG] [P|N] [param_filename]"s;
3639

37-
std::pair<Sign, std::string> cmd_line(int argc, char *argv[]) {
40+
struct CommandLineOptions {
3841
Sign sign = Sign::POS;
39-
if (argc >= 2) {
40-
const auto first_char = std::toupper(static_cast<unsigned char>(argv[1][0]));
41-
switch (first_char) {
42-
case 'P': sign = Sign::POS; break;
43-
case 'N': sign = Sign::NEG; break;
44-
default: std::cerr << usage << std::endl; exit(1);
42+
std::string param_fn = "param";
43+
std::optional<double> flat_gamma;
44+
};
45+
46+
auto uppercase(std::string text) {
47+
for (auto &ch : text) { ch = static_cast<char>(std::toupper(static_cast<unsigned char>(ch))); }
48+
return text;
49+
}
50+
51+
bool parse_sign_arg(const std::string &arg, Sign &sign) {
52+
const auto token = uppercase(arg);
53+
if (token == "P" || token == "POS" || token == "POSITIVE") {
54+
sign = Sign::POS;
55+
return true;
56+
}
57+
if (token == "N" || token == "NEG" || token == "NEGATIVE") {
58+
sign = Sign::NEG;
59+
return true;
60+
}
61+
return false;
62+
}
63+
64+
double parse_flat_gamma(const std::string &value) {
65+
std::size_t parsed = 0;
66+
double gamma = 0.0;
67+
try {
68+
gamma = std::stod(value, &parsed);
69+
} catch (const std::exception &) {
70+
throw std::invalid_argument("--flat expects a positive finite number.");
71+
}
72+
if (parsed != value.size() || !(std::isfinite(gamma) && gamma > 0.0)) {
73+
throw std::invalid_argument("--flat expects a positive finite number.");
74+
}
75+
return gamma;
76+
}
77+
78+
CommandLineOptions cmd_line(int argc, char *argv[]) {
79+
CommandLineOptions options;
80+
bool sign_set = false;
81+
bool param_set = false;
82+
83+
for (int i = 1; i < argc; i++) {
84+
const std::string arg = argv[i];
85+
if (arg == "-h" || arg == "--help") {
86+
std::cout << usage << std::endl;
87+
exit(EXIT_SUCCESS);
88+
}
89+
if (arg == "--flat") {
90+
if (options.flat_gamma) { throw std::invalid_argument("--flat specified more than once.\n" + usage); }
91+
if (i + 1 >= argc) { throw std::invalid_argument("Missing value for --flat.\n" + usage); }
92+
options.flat_gamma = parse_flat_gamma(argv[++i]);
93+
continue;
94+
}
95+
if (arg.starts_with("--flat=")) {
96+
if (options.flat_gamma) { throw std::invalid_argument("--flat specified more than once.\n" + usage); }
97+
options.flat_gamma = parse_flat_gamma(arg.substr(7));
98+
continue;
4599
}
100+
if (!arg.empty() && arg[0] == '-') { throw std::invalid_argument("Unknown option: " + arg + "\n" + usage); }
101+
102+
Sign parsed_sign;
103+
if (parse_sign_arg(arg, parsed_sign)) {
104+
if (sign_set) { throw std::invalid_argument("Sign specified more than once.\n" + usage); }
105+
options.sign = parsed_sign;
106+
sign_set = true;
107+
continue;
108+
}
109+
110+
if (param_set) { throw std::invalid_argument("Unexpected argument: " + arg + "\n" + usage); }
111+
options.param_fn = arg;
112+
param_set = true;
46113
}
47-
const std::string param_fn = argc == 3 ? std::string(argv[2]) : "param";
48-
std::cout << "# ++ " << (sign == Sign::POS ? "POSITIVE" : "NEGATIVE") << std::endl;
49-
return {sign, param_fn};
114+
115+
std::cout << "# ++ " << (options.sign == Sign::POS ? "POSITIVE" : "NEGATIVE") << std::endl;
116+
return options;
50117
}
51118

52119
int main(int argc, char *argv[]) {
53120
try {
54121
const clock_t start_clock = clock();
55122
about();
56123
help(argc, argv, usage);
57-
const auto [sign, param_fn] = cmd_line(argc, argv);
58-
Params P(param_fn);
59-
Adapt calc(P, sign);
124+
const auto options = cmd_line(argc, argv);
125+
Params P(options.param_fn);
126+
Adapt calc(P, options.sign, options.flat_gamma);
60127
calc.run();
61128
const clock_t end_clock = clock();
62129
std::cout << "# Elapsed " << double(end_clock - start_clock) / CLOCKS_PER_SEC << " s" << std::endl;

tools/adapt/adapt.hpp

Lines changed: 9 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -15,6 +15,7 @@
1515
#include <algorithm>
1616
#include <functional>
1717
#include <stdexcept>
18+
#include <optional>
1819

1920
using namespace std;
2021
using namespace std::string_literals;
@@ -56,6 +57,7 @@ class Adapt {
5657
int max_iter = 10; // Maximum number of iterations in the secant method
5758
double max_abs = 100.0; // Maximum value of |f(x)|.
5859
double bandrescale = 1.0; // Rescale the input data by this scale factor
60+
std::optional<double> flat_gamma; // Constant hybridisation supplied on the command line.
5961
bool adapt; // If adapt=false --> g(x)=1.
6062
bool hardgap;
6163
double boundary;
@@ -256,8 +258,12 @@ class Adapt {
256258
return term1 - term2;
257259
}
258260
void load_init_rho() {
259-
std::string rhofn = P.Pstr("dos", "Delta.dat");
260-
vecrho = load_rho(rhofn, sign);
261+
if (flat_gamma) {
262+
vecrho = flat_rho(*flat_gamma, sign);
263+
} else {
264+
std::string rhofn = P.Pstr("dos", "Delta.dat");
265+
vecrho = load_rho(rhofn, sign);
266+
}
261267
add_zero_point(vecrho);
262268
rescalevecxy(vecrho, 1.0/bandrescale, bandrescale);
263269
minmaxvec(vecrho, "rho");
@@ -310,7 +316,7 @@ class Adapt {
310316
factor0 = 1.0 + P.P("secant_factor", 1e-7);
311317
max_iter = P.Pint("secant_max_iter", 10);
312318
}
313-
Adapt(const Params &P_, const Sign &sign_) : P(P_), sign(sign_) {
319+
Adapt(const Params &P_, const Sign &sign_, std::optional<double> flat_gamma_ = std::nullopt) : P(P_), sign(sign_), flat_gamma(flat_gamma_) {
314320
set_parameters();
315321
report_parameters();
316322
}

tools/adapt/load.hpp

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@
1212
#include <stdexcept>
1313
#include <cctype>
1414
#include <cstdlib>
15+
#include <cmath>
1516

1617
#include "../common/tabulated.hpp"
1718

@@ -52,6 +53,14 @@ inline Vec load_rho(const std::string &filename, const Sign sign) {
5253
return vecrho;
5354
}
5455

56+
inline Vec flat_rho(const double gamma, const Sign sign) {
57+
if (!(std::isfinite(gamma) && gamma > 0.0))
58+
throw std::invalid_argument("Flat hybridisation Gamma must be a positive finite number.");
59+
Vec vecrho{{1e-10, gamma}, {1.0, gamma}};
60+
NRG::Tools::print_interval("--flat", sign == Sign::POS ? "POS" : "NEG", vecrho);
61+
return vecrho;
62+
}
63+
5564
inline void save(const std::string &fn, const Vec &v) {
5665
std::ofstream F(fn.c_str());
5766
if (!F)

0 commit comments

Comments
 (0)