Line data Source code
1 : // SPDX-FileCopyrightText: 2024 PairInteraction Developers
2 : // SPDX-License-Identifier: LGPL-3.0-or-later
3 :
4 : #include "pairinteraction/database/Database.hpp"
5 :
6 : #include "pairinteraction/basis/BasisAtom.hpp"
7 : #include "pairinteraction/database/AtomDescriptionByParameters.hpp"
8 : #include "pairinteraction/database/AtomDescriptionByRanges.hpp"
9 : #include "pairinteraction/database/GitHubDownloader.hpp"
10 : #include "pairinteraction/database/ParquetManager.hpp"
11 : #include "pairinteraction/enums/OperatorType.hpp"
12 : #include "pairinteraction/ket/KetAtom.hpp"
13 : #include "pairinteraction/ket/KetNotUniqueError.hpp"
14 : #include "pairinteraction/utils/TaskControl.hpp"
15 : #include "pairinteraction/utils/hash.hpp"
16 : #include "pairinteraction/utils/ket_id.hpp"
17 : #include "pairinteraction/utils/paths.hpp"
18 : #include "pairinteraction/utils/streamed.hpp"
19 :
20 : #include <algorithm>
21 : #include <cmath>
22 : #include <cpptrace/cpptrace.hpp>
23 : #include <duckdb.hpp>
24 : #include <fmt/core.h>
25 : #include <fmt/format.h>
26 : #include <fmt/ranges.h>
27 : #include <fstream>
28 : #include <iterator>
29 : #include <nlohmann/json.hpp>
30 : #include <oneapi/tbb.h>
31 : #include <spdlog/spdlog.h>
32 : #include <string>
33 : #include <system_error>
34 : #include <unordered_map>
35 : #include <unordered_set>
36 : #include <utility>
37 :
38 : namespace pairinteraction {
39 :
40 : namespace {
41 :
42 : // Absolute tolerance by which quantum number ranges are widened so that states are not dropped
43 : // because of the numerical noise of the (expectation) values stored in the database. This noise
44 : // stems from the expectation value calculation as well as from the nullspace solver tolerances
45 : // (MQDT in rydstate).
46 : constexpr double database_numerical_precision = 1e-9;
47 :
48 542 : std::string format_quantum_number_range(const std::string &column, const Range<double> &range) {
49 0 : return fmt::format("{} BETWEEN {} AND {}", column, range.min() - database_numerical_precision,
50 1084 : range.max() + database_numerical_precision);
51 : }
52 :
53 406 : std::string format_expectation_value_range(const std::string &value_column,
54 : const std::string &std_column,
55 : const Range<double> &range,
56 : double standard_deviation_factor) {
57 0 : return fmt::format("{} BETWEEN {}-{}*{}-{} AND {}+{}*{}+{}", value_column, range.min(),
58 : standard_deviation_factor, std_column, database_numerical_precision,
59 406 : range.max(), standard_deviation_factor, std_column,
60 812 : database_numerical_precision);
61 : }
62 :
63 : // Find the index of the result column with the given name.
64 4367 : size_t get_column_index(const std::vector<std::string> &names, const std::string &name) {
65 4367 : auto it = std::find(names.begin(), names.end(), name);
66 4367 : if (it == names.end()) {
67 0 : throw std::runtime_error("Missing database column '" + name + "'.");
68 : }
69 4367 : return static_cast<size_t>(std::distance(names.begin(), it));
70 : }
71 :
72 2335 : bool is_integer_like(duckdb::LogicalTypeId id) {
73 2335 : switch (id) {
74 903 : case duckdb::LogicalTypeId::BOOLEAN:
75 : case duckdb::LogicalTypeId::TINYINT:
76 : case duckdb::LogicalTypeId::SMALLINT:
77 : case duckdb::LogicalTypeId::INTEGER:
78 : case duckdb::LogicalTypeId::BIGINT:
79 : case duckdb::LogicalTypeId::UTINYINT:
80 : case duckdb::LogicalTypeId::USMALLINT:
81 : case duckdb::LogicalTypeId::UINTEGER:
82 : case duckdb::LogicalTypeId::UBIGINT:
83 903 : return true;
84 1432 : default:
85 1432 : return false;
86 : }
87 : }
88 :
89 : template <typename T>
90 1532439 : T get_entry_as(duckdb::Vector &vector, const duckdb::LogicalType &type, size_t row) {
91 1532439 : switch (type.id()) {
92 1388120 : case duckdb::LogicalTypeId::DOUBLE:
93 1388120 : return static_cast<T>(duckdb::FlatVector::GetData<double>(vector)[row]);
94 0 : case duckdb::LogicalTypeId::FLOAT:
95 0 : return static_cast<T>(duckdb::FlatVector::GetData<float>(vector)[row]);
96 0 : case duckdb::LogicalTypeId::BOOLEAN:
97 0 : return static_cast<T>(duckdb::FlatVector::GetData<bool>(vector)[row] ? 1 : 0);
98 0 : case duckdb::LogicalTypeId::TINYINT:
99 0 : return static_cast<T>(duckdb::FlatVector::GetData<int8_t>(vector)[row]);
100 0 : case duckdb::LogicalTypeId::SMALLINT:
101 0 : return static_cast<T>(duckdb::FlatVector::GetData<int16_t>(vector)[row]);
102 0 : case duckdb::LogicalTypeId::INTEGER:
103 0 : return static_cast<T>(duckdb::FlatVector::GetData<int32_t>(vector)[row]);
104 144318 : case duckdb::LogicalTypeId::BIGINT:
105 144318 : return static_cast<T>(duckdb::FlatVector::GetData<int64_t>(vector)[row]);
106 0 : case duckdb::LogicalTypeId::UTINYINT:
107 0 : return static_cast<T>(duckdb::FlatVector::GetData<uint8_t>(vector)[row]);
108 0 : case duckdb::LogicalTypeId::USMALLINT:
109 0 : return static_cast<T>(duckdb::FlatVector::GetData<uint16_t>(vector)[row]);
110 0 : case duckdb::LogicalTypeId::UINTEGER:
111 0 : return static_cast<T>(duckdb::FlatVector::GetData<uint32_t>(vector)[row]);
112 0 : case duckdb::LogicalTypeId::UBIGINT:
113 0 : return static_cast<T>(duckdb::FlatVector::GetData<uint64_t>(vector)[row]);
114 2 : default:
115 2 : throw std::runtime_error("Cannot read database column of type " + type.ToString() +
116 : " as a quantum number.");
117 : }
118 : }
119 :
120 1484596 : double get_entry_as_double(duckdb::Vector &vector, const duckdb::LogicalType &type, size_t row) {
121 1484596 : return get_entry_as<double>(vector, type, row);
122 : }
123 :
124 47842 : int64_t get_entry_as_int64(duckdb::Vector &vector, const duckdb::LogicalType &type, size_t row) {
125 47842 : return get_entry_as<int64_t>(vector, type, row);
126 : }
127 :
128 : struct QuantumNumbers {
129 : std::unordered_map<std::string, double> values;
130 : std::unordered_map<std::string, double> stds;
131 : };
132 :
133 47842 : QuantumNumbers get_quantum_numbers_from_row(duckdb::DataChunk &chunk,
134 : const std::vector<duckdb::LogicalType> &types,
135 : const std::vector<std::string> &names,
136 : const std::unordered_set<std::string> &excluded_columns,
137 : size_t row) {
138 47842 : QuantumNumbers quantum_numbers;
139 1626109 : for (size_t col = 0; col < names.size(); ++col) {
140 1578265 : const std::string &name = names[col];
141 1578265 : if (excluded_columns.contains(name)) {
142 143526 : continue;
143 : }
144 1434735 : double value = get_entry_as_double(chunk.data[col], types[col], row);
145 1434731 : if (name.starts_with("std_")) {
146 574101 : quantum_numbers.stds[name.substr(4)] = value;
147 860629 : } else if (name.starts_with("exp_")) {
148 574101 : quantum_numbers.values[name.substr(4)] = value;
149 : } else {
150 286536 : quantum_numbers.values[name] = value;
151 : }
152 : }
153 47842 : return quantum_numbers;
154 0 : }
155 :
156 2487 : bool must_match_exactly(const std::unordered_map<std::string, duckdb::LogicalTypeId> &columns,
157 : const std::string &column) {
158 2487 : if (column == "f" || column == "m" || column == "parity") {
159 152 : return true;
160 : }
161 2335 : auto it = columns.find(column);
162 2335 : return it != columns.end() && is_integer_like(it->second);
163 : }
164 : } // namespace
165 :
166 47840 : void ensure_consistent_quantum_numbers(double quantum_number_f, double quantum_number_m) {
167 47840 : if (2 * quantum_number_m != std::rint(2 * quantum_number_m)) {
168 0 : throw std::runtime_error("The quantum number m must be an integer or half-integer.");
169 : }
170 47840 : if (2 * quantum_number_f != std::rint(2 * quantum_number_f)) {
171 0 : throw std::runtime_error("The quantum number f must be an integer or half-integer.");
172 : }
173 47840 : if (quantum_number_f + quantum_number_m != std::rint(quantum_number_f + quantum_number_m)) {
174 0 : throw std::invalid_argument(
175 0 : "The quantum numbers f and m must be both either integers or half-integers.");
176 : }
177 47840 : if (std::abs(quantum_number_m) > quantum_number_f) {
178 1 : throw std::invalid_argument(
179 2 : "The absolute value of the quantum number m must be less than or equal to f.");
180 : }
181 47839 : }
182 :
183 0 : Database::Database() : Database(default_download_missing) {}
184 :
185 0 : Database::Database(bool download_missing)
186 0 : : Database(download_missing, default_use_cache, default_database_dir) {}
187 :
188 0 : Database::Database(std::filesystem::path database_dir)
189 0 : : Database(default_download_missing, default_use_cache, std::move(database_dir)) {}
190 :
191 3 : Database::Database(bool download_missing, bool use_cache, std::filesystem::path database_dir)
192 3 : : download_missing_(download_missing), use_cache_(use_cache),
193 3 : database_dir_(std::move(database_dir)), db(std::make_unique<duckdb::DuckDB>(nullptr)),
194 15 : con(std::make_unique<duckdb::Connection>(*db)) {
195 :
196 3 : if (database_dir_.empty()) {
197 0 : database_dir_ = default_database_dir;
198 : }
199 :
200 : // Ensure the database directory exists
201 3 : if (!std::filesystem::exists(database_dir_)) {
202 0 : std::filesystem::create_directories(database_dir_);
203 : }
204 3 : database_dir_ = std::filesystem::canonical(database_dir_);
205 3 : if (!std::filesystem::is_directory(database_dir_)) {
206 0 : throw std::filesystem::filesystem_error("Cannot access database", database_dir_.string(),
207 0 : std::make_error_code(std::errc::not_a_directory));
208 : }
209 3 : SPDLOG_INFO("Using database directory: {}", database_dir_.string());
210 :
211 : // Ensure that the config directory exists
212 3 : std::filesystem::path configdir = paths::get_config_directory();
213 3 : if (!std::filesystem::exists(configdir)) {
214 1 : std::filesystem::create_directories(configdir);
215 2 : } else if (!std::filesystem::is_directory(configdir)) {
216 0 : throw std::filesystem::filesystem_error("Cannot access config directory ",
217 0 : configdir.string(),
218 0 : std::make_error_code(std::errc::not_a_directory));
219 : }
220 :
221 : // Read in the database_repo_paths if a config file exists, otherwise use the default and
222 : // write it to the config file
223 3 : std::filesystem::path configfile = configdir / "database.json";
224 3 : std::string database_repo_host;
225 3 : std::vector<std::string> database_repo_paths;
226 3 : if (std::filesystem::exists(configfile)) {
227 2 : std::ifstream file(configfile);
228 2 : nlohmann::json doc = nlohmann::json::parse(file, nullptr, false);
229 :
230 4 : if (!doc.is_discarded() && doc.contains("hash") && doc.contains("database_repo_host") &&
231 2 : doc.contains("database_repo_paths")) {
232 2 : database_repo_host = doc["database_repo_host"].get<std::string>();
233 2 : database_repo_paths = doc["database_repo_paths"].get<std::vector<std::string>>();
234 :
235 : // If the values are not equal to the default values but the hash is consistent (i.e.,
236 : // the user has not changed anything manually), clear the values so that they can be
237 : // updated
238 4 : if (database_repo_host != default_database_repo_host ||
239 2 : database_repo_paths != default_database_repo_paths) {
240 0 : std::size_t seed = 0;
241 0 : utils::hash_combine(seed, database_repo_paths);
242 0 : utils::hash_combine(seed, database_repo_host);
243 0 : if (seed == doc["hash"].get<std::size_t>()) {
244 0 : database_repo_host.clear();
245 0 : database_repo_paths.clear();
246 : } else {
247 0 : SPDLOG_INFO("The database repository host and paths have been changed "
248 : "manually. Thus, they will not be updated automatically. To reset "
249 : "them, delete the file '{}'.",
250 : configfile.string());
251 : }
252 : }
253 : }
254 2 : }
255 :
256 : // Read in and store the default values if necessary
257 3 : if (database_repo_host.empty() || database_repo_paths.empty()) {
258 2 : SPDLOG_INFO("Updating the database repository host and paths:");
259 :
260 1 : database_repo_host = default_database_repo_host;
261 1 : database_repo_paths = default_database_repo_paths;
262 1 : std::ofstream file(configfile);
263 1 : nlohmann::json doc;
264 :
265 1 : SPDLOG_INFO("* New host: {}", default_database_repo_host);
266 2 : SPDLOG_INFO("* New paths: {}", fmt::join(default_database_repo_paths, ", "));
267 :
268 1 : doc["database_repo_host"] = default_database_repo_host;
269 1 : doc["database_repo_paths"] = database_repo_paths;
270 :
271 1 : std::size_t seed = 0;
272 1 : utils::hash_combine(seed, default_database_repo_paths);
273 1 : utils::hash_combine(seed, default_database_repo_host);
274 1 : doc["hash"] = seed;
275 :
276 1 : file << doc.dump(4);
277 1 : }
278 :
279 : // Limit the memory usage of duckdb's buffer manager
280 : {
281 3 : auto result = con->Query("PRAGMA max_memory = '8GB';");
282 3 : if (result->HasError()) {
283 0 : throw cpptrace::runtime_error("Error setting the memory limit: " + result->GetError());
284 : }
285 3 : }
286 :
287 : // Instantiate a database manager that provides access to database tables. If a table
288 : // is outdated/not available locally, it will be downloaded if download_missing_ is true.
289 3 : if (!download_missing_) {
290 3 : database_repo_paths.clear();
291 : }
292 3 : downloader = std::make_unique<GitHubDownloader>();
293 6 : manager = std::make_unique<ParquetManager>(database_dir_, *downloader, database_repo_paths,
294 6 : *con, use_cache_);
295 3 : manager->scan_local();
296 3 : manager->scan_remote();
297 :
298 : // Print versions of tables
299 3 : std::istringstream iss(manager->get_versions_info());
300 36 : for (std::string line; std::getline(iss, line);) {
301 33 : SPDLOG_INFO(line);
302 3 : }
303 9 : }
304 :
305 2 : Database::~Database() = default;
306 :
307 2187 : const Database::column_info_t &Database::get_states_column_info(const std::string &species) {
308 2187 : std::string table_path = manager->get_path(species, "states");
309 2187 : if (auto it = column_info_cache.find(table_path); it != column_info_cache.end()) {
310 2177 : return it->second;
311 : }
312 20 : auto result = con->Query(fmt::format(R"(SELECT * FROM '{}' LIMIT 0)", table_path));
313 10 : if (result->HasError()) {
314 0 : throw cpptrace::runtime_error("Error querying the database columns: " + result->GetError());
315 : }
316 10 : column_info_t columns;
317 10 : columns.reserve(result->names.size());
318 320 : for (size_t col = 0; col < result->names.size(); ++col) {
319 310 : columns.emplace(result->names[col], result->types[col].id());
320 : }
321 :
322 : // Every states table must provide the columns required for constructing kets. All other
323 : // columns are treated generically as quantum numbers, whatever they are named.
324 50 : for (const auto &required : {"id", "energy", "f", "parity"}) {
325 40 : if (!columns.contains(required)) {
326 0 : throw std::runtime_error(
327 0 : fmt::format("The table of states of species '{}' is missing the required column "
328 : "'{}'.",
329 0 : species, required));
330 : }
331 : }
332 40 : for (const auto &forbidden : {"m", "exp_m", "std_m"}) {
333 30 : if (columns.contains(forbidden)) {
334 0 : throw std::runtime_error(
335 0 : fmt::format("The table of states of species '{}' must not contain a column '{}', "
336 : "since the quantum number m is generated by the software itself.",
337 0 : species, forbidden));
338 : }
339 : }
340 :
341 : // If another thread inserted the same entry concurrently, insert keeps the existing value and
342 : // returns an iterator to it, so the reference stays valid (elements are never erased).
343 10 : return column_info_cache.insert({table_path, std::move(columns)}).first->second;
344 2187 : }
345 :
346 520 : std::shared_ptr<const KetAtom> Database::get_ket(const std::string &species,
347 : const AtomDescriptionByParameters &description) {
348 : // Check that the specifications are valid
349 520 : if (!description.quantum_numbers.contains("m")) {
350 0 : throw std::invalid_argument("The quantum number m must be specified.");
351 : }
352 2581 : for (const auto &[name, value] : description.quantum_numbers) {
353 2061 : if ((name == "f" || name == "m") && 2 * value != std::rint(2 * value)) {
354 0 : throw std::invalid_argument("The quantum number " + name +
355 0 : " must be an integer or half-integer.");
356 : }
357 2061 : if (name == "f" && value < 0) {
358 0 : throw std::invalid_argument("The quantum number " + name + " must be positive.");
359 : }
360 : }
361 :
362 520 : const auto &columns = get_states_column_info(species);
363 :
364 : // Describe the state. Integer-valued quantum numbers and (f, m, parity) are matched exactly,
365 : // while all other quantum numbers are matched within a +-0.5 window (they can deviate from the
366 : // requested value, e.g. expectation values in MQDT). The result is ordered by the distance to
367 : // the requested values, so the nearest state is returned.
368 520 : std::string where;
369 520 : std::string where_separator;
370 520 : std::string orderby;
371 520 : std::string orderby_separator;
372 520 : if (description.energy.has_value()) {
373 : // The following condition derives from demanding that quantum number n that corresponds to
374 : // the energy "E_n = -1/(2*n^2)" is not off by more than 1 from the actual quantum number n,
375 : // i.e., "sqrt(-1/(2*E_n)) - sqrt(-1/(2*E_{n-1})) = 1"
376 0 : double n_from_energy = std::sqrt(-1 / (2 * description.energy.value()));
377 0 : where += where_separator +
378 0 : fmt::format("SQRT(-1/(2*energy)) BETWEEN {} AND {}", n_from_energy - 0.5,
379 0 : n_from_energy + 0.5);
380 0 : where_separator = " AND ";
381 0 : orderby += orderby_separator + fmt::format("(SQRT(-1/(2*energy)) - {})^2", n_from_energy);
382 0 : orderby_separator = " + ";
383 : }
384 2579 : for (const auto &[name, value] : description.quantum_numbers) {
385 2060 : if (name == "m") {
386 520 : continue; // m is encoded into the id, not stored as a queryable column
387 : }
388 1540 : std::string column = columns.contains("exp_" + name) ? "exp_" + name : name;
389 1540 : if (!columns.contains(column)) {
390 2 : throw std::invalid_argument(
391 2 : fmt::format("The quantum number '{}' is not stored in the database table for "
392 : "species '{}'.",
393 2 : name, species));
394 : }
395 1539 : double tolerance = must_match_exactly(columns, column) ? 0.0 : 0.5;
396 1539 : where += where_separator +
397 4617 : fmt::format("{} BETWEEN {} AND {}", column, value - tolerance, value + tolerance);
398 1539 : where_separator = " AND ";
399 3078 : orderby += orderby_separator + fmt::format("({} - {})^2", column, value);
400 1539 : orderby_separator = " + ";
401 1540 : }
402 519 : if (where_separator.empty()) {
403 0 : where += "FALSE";
404 : }
405 519 : if (orderby_separator.empty()) {
406 0 : orderby += "id";
407 : }
408 :
409 : // Ask the database for the described state
410 519 : set_task_status("Loading atomic ket from database...");
411 1038 : auto result = con->Query(fmt::format(
412 : R"(SELECT *, {} AS order_val FROM '{}' WHERE {} ORDER BY order_val ASC LIMIT 2)", orderby,
413 1557 : manager->get_path(species, "states"), where));
414 :
415 519 : if (result->HasError()) {
416 0 : throw cpptrace::runtime_error("Error querying the database: " + result->GetError());
417 : }
418 :
419 519 : if (result->RowCount() == 0) {
420 5 : throw std::invalid_argument("No state found.");
421 : }
422 :
423 : // Get the first chunk of the results (the first chunk is sufficient as we need two rows at
424 : // most). Every column except energy, id and the synthetic order_val is treated as a quantum
425 : // number; m is not a database column and is injected from the description.
426 514 : const auto &types = result->types;
427 514 : const auto &names = result->names;
428 2056 : const std::unordered_set<std::string> excluded_columns = {"energy", "id", "order_val"};
429 514 : auto chunk = result->Fetch();
430 :
431 514 : size_t energy_column = get_column_index(names, "energy");
432 514 : size_t id_column = get_column_index(names, "id");
433 514 : double quantum_number_m = description.quantum_numbers.at("m");
434 :
435 515 : auto make_ket = [&](size_t row) {
436 : auto quantum_numbers =
437 515 : get_quantum_numbers_from_row(*chunk, types, names, excluded_columns, row);
438 515 : quantum_numbers.values["m"] = quantum_number_m;
439 515 : double energy = get_entry_as_double(chunk->data[energy_column], types[energy_column], row);
440 1030 : auto id = utils::encode_as_ket_id({.id = static_cast<size_t>(get_entry_as_int64(
441 515 : chunk->data[id_column], types[id_column], row)),
442 515 : .m = quantum_number_m});
443 515 : return KetAtom(typename KetAtom::Private(), energy, species,
444 1030 : std::move(quantum_numbers.values), std::move(quantum_numbers.stds), *this,
445 1545 : id);
446 515 : };
447 :
448 : // Check that the ket is uniquely specified. If not, throw a KetNotUniqueError.
449 514 : if (chunk->size() > 1) {
450 7 : size_t order_val_column = get_column_index(names, "order_val");
451 : auto order_val_0 =
452 7 : get_entry_as_double(chunk->data[order_val_column], types[order_val_column], 0);
453 : auto order_val_1 =
454 7 : get_entry_as_double(chunk->data[order_val_column], types[order_val_column], 1);
455 :
456 7 : if (order_val_1 - order_val_0 <= order_val_0) {
457 4 : throw KetNotUniqueError({
458 0 : std::make_shared<const KetAtom>(make_ket(0)),
459 1 : std::make_shared<const KetAtom>(make_ket(1)),
460 5 : });
461 : }
462 : }
463 :
464 : // Construct the state
465 513 : auto ket = make_ket(0);
466 :
467 : // Check database consistency
468 515 : ensure_consistent_quantum_numbers(ket.get_quantum_number("f"), quantum_number_m);
469 :
470 1024 : return std::make_shared<const KetAtom>(std::move(ket));
471 1585 : }
472 :
473 : template <typename Scalar>
474 : std::shared_ptr<const BasisAtom<Scalar>>
475 1667 : Database::get_basis(const std::string &species, const AtomDescriptionByRanges &description,
476 : const std::vector<size_t> &additional_ket_ids) {
477 : // The quantum number m is restricted separately because it is generated by UNNEST below.
478 3334 : auto range_quantum_number_m = [&description]() {
479 1667 : auto it = description.quantum_number_ranges.find("m");
480 3334 : return it != description.quantum_number_ranges.end() ? it->second : Range<double>{};
481 1667 : }();
482 :
483 1667 : const auto &columns = get_states_column_info(species);
484 :
485 : // Describe the states by all restrictions that do not involve the quantum number m
486 1667 : std::string where = "(";
487 1667 : std::string separator;
488 1667 : if (description.range_energy.is_finite()) {
489 57 : where += separator +
490 57 : fmt::format("energy BETWEEN {} AND {}", description.range_energy.min(),
491 57 : description.range_energy.max());
492 57 : separator = " AND ";
493 : }
494 3446 : for (const auto &[name, range] : description.quantum_number_ranges) {
495 949 : if (name == "m" || !range.is_finite()) {
496 119 : continue;
497 : }
498 830 : std::string exp_column = "exp_" + name;
499 830 : std::string std_column = "std_" + name;
500 830 : if (columns.contains(exp_column) && columns.contains(std_column)) {
501 406 : where += separator +
502 : format_expectation_value_range(
503 : exp_column, std_column, range,
504 406 : description.quantum_number_standard_deviation_factor);
505 424 : } else if (columns.contains(name)) {
506 423 : where += separator + format_quantum_number_range(name, range);
507 : } else {
508 2 : throw std::invalid_argument(
509 : fmt::format("The quantum number '{}' is not stored in the database table for "
510 : "species '{}'.",
511 : name, species));
512 : }
513 829 : separator = " AND ";
514 : }
515 1666 : if (separator.empty()) {
516 : // If the description contains no restrictions at all, it describes no states
517 1247 : where += range_quantum_number_m.is_finite() ? "TRUE" : "FALSE";
518 : }
519 1666 : where += ")";
520 :
521 : // Describe the restriction of the quantum number m
522 1666 : std::string where_m = "(";
523 1666 : if (range_quantum_number_m.is_finite()) {
524 119 : where_m += format_quantum_number_range("m", range_quantum_number_m);
525 : } else {
526 1547 : where_m += "TRUE";
527 : }
528 1666 : where_m += ")";
529 :
530 : // Create a table containing the described states
531 1666 : std::string canonical_basis_id;
532 : {
533 1666 : auto result = con->Query(R"(SELECT UUID()::varchar)");
534 1666 : if (result->HasError()) {
535 0 : throw cpptrace::runtime_error("Error selecting canonical_basis_id: " +
536 0 : result->GetError());
537 : }
538 1666 : canonical_basis_id =
539 3332 : duckdb::FlatVector::GetData<duckdb::string_t>(result->Fetch()->data[0])[0].GetString();
540 1666 : }
541 : {
542 1666 : set_task_status("Selecting atomic basis states...");
543 4998 : auto result = con->Query(fmt::format(
544 : R"(CREATE TEMP TABLE '{}' AS SELECT *, id*{}+(2*m+{})::bigint AS ketid FROM (
545 : SELECT *,
546 : UNNEST(list_transform(generate_series(0,(2*f)::bigint),
547 : x -> x::double-f)) AS m FROM (
548 : SELECT * FROM '{}' WHERE {}
549 : )
550 : ) WHERE {})",
551 : canonical_basis_id, utils::KET_ID_STRIDE, utils::M_OFFSET,
552 3332 : manager->get_path(species, "states"), where, where_m));
553 :
554 1666 : if (result->HasError()) {
555 0 : throw cpptrace::runtime_error("Error creating table: " + result->GetError());
556 : }
557 1666 : }
558 : {
559 3332 : auto result = con->Query(
560 : fmt::format(R"(ALTER TABLE '{}' ADD PRIMARY KEY (ketid))", canonical_basis_id));
561 :
562 1666 : if (result->HasError()) {
563 0 : throw cpptrace::runtime_error("Error adding primary key: " + result->GetError());
564 : }
565 1666 : }
566 :
567 : // Add the additional kets to the table if they are not already contained in it
568 1666 : if (!additional_ket_ids.empty()) {
569 1321 : set_task_status("Selecting additional kets...");
570 1321 : std::vector<std::string> values;
571 1321 : values.reserve(additional_ket_ids.size());
572 3527 : for (size_t ket_id : additional_ket_ids) {
573 2206 : auto [id, m] = utils::decode_from_ket_id(ket_id);
574 4412 : values.push_back(fmt::format("({}, {}::double, {})", id, m, ket_id));
575 : }
576 5284 : auto result = con->Query(fmt::format(
577 : R"(INSERT OR IGNORE INTO '{0}'
578 : SELECT s.*, v.m, v.ketid FROM '{1}' AS s
579 : JOIN (VALUES {2}) AS v(id, m, ketid) ON s.id = v.id)",
580 3963 : canonical_basis_id, manager->get_path(species, "states"), fmt::join(values, ",")));
581 :
582 1321 : if (result->HasError()) {
583 0 : throw cpptrace::runtime_error("Error adding additional kets: " + result->GetError());
584 : }
585 1321 : }
586 :
587 : // Ask the table for the extreme values of the quantum numbers
588 : {
589 1666 : set_task_status("Validating atomic basis coverage...");
590 :
591 : // Collect the finite quantum-number ranges to validate against the loaded basis. Energy is
592 : // handled separately because it is compared via its corresponding effective quantum number.
593 : struct CoverageCheck {
594 : std::string name;
595 : std::string column;
596 : double tolerance{}; // allowed slack between the requested and the available range
597 : Range<double> range;
598 : };
599 1666 : std::vector<CoverageCheck> checks;
600 2614 : for (const auto &[name, range] : description.quantum_number_ranges) {
601 948 : if (!range.is_finite()) {
602 0 : continue;
603 : }
604 948 : bool is_expectation_value =
605 948 : columns.contains("exp_" + name) && columns.contains("std_" + name);
606 948 : std::string column = is_expectation_value ? "exp_" + name : name;
607 948 : double tolerance = must_match_exactly(columns, column) ? 0.0 : 1.0;
608 948 : checks.push_back({name, column, tolerance, range});
609 : }
610 :
611 1666 : std::string select;
612 1666 : std::string separator;
613 1666 : if (description.range_energy.is_finite()) {
614 57 : select += separator + "MIN(energy) AS min_energy, MAX(energy) AS max_energy";
615 57 : separator = ", ";
616 : }
617 2614 : for (const auto &check : checks) {
618 948 : select += separator +
619 948 : fmt::format("MIN({0}) AS min_{1}, MAX({0}) AS max_{1}", check.column, check.name);
620 948 : separator = ", ";
621 : }
622 :
623 1666 : if (!separator.empty()) {
624 838 : auto result =
625 : con->Query(fmt::format(R"(SELECT {} FROM '{}')", select, canonical_basis_id));
626 :
627 419 : if (result->HasError()) {
628 0 : throw cpptrace::runtime_error("Error querying the database: " + result->GetError());
629 : }
630 :
631 419 : auto chunk = result->Fetch();
632 419 : const auto &types = result->types;
633 :
634 2429 : for (size_t i = 0; i < chunk->ColumnCount(); i++) {
635 2010 : if (duckdb::FlatVector::IsNull(chunk->data[i], 0)) {
636 0 : throw std::invalid_argument("No state found.");
637 : }
638 : }
639 :
640 419 : size_t idx = 0;
641 419 : if (description.range_energy.is_finite()) {
642 57 : auto min_energy = get_entry_as_double(chunk->data[idx], types[idx], 0);
643 57 : idx++;
644 114 : if (std::sqrt(-1 / (2 * min_energy)) - 1 >
645 57 : std::sqrt(-1 / (2 * description.range_energy.min()))) {
646 0 : SPDLOG_DEBUG("No state found with the requested minimum energy. Requested: {}, "
647 : "found: {}.",
648 : description.range_energy.min(), min_energy);
649 : }
650 57 : auto max_energy = get_entry_as_double(chunk->data[idx], types[idx], 0);
651 57 : idx++;
652 114 : if (std::sqrt(-1 / (2 * max_energy)) + 1 <
653 57 : std::sqrt(-1 / (2 * description.range_energy.max()))) {
654 0 : SPDLOG_DEBUG("No state found with the requested maximum energy. Requested: {}, "
655 : "found: {}.",
656 : description.range_energy.max(), max_energy);
657 : }
658 : }
659 1367 : for (const auto &check : checks) {
660 948 : auto min_value = get_entry_as_double(chunk->data[idx], types[idx], 0);
661 948 : idx++;
662 948 : if (min_value - check.tolerance > check.range.min()) {
663 144 : SPDLOG_DEBUG("No state found with the requested minimum quantum number {}. "
664 : "Requested: {}, found: {}.",
665 : check.name, check.range.min(), min_value);
666 : }
667 948 : auto max_value = get_entry_as_double(chunk->data[idx], types[idx], 0);
668 948 : idx++;
669 948 : if (max_value + check.tolerance < check.range.max()) {
670 52 : SPDLOG_DEBUG("No state found with the requested maximum quantum number {}. "
671 : "Requested: {}, found: {}.",
672 : check.name, check.range.max(), max_value);
673 : }
674 : }
675 419 : }
676 1666 : }
677 :
678 : // Ask the table for the described states
679 1666 : set_task_status("Loading atomic basis states...");
680 3332 : auto result =
681 : con->Query(fmt::format(R"(SELECT * FROM '{}' ORDER BY ketid ASC)", canonical_basis_id));
682 :
683 1666 : if (result->HasError()) {
684 0 : throw cpptrace::runtime_error("Error querying the database: " + result->GetError());
685 : }
686 :
687 1666 : if (result->RowCount() == 0) {
688 0 : throw std::invalid_argument("No state found.");
689 : }
690 :
691 : // Construct the states. Every column except energy and the raw/encoded id is treated as a
692 : // quantum number ("id" is the raw states-table id, "ketid" the m-encoded id of the basis
693 : // state).
694 1666 : const auto &types = result->types;
695 1666 : const auto &names = result->names;
696 6664 : const std::unordered_set<std::string> excluded_columns = {"energy", "id", "ketid"};
697 1666 : size_t energy_column = get_column_index(names, "energy");
698 1666 : size_t ketid_column = get_column_index(names, "ketid");
699 :
700 1666 : std::vector<std::shared_ptr<const KetAtom>> kets;
701 1666 : kets.reserve(result->RowCount());
702 1666 : double last_energy = std::numeric_limits<double>::lowest();
703 :
704 3332 : for (auto chunk = result->Fetch(); chunk; chunk = result->Fetch()) {
705 1666 : set_task_status("Constructing atomic basis...");
706 :
707 48993 : for (size_t i = 0; i < chunk->size(); i++) {
708 94654 : auto quantum_numbers =
709 47327 : get_quantum_numbers_from_row(*chunk, types, names, excluded_columns, i);
710 47327 : double energy =
711 47327 : get_entry_as_double(chunk->data[energy_column], types[energy_column], i);
712 47327 : auto id = static_cast<size_t>(
713 47327 : get_entry_as_int64(chunk->data[ketid_column], types[ketid_column], i));
714 :
715 : // Check database consistency
716 47327 : ensure_consistent_quantum_numbers(quantum_numbers.values.at("f"),
717 94654 : quantum_numbers.values.at("m"));
718 47327 : if (energy < last_energy) {
719 0 : throw std::runtime_error("The states are not sorted by energy.");
720 : }
721 47327 : last_energy = energy;
722 :
723 : // Append a new state
724 47327 : kets.push_back(std::make_shared<const KetAtom>(
725 94654 : typename KetAtom::Private(), energy, species, std::move(quantum_numbers.values),
726 47327 : std::move(quantum_numbers.stds), *this, id));
727 : }
728 : }
729 :
730 0 : return std::make_shared<const BasisAtom<Scalar>>(typename BasisAtom<Scalar>::Private(),
731 1666 : std::move(kets), std::move(canonical_basis_id),
732 3332 : *this);
733 5948 : }
734 :
735 : template <typename Scalar>
736 5888 : Eigen::SparseMatrix<Scalar, Eigen::RowMajor> Database::get_matrix_elements_in_canonical_basis(
737 : std::shared_ptr<const BasisAtom<Scalar>> initial_basis,
738 : std::shared_ptr<const BasisAtom<Scalar>> final_basis, OperatorType type, int q) {
739 : using real_t = typename traits::NumTraits<Scalar>::real_t;
740 : using cached_matrix_ptr_t = std::shared_ptr<const cached_matrix_t>;
741 :
742 5888 : if (&initial_basis->get_database() != this || &final_basis->get_database() != this) {
743 1 : throw std::invalid_argument(
744 : "The initial and final bases must belong to the Database instance used for the "
745 : "matrix element calculation.");
746 : }
747 5916 : if (initial_basis->get_species() != final_basis->get_species()) {
748 0 : throw std::invalid_argument(fmt::format(
749 : "The initial and final bases must have the same species, but got '{}' and '{}'.",
750 0 : initial_basis->get_species(), final_basis->get_species()));
751 : }
752 :
753 5909 : std::string specifier;
754 5971 : int kappa{};
755 5971 : switch (type) {
756 129 : case OperatorType::ELECTRIC_MONOPOLE:
757 129 : specifier = "matrix_elements_m";
758 129 : kappa = 0;
759 129 : break;
760 4023 : case OperatorType::ELECTRIC_DIPOLE:
761 4023 : specifier = "matrix_elements_d";
762 4062 : kappa = 1;
763 4062 : break;
764 499 : case OperatorType::ELECTRIC_QUADRUPOLE:
765 499 : specifier = "matrix_elements_q";
766 499 : kappa = 2;
767 499 : break;
768 95 : case OperatorType::ELECTRIC_QUADRUPOLE_ZERO:
769 95 : specifier = "matrix_elements_q0";
770 94 : kappa = 0;
771 94 : break;
772 0 : case OperatorType::ELECTRIC_OCTUPOLE:
773 0 : specifier = "matrix_elements_o";
774 0 : kappa = 3;
775 0 : break;
776 56 : case OperatorType::MAGNETIC_DIPOLE:
777 56 : specifier = "matrix_elements_mu";
778 57 : kappa = 1;
779 57 : break;
780 4 : case OperatorType::ENERGY:
781 4 : specifier = "energy";
782 4 : kappa = 0;
783 4 : break;
784 1137 : case OperatorType::IDENTITY:
785 1137 : specifier = "identity";
786 1137 : kappa = 0;
787 1137 : break;
788 28 : default:
789 28 : throw std::invalid_argument("Unknown operator type.");
790 : }
791 :
792 5982 : std::string canonical_basis_id_initial = initial_basis->get_canonical_basis_id();
793 5952 : std::string canonical_basis_id_final = final_basis->get_canonical_basis_id();
794 5968 : std::string cache_key = fmt::format("{}_{}_{}_{}", specifier, q, canonical_basis_id_initial,
795 : canonical_basis_id_final);
796 5973 : auto &matrix_elements_cache = get_matrix_elements_cache();
797 5885 : std::promise<cached_matrix_ptr_t> matrix_promise;
798 5907 : auto [cache_it, inserted] =
799 : matrix_elements_cache.insert({cache_key, matrix_promise.get_future().share()});
800 :
801 5884 : if (inserted) {
802 : try {
803 1775 : Eigen::Index num_rows = final_basis->get_number_of_kets();
804 1774 : Eigen::Index num_cols = initial_basis->get_number_of_kets();
805 :
806 1774 : std::vector<int> outerIndexPtr;
807 1774 : std::vector<int> innerIndices;
808 1774 : std::vector<real_t> values;
809 :
810 : // Check that the specifications are valid
811 1774 : if (std::abs(q) > kappa) {
812 0 : throw std::invalid_argument("Invalid q.");
813 : }
814 :
815 : // Ask the database for the operator
816 1774 : set_task_status("Loading matrix elements from database...");
817 1777 : std::string species = initial_basis->get_species();
818 1777 : duckdb::unique_ptr<duckdb::MaterializedQueryResult> result;
819 1777 : if (specifier == "identity") {
820 1972 : result = con->Query(fmt::format(
821 : R"(SELECT s2.ketid AS row, s1.ketid AS col, 1.0::DOUBLE AS val
822 : FROM '{}' AS s1
823 : INNER JOIN '{}' AS s2 ON s1.ketid = s2.ketid
824 : ORDER BY row ASC)",
825 : canonical_basis_id_initial, canonical_basis_id_final));
826 790 : } else if (specifier == "energy") {
827 6 : result = con->Query(fmt::format(
828 : R"(SELECT s2.ketid AS row, s1.ketid AS col, s1.energy AS val
829 : FROM '{}' AS s1
830 : INNER JOIN '{}' AS s2 ON s1.ketid = s2.ketid
831 : ORDER BY row ASC)",
832 : canonical_basis_id_initial, canonical_basis_id_final));
833 : } else {
834 3151 : result = con->Query(fmt::format(
835 : R"(WITH s1 AS (
836 : SELECT id, f, m, ketid FROM '{}'
837 : ),
838 : s2 AS (
839 : SELECT id, f, m, ketid FROM '{}'
840 : ),
841 : b AS (
842 : SELECT MIN(f) AS min_f, MAX(f) AS max_f,
843 : MIN(id) AS min_id, MAX(id) AS max_id
844 : FROM (SELECT f, id FROM s1 UNION ALL SELECT f, id FROM s2)
845 : ),
846 : w_filtered AS (
847 : SELECT *
848 : FROM '{}'
849 : WHERE kappa = {} AND q = {} AND
850 : f_initial BETWEEN (SELECT min_f FROM b) AND (SELECT max_f FROM b) AND
851 : f_final BETWEEN (SELECT min_f FROM b) AND (SELECT max_f FROM b)
852 : ),
853 : e_filtered AS (
854 : SELECT *
855 : FROM '{}'
856 : WHERE
857 : id_initial BETWEEN (SELECT min_id FROM b) AND (SELECT max_id FROM b) AND
858 : id_final BETWEEN (SELECT min_id FROM b) AND (SELECT max_id FROM b)
859 : )
860 : SELECT
861 : s2.ketid AS row,
862 : s1.ketid AS col,
863 : e.val*w.val AS val
864 : FROM e_filtered AS e
865 : JOIN s1 ON e.id_initial = s1.id
866 : JOIN s2 ON e.id_final = s2.id
867 : JOIN w_filtered AS w ON
868 : w.f_initial = s1.f AND w.m_initial = s1.m AND
869 : w.f_final = s2.f AND w.m_final = s2.m
870 : ORDER BY row ASC, col ASC)",
871 : canonical_basis_id_initial, canonical_basis_id_final,
872 1575 : manager->get_path("misc", "wigner"), kappa, q,
873 : manager->get_path(species, specifier)));
874 : }
875 :
876 1777 : if (result->HasError()) {
877 0 : throw cpptrace::runtime_error("Error querying the database: " + result->GetError());
878 : }
879 :
880 : // Check the types of the columns
881 1777 : const auto &types = result->types;
882 1777 : const auto &labels = result->names;
883 7108 : const std::vector<duckdb::LogicalType> ref_types = {duckdb::LogicalType::BIGINT,
884 : duckdb::LogicalType::BIGINT,
885 : duckdb::LogicalType::DOUBLE};
886 7108 : for (size_t i = 0; i < types.size(); i++) {
887 5331 : if (types[i] != ref_types[i]) {
888 0 : throw std::runtime_error("Wrong type for '" + labels[i] + "'.");
889 : }
890 : }
891 :
892 1777 : set_task_status("Constructing matrix elements...");
893 :
894 : // Construct the matrix
895 1777 : int num_entries = static_cast<int>(result->RowCount());
896 1777 : outerIndexPtr.reserve(num_rows + 1);
897 1777 : innerIndices.reserve(num_entries);
898 1777 : values.reserve(num_entries);
899 :
900 1777 : int last_row = -1;
901 :
902 3723 : for (auto chunk = result->Fetch(); chunk; chunk = result->Fetch()) {
903 1946 : auto *chunk_row = duckdb::FlatVector::GetData<int64_t>(chunk->data[0]);
904 1946 : auto *chunk_col = duckdb::FlatVector::GetData<int64_t>(chunk->data[1]);
905 1946 : auto *chunk_val = duckdb::FlatVector::GetData<double>(chunk->data[2]);
906 :
907 878163 : for (size_t i = 0; i < chunk->size(); i++) {
908 876217 : int row = final_basis->get_ket_index_from_id(chunk_row[i]);
909 876217 : if (row != last_row) {
910 58970 : if (row < last_row) {
911 0 : throw std::runtime_error("The rows are not sorted.");
912 : }
913 124143 : for (; last_row < row; last_row++) {
914 65173 : outerIndexPtr.push_back(static_cast<int>(innerIndices.size()));
915 : }
916 : }
917 876217 : innerIndices.push_back(initial_basis->get_ket_index_from_id(chunk_col[i]));
918 876217 : values.push_back(chunk_val[i]);
919 : }
920 : }
921 :
922 6685 : for (; last_row < num_rows + 1; last_row++) {
923 4908 : outerIndexPtr.push_back(static_cast<int>(innerIndices.size()));
924 : }
925 :
926 1777 : Eigen::Map<const cached_matrix_t> matrix_map(num_rows, num_cols, values.size(),
927 1777 : outerIndexPtr.data(), innerIndices.data(),
928 1777 : values.data());
929 :
930 1777 : auto cached_matrix = std::make_shared<const cached_matrix_t>(matrix_map);
931 1777 : matrix_promise.set_value(std::move(cached_matrix));
932 1777 : } catch (...) {
933 0 : matrix_promise.set_exception(std::current_exception());
934 0 : throw;
935 : }
936 : }
937 :
938 5886 : set_task_status("Returning matrix elements in canonical basis...");
939 :
940 11995 : return cache_it->second.get()->template cast<Scalar>();
941 7747 : }
942 :
943 1266 : int Database::get_minimal_kappa(const std::string &species) {
944 : // The monopole operator is only stored for charged species. Thus, if the table exists, the
945 : // species is charged and the monopole operator has to be taken into account.
946 1266 : return manager->has_table(species, "matrix_elements_m") ? 0 : 1;
947 : }
948 :
949 4 : bool Database::get_download_missing() const { return download_missing_; }
950 :
951 2 : bool Database::get_use_cache() const { return use_cache_; }
952 :
953 10 : std::filesystem::path Database::get_database_dir() const { return database_dir_; }
954 :
955 0 : std::string Database::get_versions_info() const { return manager->get_versions_info(); }
956 :
957 5985 : Database::matrix_elements_cache_t &Database::get_matrix_elements_cache() {
958 5985 : static matrix_elements_cache_t matrix_elements_cache;
959 5946 : return matrix_elements_cache;
960 : }
961 :
962 63 : Database &Database::get_global_instance() {
963 126 : return get_global_instance_without_checks(default_download_missing, default_use_cache,
964 126 : default_database_dir);
965 : }
966 :
967 0 : Database &Database::get_global_instance(bool download_missing) {
968 0 : Database &database = get_global_instance_without_checks(download_missing, default_use_cache,
969 : default_database_dir);
970 0 : if (download_missing != database.download_missing_) {
971 0 : throw std::invalid_argument(
972 0 : "The 'download_missing' argument must not change between calls to the method.");
973 : }
974 0 : return database;
975 : }
976 :
977 0 : Database &Database::get_global_instance(std::filesystem::path database_dir) {
978 0 : if (database_dir.empty()) {
979 0 : database_dir = default_database_dir;
980 : }
981 0 : Database &database = get_global_instance_without_checks(default_download_missing,
982 : default_use_cache, database_dir);
983 0 : if (!std::filesystem::exists(database_dir) ||
984 0 : std::filesystem::canonical(database_dir) != database.database_dir_) {
985 0 : throw std::invalid_argument(
986 0 : "The 'database_dir' argument must not change between calls to the method.");
987 : }
988 0 : return database;
989 : }
990 :
991 1 : Database &Database::get_global_instance(bool download_missing, bool use_cache,
992 : std::filesystem::path database_dir) {
993 1 : if (database_dir.empty()) {
994 0 : database_dir = default_database_dir;
995 : }
996 : Database &database =
997 1 : get_global_instance_without_checks(download_missing, use_cache, database_dir);
998 1 : if (download_missing != database.download_missing_ || use_cache != database.use_cache_ ||
999 3 : !std::filesystem::exists(database_dir) ||
1000 2 : std::filesystem::canonical(database_dir) != database.database_dir_) {
1001 0 : throw std::invalid_argument(
1002 : "The 'download_missing', 'use_cache' and 'database_dir' arguments must not "
1003 0 : "change between calls to the method.");
1004 : }
1005 1 : return database;
1006 : }
1007 :
1008 64 : Database &Database::get_global_instance_without_checks(bool download_missing, bool use_cache,
1009 : std::filesystem::path database_dir) {
1010 64 : static Database database(download_missing, use_cache, std::move(database_dir));
1011 64 : return database;
1012 : }
1013 :
1014 : struct database_dir_noexcept : std::filesystem::path {
1015 2 : explicit database_dir_noexcept() noexcept try : std
1016 2 : ::filesystem::path(paths::get_cache_directory() / "database") {}
1017 0 : catch (...) {
1018 0 : SPDLOG_ERROR("Error getting the PairInteraction cache directory.");
1019 0 : std::terminate();
1020 2 : }
1021 : };
1022 :
1023 : const std::filesystem::path Database::default_database_dir = database_dir_noexcept();
1024 :
1025 : // Explicit instantiations
1026 : // NOLINTBEGIN(bugprone-macro-parentheses, cppcoreguidelines-macro-usage)
1027 : #define INSTANTIATE_GETTERS(SCALAR) \
1028 : template std::shared_ptr<const BasisAtom<SCALAR>> Database::get_basis<SCALAR>( \
1029 : const std::string &species, const AtomDescriptionByRanges &description, \
1030 : const std::vector<size_t> &additional_ket_ids); \
1031 : template Eigen::SparseMatrix<SCALAR, Eigen::RowMajor> \
1032 : Database::get_matrix_elements_in_canonical_basis<SCALAR>( \
1033 : std::shared_ptr<const BasisAtom<SCALAR>> initial_basis, \
1034 : std::shared_ptr<const BasisAtom<SCALAR>> final_basis, OperatorType type, int q);
1035 : // NOLINTEND(bugprone-macro-parentheses, cppcoreguidelines-macro-usage)
1036 :
1037 : INSTANTIATE_GETTERS(double)
1038 : INSTANTIATE_GETTERS(std::complex<double>)
1039 :
1040 : #undef INSTANTIATE_GETTERS
1041 : } // namespace pairinteraction
|