From 77dddafc7ac710843a0d1568d2d481485eafe89c Mon Sep 17 00:00:00 2001 From: mfcmquintela Date: Fri, 13 Feb 2026 10:16:57 +0100 Subject: [PATCH] Energy cutoff for exporting BSE eigensystem Implements a CLI flag "-t" or "--encut" that enhances the previous --states flag for setting a number of states to export by exporting the eigensystem up to an energy set via this flag. Energy cutoff takes priority over number of states cutoff, making it easier to export the states needed to study optical response up to a specific frequency. --- include/xatu/Result.hpp | 6 +++--- include/xatu/Result.tpp | 38 ++++++++++++++++++++++++++++---------- include/xatu/utils.hpp | 12 +++++++++--- main/xatu.cpp | 10 ++++++---- 4 files changed, 46 insertions(+), 20 deletions(-) diff --git a/include/xatu/Result.hpp b/include/xatu/Result.hpp index 9161cb7..4a5a11b 100755 --- a/include/xatu/Result.hpp +++ b/include/xatu/Result.hpp @@ -53,9 +53,9 @@ class Result { void writePhase(int, FILE*); void writeExtendedPhase(const arma::cx_vec&, FILE*); void writeExtendedPhase(int, FILE*); - void writeEigenvalues(FILE*, int n = 0); - void writeStates(FILE*, int n = 0); - void writeSpin(int, FILE*); + void writeEigenvalues(FILE*, int n = 0, double encut = 0.0); + void writeStates(FILE*, int n = 0, double encut = 0.0); + void writeSpin(int, double, FILE*); void writeRealspaceAmplitude(int, int, const arma::rowvec&, FILE*, int); virtual void writeRealspaceAmplitude(const arma::cx_vec&, int, const arma::rowvec&, FILE*, int) = 0; virtual void writeAbsorptionSpectrum() = 0; diff --git a/include/xatu/Result.tpp b/include/xatu/Result.tpp index 09f9939..e4e19e6 100755 --- a/include/xatu/Result.tpp +++ b/include/xatu/Result.tpp @@ -401,9 +401,15 @@ void Result::writeRealspaceAmplitude(int stateindex, int holeIndex, * @return void */ template -void Result::writeEigenvalues(FILE* textfile, int n){ +void Result::writeEigenvalues(FILE* textfile, int n, double encut){ - if(n > exciton->excitonbasisdim || n < 0){ + int newn = n; + + if (encut != 0.0){ + newn = (int) arma::abs(eigval - encut).index_min() + 1; + } + + if(newn > exciton->excitonbasisdim || newn < 0){ throw std::invalid_argument("Optional argument n must be a positive integer equal or below basisdim"); } @@ -411,7 +417,7 @@ void Result::writeEigenvalues(FILE* textfile, int n){ fprintf(textfile, "%d\n", exciton->excitonbasisdim); // second line: how many eigenvalues will follow - int maxEigval = (n == 0) ? exciton->excitonbasisdim : n; + int maxEigval = (newn == 0) ? exciton->excitonbasisdim : newn; fprintf(textfile, "%d\n", maxEigval); // then one eigenvalue per line @@ -427,8 +433,14 @@ void Result::writeEigenvalues(FILE* textfile, int n){ * @return void */ template -void Result::writeStates(FILE* textfile, int n){ - if(n > exciton->excitonbasisdim || n < 0){ +void Result::writeStates(FILE* textfile, int n, double encut){ + int newn = n; + + if (encut != 0.0){ + newn = (int) arma::abs(eigval - encut).index_min() + 1; + } + + if(newn > exciton->excitonbasisdim || newn < 0){ throw std::invalid_argument("Optional argument n must be a positive integer equal or below basisdim"); } // First write basis @@ -442,7 +454,7 @@ void Result::writeStates(FILE* textfile, int n){ kpoint(0), kpoint(1), kpoint(2), v, c); } - int nstates = (n == 0) ? exciton->excitonbasisdim : n; + int nstates = (newn == 0) ? exciton->excitonbasisdim : newn; for(unsigned int i = 0; i < nstates; i++){ for(unsigned int j = 0; j < exciton->excitonbasisdim; j++){ fprintf(textfile, "%11.7lf\t%11.7lf\t", @@ -458,13 +470,19 @@ void Result::writeStates(FILE* textfile, int n){ * @param textfile Textfile where the spins are written. */ template -void Result::writeSpin(int n, FILE* textfile){ - - if(n > exciton->excitonbasisdim || n < 0){ +void Result::writeSpin(int n, double encut, FILE* textfile){ + + int newn = n; + + if (encut != 0.0){ + newn = (int) arma::abs(eigval - encut).index_min() + 1; + } + + if(newn > exciton->excitonbasisdim || newn < 0){ throw std::invalid_argument("Optional argument n must be a positive integer equal or below basisdim"); } - int maxState = (n == 0) ? exciton->excitonbasisdim : n; + int maxState = (newn == 0) ? exciton->excitonbasisdim : newn; fprintf(textfile, "n\tSt\tSe\tSh\n"); for(unsigned int i = 0; i < maxState; i++){ auto spin = spinX(i); diff --git a/include/xatu/utils.hpp b/include/xatu/utils.hpp index 6f69b33..e767146 100755 --- a/include/xatu/utils.hpp +++ b/include/xatu/utils.hpp @@ -36,14 +36,20 @@ bool checkIfTriangular(const arma::cx_mat&); * @param precision Number of decimals (sets degeneracy threshold) **/ template -void printEnergies(const std::unique_ptr& results, int n = 8, int precision = 6){ +void printEnergies(const std::unique_ptr& results, int n = 8, double encut = 0.0, int precision = 6){ // Print header printf("+---------------+-----------------------------+-----------------------------+\n"); printf("| N | Eigval (eV) | Degeneracy |\n"); printf("+---------------+-----------------------------+-----------------------------+\n"); - - std::vector> pairs = detectDegeneracies(results->eigval, n, precision); + + int newn = n; + + if (encut != 0.0){ + newn = (int) arma::abs(results->eigval - encut).index_min() + 1; + } + + std::vector> pairs = detectDegeneracies(results->eigval, newn, precision); int it = 1; for(auto pair : pairs){ diff --git a/main/xatu.cpp b/main/xatu.cpp index f4e43c2..1059ddf 100755 --- a/main/xatu.cpp +++ b/main/xatu.cpp @@ -22,6 +22,7 @@ int main(int argc, char* argv[]){ TCLAP::CmdLine cmd("Command line interface options of the Xatu binary. For a more detailed description, refer to the user guide or the API documentation.", ' ', "1.0"); TCLAP::ValueArg statesArg("n", "states", "Specify number of exciton states to show.", false, 8, "No. states", cmd); + TCLAP::ValueArg energycutoff("t", "encut", "Specify up to high energy to print excitons (eV).", false, 0.0, "En cutoff", cmd); TCLAP::ValueArg precisionArg("p", "precision", "Desired energy precision. Used to compute degeneracies.", false, 6, "No. decimals", cmd); TCLAP::SwitchArg spinArg("s", "spin", "Compute exciton spin and write it to file.", cmd, false); TCLAP::ValueArg dftArg("d", "dft", "Indicates that the system file is a .outp CRYSTAL file.", false, -1, "No. Fock matrices", cmd); @@ -50,6 +51,7 @@ int main(int argc, char* argv[]){ // Extract information from parsed CLI options int nstates = statesArg.getValue(); + double encut = energycutoff.getValue(); int ncells = dftArg.getValue(); int electronNum = w90Arg.getValue(); int decimals = precisionArg.getValue(); @@ -152,7 +154,7 @@ int main(int argc, char* argv[]){ cout << "| Results |" << endl; cout << "+---------------------------------------------------------------------------+" << endl; - xatu::printEnergies(results, nstates, decimals); + xatu::printEnergies(results, nstates, encut, decimals); cout << "+---------------------------------------------------------------------------+" << endl; cout << "| Output |" << endl; @@ -168,7 +170,7 @@ int main(int argc, char* argv[]){ std::cout << "Writing eigvals to file: " << filename_en << std::endl; fprintf(textfile_en, "%d\n", excitonConfig->excitonInfo.ncell); - results->writeEigenvalues(textfile_en, nstates); + results->writeEigenvalues(textfile_en, nstates, encut); fclose(textfile_en); } @@ -179,7 +181,7 @@ int main(int argc, char* argv[]){ FILE* textfile_st = fopen(filename_st.c_str(), "w"); std::cout << "Writing states to file: " << filename_st << std::endl; - results->writeStates(textfile_st, nstates); + results->writeStates(textfile_st, nstates, encut); fclose(textfile_st); } @@ -231,7 +233,7 @@ int main(int argc, char* argv[]){ FILE* textfile_spin = fopen(filename_spin.c_str(), "w"); std::cout << "Writing excitons spin fo file: " << filename_spin << std::endl; - results->writeSpin(nstates, textfile_spin); + results->writeSpin(nstates, encut, textfile_spin); } auto stop = high_resolution_clock::now();