1 /* Ergo, version 3.8, a program for linear scaling electronic structure 2 * calculations. 3 * Copyright (C) 2019 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, 4 * and Anastasia Kruchinina. 5 * 6 * This program is free software: you can redistribute it and/or modify 7 * it under the terms of the GNU General Public License as published by 8 * the Free Software Foundation, either version 3 of the License, or 9 * (at your option) any later version. 10 * 11 * This program is distributed in the hope that it will be useful, 12 * but WITHOUT ANY WARRANTY; without even the implied warranty of 13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 14 * GNU General Public License for more details. 15 * 16 * You should have received a copy of the GNU General Public License 17 * along with this program. If not, see <http://www.gnu.org/licenses/>. 18 * 19 * Primary academic reference: 20 * Ergo: An open-source program for linear-scaling electronic structure 21 * calculations, 22 * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia 23 * Kruchinina, 24 * SoftwareX 7, 107 (2018), 25 * <http://dx.doi.org/10.1016/j.softx.2018.03.005> 26 * 27 * For further information about Ergo, see <http://www.ergoscf.org>. 28 */ 29 30 /** @file GetDensFromFock.h 31 * 32 * @brief Routines for getting density matrix from a given Fock 33 * matrix. 34 * 35 * @author Anastasia Kruchinina <em>responsible</em> 36 * @author Elias Rudberg 37 */ 38 39 #ifndef GETDENSFROMFOCKHEADER 40 #define GETDENSFROMFOCKHEADER 41 42 #include "realtype.h" 43 #include "matrix_typedefs.h" 44 #include "matrix_typedefs_chtml.h" 45 #include "transform.h" 46 #include "output.h" 47 48 49 50 /** GetDensFromFock class containing parameters and functions for computing density matrix. 51 * 52 * Flags are set to undefined value by default. User should define 53 * them explicitly, otherwise exception is thrown if undefined flag is 54 * used. 55 */ 56 class GetDensFromFock 57 { 58 public: 59 60 static const int UNDEF_VALUE_UINT; 61 static const std::string NA_STRING; 62 63 static const bool BOOL_TRUE; 64 static const bool BOOL_FALSE; 65 66 void create_checkpoint(symmMatrix& Finput, /**< [in] Effective Hamiltonian matrix (written to file) */ 67 symmMatrix& F_ort_prev, /**< [in/out] 68 * Input: Previous F matrix in orthogonal basis. (written to file) 69 * Output: New F matrix in orthogonal basis ( ZT*Finput*Z ). (written to file) */ 70 generalVector *eigVecLUMO, /**< [out] LUMO eigenvector */ 71 generalVector *eigVecHOMO, /**< [out] HOMO eigenvector */ 72 std::string IDstr /**< [in] File identificator; added to the name of each file */ 73 ); 74 75 76 static void restore_from_checkpoint(GetDensFromFock& DensFromFock, /**< [out] Instance of GetDensFromFock class contatining all data for computing the density matrix */ 77 symmMatrix& Finput, /**< [out] Effective Hamiltonian matrix (written to file) */ 78 symmMatrix& F_ort_prev, /**< [out] F matrix in orthogonal basis ( ZT*Finput*Z ). (written to file) */ 79 generalVector *eigVecLUMO, /**< [out] LUMO eigenvector */ 80 generalVector *eigVecHOMO, /**< [out] HOMO eigenvector */ 81 std::string checkpoint_path, /**< [out] HOMO eigenvector */ 82 std::string IDstr, /**< [in] File identificator; added to the name of each file. */ 83 int SCF_step /**< [in] SCF step which should be restored; added to the name of each file in given SCF cycle. */ 84 ); 85 86 87 //constructor GetDensFromFock()88 GetDensFromFock() 89 { 90 do_output(LOG_CAT_INFO, LOG_AREA_DENSFROMF, "Create object from GetDensFromFock."); 91 92 // set all variables and flags to default values 93 n = UNDEF_VALUE_UINT; 94 noOfOccupiedOrbs = UNDEF_VALUE_UINT; 95 factor = UNDEF_VALUE_UINT; 96 invCholFactor_euclnorm = 0; 97 maxMul = 0; 98 plot_puri_results = BOOL_FALSE; 99 SCF_step = UNDEF_VALUE_UINT; 100 101 use_diagonalization = BOOL_FALSE; 102 use_purification = BOOL_FALSE; 103 104 electronicTemperature = 0; 105 gap_expected_lower_bound = 0; 106 eigvalueErrorLimit = 0; 107 subspaceErrorLimit = 0; 108 puri_eig_acc_factor_for_guess = 0; 109 use_diag_on_error = BOOL_FALSE; 110 use_diag_on_error_guess = BOOL_FALSE; 111 112 create_m_files = BOOL_FALSE; 113 output_homo_and_lumo_eigenvectors = BOOL_FALSE; 114 number_of_occupied_eigenvectors = 0; 115 number_of_unoccupied_eigenvectors = 0; 116 go_back_X_iter_proj_method = 0; 117 jump_over_X_iter_proj_method = 0; 118 ignore_purification_failure = BOOL_FALSE; 119 use_rand_perturbation_for_alleigsint = BOOL_FALSE; 120 use_acceleration = BOOL_FALSE; 121 use_new_stopping_criterion = BOOL_FALSE; 122 store_all_eigenvalues_to_file = BOOL_FALSE; 123 try_eigv_on_next_iteration_if_fail = BOOL_FALSE; 124 125 leavesSizeMax = UNDEF_VALUE_UINT; 126 blocksize = UNDEF_VALUE_UINT; 127 128 eigenvectors_method = NA_STRING; 129 eigenvectors_iterative_method = NA_STRING; 130 use_prev_vector_as_initial_guess = BOOL_FALSE; 131 puri_compute_eigv_in_each_iteration = BOOL_FALSE; 132 run_shift_and_square_method_on_F = BOOL_FALSE; 133 save_permuted_F_matrix_in_bin = BOOL_FALSE; 134 135 eigensolver_accuracy = 0; 136 eigensolver_maxiter = 0; 137 138 std::string stats_prefix = ""; // default value 139 140 truncationNormPurification = mat::euclNorm; // default value 141 stopCriterionNormPurification = mat::euclNorm; // default value 142 143 clean_eigs_intervals(); 144 clean_puri_stats(); 145 146 147 filenameFinput = "matrix_Finput"; 148 filenameF_ort_prev = "matrix_F_ort_prev"; 149 filenameeigVecLUMO = "vector_eigVecLUMO"; 150 filenameeigVecHOMO = "vector_eigVecHOMO"; 151 filenameOverlap = "matrix_Overlap"; 152 filenameD_ort_prev = "matrix_D_ort_prev"; 153 filenameinvCholFactor = "matrix_invCholFactor"; 154 file_for_basic_types = "basic_types"; 155 } 156 157 /** Choose which method to use for computing the density matrix from Fock matrix. 158 */ 159 int get_dens_from_fock(symmMatrix& Finput, /**< [in] Effective Hamiltonian matrix. (written to file) */ 160 symmMatrix& resultDens, /**< [out] Density matrix. (written to file) */ 161 symmMatrix& F_ort_prev /**< [in/out] 162 * Input: Previous F matrix in orthogonal basis. (written to file) 163 * Output: New F matrix in orthogonal basis ( ZT*Finput*Z ). (written to file) */ 164 ); 165 166 167 168 /** Use recursive expansion for computing the density matrix from Fock matrix. 169 */ 170 int get_dens_from_fock_sparse( 171 symmMatrix& F, /**< [in] Effective Hamiltonian matrix. (written to file) */ 172 symmMatrix& resultDens, /**< [out] Density matrix. (written to file) */ 173 symmMatrix& F_ort_prev /**< [in/out] 174 Input: Previous F matrix in orthogonal basis. (written to file) 175 Output: New F matrix in orthogonal basis ( ZT*Finput*Z ). (written to file) 176 */ 177 ); 178 179 get_computed_eigenpairs(std::vector<generalVector> & eigVecUNOCCref,std::vector<generalVector> & eigVecOCCref,std::vector<ergo_real> & eigValUNOCCref,std::vector<ergo_real> & eigValOCCref)180 void get_computed_eigenpairs( 181 std::vector<generalVector> &eigVecUNOCCref, /**< [out] Unoccupied eigenvectors */ 182 std::vector<generalVector> &eigVecOCCref, /**< [out] Occupied eigenvectors */ 183 std::vector<ergo_real> &eigValUNOCCref, /**< [out] Unoccupied eigenvalues */ 184 std::vector<ergo_real> &eigValOCCref /**< [out] Occupied eigenvalues */ 185 ) 186 { 187 eigVecUNOCCref = eigVecUNOCC; 188 eigVecOCCref = eigVecOCC; 189 eigValOCCref = eigValOCC; 190 eigValUNOCCref = eigValUNOCC; 191 } 192 193 /** Set bounds for HOMO and LUMO eigenvalues to -/+ inf, thus remove 194 * any known bounds. 195 */ clean_eigs_intervals()196 inline void clean_eigs_intervals() 197 { 198 homoInterval_Finput = intervalType(-1e22, 1e22); 199 lumoInterval_Finput = intervalType(-1e22, 1e22); 200 homoInterval_Finput = intervalType(-1e22, 1e22); 201 lumoInterval_Finput = intervalType(-1e22, 1e22); 202 203 homoInterval_F_ort_prev = intervalType(-1e22, 1e22); 204 lumoInterval_F_ort_prev = intervalType(-1e22, 1e22); 205 homoInterval_F_ort_prev = intervalType(-1e22, 1e22); 206 lumoInterval_F_ort_prev = intervalType(-1e22, 1e22); 207 } 208 209 210 211 set_SCF_step(int step)212 inline void set_SCF_step(int step /**< [in] Current SCF step */) 213 { SCF_step = step; } unset_SCF_step()214 inline void unset_SCF_step() 215 { SCF_step = UNDEF_VALUE_UINT; } 216 217 218 /** Plot figures from the recursive expansion. 219 */ 220 inline void set_generate_figures(std::string str = "" /**< [in] String added to each generated file. */) 221 { 222 if (create_m_files) 223 { 224 assert(SCF_step >= 0); 225 plot_puri_results = BOOL_TRUE; 226 plot_puri_results_str = str; 227 } 228 } 229 230 /** Do not plot figures from the recursive expansion. 231 */ unset_generate_figures()232 inline void unset_generate_figures() 233 { 234 plot_puri_results = BOOL_FALSE; 235 plot_puri_results_str = ""; 236 } 237 set_general_params(const int n_,mat::SizesAndBlocks const & matrixSizesAndBlocks_)238 inline void set_general_params(const int n_, /**< [in] Number of basis functions. */ 239 mat::SizesAndBlocks const& matrixSizesAndBlocks_ /**< [in] Matrix library parameters. */ 240 ) 241 { 242 assert(n_ >= 1); 243 n = n_; 244 matrixSizesAndBlocks = matrixSizesAndBlocks_; 245 } 246 set_cht_matrix_params(const int leavesSizeMax_,const int blocksize_)247 inline void set_cht_matrix_params(const int leavesSizeMax_, /**< [in] CHTMatrix library parameter leavesSizeMax. */ 248 const int blocksize_ /**< [in] CHTMatrix library parameter blocksize. */ 249 ) 250 { 251 assert(leavesSizeMax_ >= 1); 252 assert(blocksize_ >= 1); 253 leavesSizeMax = leavesSizeMax_; 254 blocksize = blocksize_; 255 } 256 get_SizesAndBlocks(mat::SizesAndBlocks & matrixSizesAndBlocks_)257 inline void get_SizesAndBlocks(mat::SizesAndBlocks& matrixSizesAndBlocks_ /**< [out] Matrix library parameters. */ 258 ) const 259 { 260 matrixSizesAndBlocks_ = matrixSizesAndBlocks; 261 } 262 263 /** Set truncation norm used in the recursive expansion. 264 * Possible norms: spectral, Frobenius or mixed. 265 */ set_truncationNormPurification(mat::normType const truncationNormPurification_)266 inline void set_truncationNormPurification(mat::normType const truncationNormPurification_ /**< [in] Norm used in truncation. */) 267 { truncationNormPurification = truncationNormPurification_; } 268 269 /** Set stopping criterion norm used in the recursive expansion. 270 * Possible norms: spectral, Frobenius or mixed. 271 */ set_stopCriterionNormPurification(mat::normType const stopCriterionNormPurification_)272 inline void set_stopCriterionNormPurification(mat::normType const stopCriterionNormPurification_ /**< [in] Norm used in the stopping criterion. */) 273 { stopCriterionNormPurification = stopCriterionNormPurification_; } 274 275 do_restricted_calculations()276 inline void do_restricted_calculations() 277 { factor = 2; } 278 do_unrestricted_calculations()279 inline void do_unrestricted_calculations() 280 { factor = 1; } 281 set_no_occupied_orbs(int noOfOccupiedOrbs_)282 inline void set_no_occupied_orbs(int noOfOccupiedOrbs_) 283 { 284 assert(noOfOccupiedOrbs_ >= 0); 285 noOfOccupiedOrbs = noOfOccupiedOrbs_; 286 } 287 clean_puri_stats()288 inline void clean_puri_stats() 289 { puri_stats.clear(); } 290 291 set_invCholFactor(triangMatrix const & invCholFactor_,ergo_real invCholFactor_euclnorm_)292 inline void set_invCholFactor(triangMatrix const& invCholFactor_, 293 ergo_real invCholFactor_euclnorm_) 294 { 295 invCholFactor = invCholFactor_; 296 assert(invCholFactor_euclnorm_ >= 0); 297 invCholFactor_euclnorm = invCholFactor_euclnorm_; 298 } 299 set_gap_expected_lower_bound(ergo_real gap_expected_lower_bound_)300 inline void set_gap_expected_lower_bound(ergo_real gap_expected_lower_bound_) 301 { 302 assert(gap_expected_lower_bound_ >= 0); 303 gap_expected_lower_bound = gap_expected_lower_bound_; 304 } 305 306 /** Set maximum allowed number of iterations in recursive expansion. 307 */ set_purification_maxmul(ergo_real purification_maxmul_)308 inline void set_purification_maxmul(ergo_real purification_maxmul_) 309 { 310 assert(purification_maxmul_ > 0); 311 maxMul = purification_maxmul_; 312 } 313 set_number_of_eigenvectors_to_compute(int occ,int unocc)314 inline void set_number_of_eigenvectors_to_compute(int occ, int unocc) 315 { 316 number_of_occupied_eigenvectors = occ; 317 number_of_unoccupied_eigenvectors = unocc; 318 } 319 set_projection_method_params(int go_back,int step)320 inline void set_projection_method_params(int go_back, int step) 321 { 322 go_back_X_iter_proj_method = go_back; 323 jump_over_X_iter_proj_method = step; 324 } 325 326 /**** SET/UNSET SECTION *****/ 327 get_purification_create_m_files()328 inline int get_purification_create_m_files() const 329 { return create_m_files == BOOL_TRUE; } set_purification_create_m_files()330 inline void set_purification_create_m_files() 331 { create_m_files = BOOL_TRUE; } unset_purification_create_m_files()332 inline void unset_purification_create_m_files() 333 { create_m_files = BOOL_FALSE; } 334 335 336 get_output_homo_and_lumo_eigenvectors()337 inline int get_output_homo_and_lumo_eigenvectors() const 338 { return output_homo_and_lumo_eigenvectors == BOOL_TRUE; } set_output_homo_and_lumo_eigenvectors()339 inline void set_output_homo_and_lumo_eigenvectors() 340 { output_homo_and_lumo_eigenvectors = BOOL_TRUE; } unset_output_homo_and_lumo_eigenvectors()341 inline void unset_output_homo_and_lumo_eigenvectors() 342 { output_homo_and_lumo_eigenvectors = BOOL_FALSE; } 343 344 get_purification_ignore_failure()345 inline int get_purification_ignore_failure() const 346 { return ignore_purification_failure == BOOL_TRUE; } set_purification_ignore_failure()347 inline void set_purification_ignore_failure() 348 { ignore_purification_failure = BOOL_TRUE; } unset_purification_ignore_failure()349 inline void unset_purification_ignore_failure() 350 { ignore_purification_failure = BOOL_FALSE; } 351 352 get_use_rand_perturbation_for_alleigsint()353 inline int get_use_rand_perturbation_for_alleigsint() const 354 { return use_rand_perturbation_for_alleigsint == BOOL_TRUE; } set_purification_use_rand_perturbation_for_alleigsint()355 inline void set_purification_use_rand_perturbation_for_alleigsint() 356 { use_rand_perturbation_for_alleigsint = BOOL_TRUE; } unset_purification_use_rand_perturbation_for_alleigsint()357 inline void unset_purification_use_rand_perturbation_for_alleigsint() 358 { use_rand_perturbation_for_alleigsint = BOOL_FALSE; } 359 get_use_diagonalization()360 inline int get_use_diagonalization() const 361 { return use_diagonalization == BOOL_TRUE; } set_use_diagonalization()362 inline void set_use_diagonalization() 363 { use_diagonalization = BOOL_TRUE; } unset_use_diagonalization()364 inline void unset_use_diagonalization() 365 { use_diagonalization = BOOL_FALSE; } 366 367 get_use_purification()368 inline int get_use_purification() const 369 { return use_purification == BOOL_TRUE; } set_use_purification()370 inline void set_use_purification() 371 { use_purification = BOOL_TRUE; } unset_use_purification()372 inline void unset_use_purification() 373 { use_purification = BOOL_FALSE; } 374 375 get_use_diag_on_error_guess()376 inline int get_use_diag_on_error_guess() const 377 { return use_diag_on_error_guess == BOOL_TRUE; } set_use_diag_on_error_guess()378 inline void set_use_diag_on_error_guess() 379 { use_diag_on_error_guess = BOOL_TRUE; } unset_use_diag_on_error_guess()380 inline void unset_use_diag_on_error_guess() 381 { use_diag_on_error_guess = BOOL_FALSE; } 382 383 get_use_diag_on_error()384 inline int get_use_diag_on_error() const 385 { return use_diag_on_error == BOOL_TRUE; } set_use_diag_on_error()386 inline void set_use_diag_on_error() 387 { use_diag_on_error = BOOL_TRUE; } unset_use_diag_on_error()388 inline void unset_use_diag_on_error() 389 { use_diag_on_error = BOOL_FALSE; } 390 391 get_stats_prefix()392 inline std::string get_stats_prefix() const 393 { return stats_prefix; } set_stats_prefix(std::string stats_prefix_)394 inline void set_stats_prefix(std::string stats_prefix_) 395 { stats_prefix = stats_prefix_; } unset_stats_prefix()396 inline void unset_stats_prefix() 397 { stats_prefix = ""; } 398 399 get_use_acceleration()400 inline int get_use_acceleration() const 401 { return use_acceleration == BOOL_TRUE; } set_use_acceleration()402 inline void set_use_acceleration() 403 { use_acceleration = BOOL_TRUE; } unset_use_acceleration()404 inline void unset_use_acceleration() 405 { use_acceleration = BOOL_FALSE; } 406 get_use_new_stopping_criterion()407 inline int get_use_new_stopping_criterion() const 408 { return use_new_stopping_criterion == BOOL_TRUE; } set_use_new_stopping_criterion()409 inline void set_use_new_stopping_criterion() 410 { use_new_stopping_criterion = BOOL_TRUE; } unset_use_new_stopping_criterion()411 inline void unset_use_new_stopping_criterion() 412 { use_new_stopping_criterion = BOOL_FALSE; } 413 get_store_all_eigenvalues_to_file()414 inline int get_store_all_eigenvalues_to_file() const 415 { return store_all_eigenvalues_to_file == BOOL_TRUE; } set_store_all_eigenvalues_to_file()416 inline void set_store_all_eigenvalues_to_file() 417 { store_all_eigenvalues_to_file = BOOL_TRUE; } unset_store_all_eigenvalues_to_file()418 inline void unset_store_all_eigenvalues_to_file() 419 { store_all_eigenvalues_to_file = BOOL_FALSE; } 420 get_save_permuted_F_matrix_in_bin()421 inline int get_save_permuted_F_matrix_in_bin() 422 { return save_permuted_F_matrix_in_bin == BOOL_TRUE; } set_save_permuted_F_matrix_in_bin()423 inline void set_save_permuted_F_matrix_in_bin() 424 { save_permuted_F_matrix_in_bin = BOOL_TRUE; } unset_save_permuted_F_matrix_in_bin()425 inline void unset_save_permuted_F_matrix_in_bin() 426 { save_permuted_F_matrix_in_bin = BOOL_FALSE; } 427 428 429 get_puri_compute_eigv_in_each_iteration()430 inline int get_puri_compute_eigv_in_each_iteration() 431 { return puri_compute_eigv_in_each_iteration == BOOL_TRUE; } set_puri_compute_eigv_in_each_iteration()432 inline void set_puri_compute_eigv_in_each_iteration() 433 { puri_compute_eigv_in_each_iteration = BOOL_TRUE; } unset_puri_compute_eigv_in_each_iteration()434 inline void unset_puri_compute_eigv_in_each_iteration() 435 { puri_compute_eigv_in_each_iteration = BOOL_FALSE; } 436 437 get_run_shift_and_square_method_on_F()438 inline int get_run_shift_and_square_method_on_F() 439 { return run_shift_and_square_method_on_F == BOOL_TRUE; } set_run_shift_and_square_method_on_F()440 inline void set_run_shift_and_square_method_on_F() 441 { run_shift_and_square_method_on_F = BOOL_TRUE; } unset_run_shift_and_square_method_on_F()442 inline void unset_run_shift_and_square_method_on_F() 443 { run_shift_and_square_method_on_F = BOOL_FALSE; } 444 445 get_try_eigv_on_next_iteration_if_fail()446 inline int get_try_eigv_on_next_iteration_if_fail() 447 { return try_eigv_on_next_iteration_if_fail == BOOL_TRUE; } set_try_eigv_on_next_iteration_if_fail()448 inline void set_try_eigv_on_next_iteration_if_fail() 449 { try_eigv_on_next_iteration_if_fail = BOOL_TRUE; } unset_try_eigv_on_next_iteration_if_fail()450 inline void unset_try_eigv_on_next_iteration_if_fail() 451 { try_eigv_on_next_iteration_if_fail = BOOL_FALSE; } 452 453 get_use_prev_vector_as_initial_guess()454 inline int get_use_prev_vector_as_initial_guess() 455 { return use_prev_vector_as_initial_guess == BOOL_TRUE; } set_use_prev_vector_as_initial_guess()456 inline void set_use_prev_vector_as_initial_guess() 457 { use_prev_vector_as_initial_guess = BOOL_TRUE; } unset_use_prev_vector_as_initial_guess()458 inline void unset_use_prev_vector_as_initial_guess() 459 { use_prev_vector_as_initial_guess = BOOL_FALSE; } 460 461 set_diagonalization_params(ergo_real electronicTemperature_,symmMatrix & overlapMatrix_)462 inline void set_diagonalization_params(ergo_real electronicTemperature_, 463 symmMatrix& overlapMatrix_) 464 { 465 set_overlapMatrix(overlapMatrix_); 466 assert(electronicTemperature_ >= 0); 467 electronicTemperature = electronicTemperature_; 468 } 469 set_overlapMatrix(symmMatrix & overlapMatrix_)470 inline void set_overlapMatrix(symmMatrix& overlapMatrix_) 471 { overlapMatrix = overlapMatrix_; } 472 473 474 inline void set_purification_limits(ergo_real subspaceErrorLimit_, 475 ergo_real eigvalueErrorLimit_ = 0, 476 ergo_real puri_eig_acc_factor_for_guess = 0) 477 { 478 set_eigvalueErrorLimit(eigvalueErrorLimit_); 479 set_subspaceErrorLimit(subspaceErrorLimit_); 480 set_puri_eig_acc_factor_for_guess(puri_eig_acc_factor_for_guess); 481 } 482 483 /** Set maximum allowed error in eigenvalues of the density matrix. 484 */ set_eigvalueErrorLimit(ergo_real eigvalueErrorLimit_)485 inline void set_eigvalueErrorLimit(ergo_real eigvalueErrorLimit_) 486 { eigvalueErrorLimit = eigvalueErrorLimit_; } 487 488 /** Set maximum allowed error in invariant subspaces of the density 489 * matrix. 490 */ set_subspaceErrorLimit(ergo_real subspaceErrorLimit_)491 inline void set_subspaceErrorLimit(ergo_real subspaceErrorLimit_) 492 { subspaceErrorLimit = subspaceErrorLimit_; } 493 494 /** Set puri_eig_acc_factor_for_guess parameter. 495 * 496 * Obsolete parameter needed for the old stopping criterion for 497 * creating the initial guess. 498 */ set_puri_eig_acc_factor_for_guess(ergo_real puri_eig_acc_factor_for_guess_)499 inline void set_puri_eig_acc_factor_for_guess(ergo_real puri_eig_acc_factor_for_guess_) 500 { puri_eig_acc_factor_for_guess = puri_eig_acc_factor_for_guess_; } 501 502 503 504 // get some results from the purification 505 get_result_entropy_term()506 ergo_real get_result_entropy_term() const 507 { return resultEntropyTerm; } 508 get_puri_stats(std::map<std::string,double> & puri_stats_)509 inline void get_puri_stats(std::map<std::string, double>& puri_stats_) const 510 { puri_stats_ = puri_stats; } 511 512 513 514 // Fprev is effective Hamiltonian matrix (=Finput) 515 // Intervals contain the homo and lumo eigenvalues of Fprev set_eigs_Fprev(intervalType & homoInterval_Finput_,intervalType & lumoInterval_Finput_)516 inline void set_eigs_Fprev(intervalType& homoInterval_Finput_, 517 intervalType& lumoInterval_Finput_) 518 { 519 homoInterval_Finput = intervalType(homoInterval_Finput_); 520 lumoInterval_Finput = intervalType(lumoInterval_Finput_); 521 } 522 get_eigs_Fprev(intervalType & homoInterval_Finput_,intervalType & lumoInterval_Finput_)523 inline void get_eigs_Fprev(intervalType& homoInterval_Finput_, 524 intervalType& lumoInterval_Finput_) const 525 { 526 homoInterval_Finput_ = intervalType(homoInterval_Finput); 527 lumoInterval_Finput_ = intervalType(lumoInterval_Finput); 528 } 529 530 // F_ort_prev is matrix in orthogonal basis 531 // Intervals contain the homo and lumo eigenvalues of F_ort_prev set_eigs_F_ort_prev(intervalType & homoInterval_F_ort_prev_,intervalType & lumoInterval_F_ort_prev_)532 inline void set_eigs_F_ort_prev(intervalType& homoInterval_F_ort_prev_, 533 intervalType& lumoInterval_F_ort_prev_) 534 { 535 homoInterval_F_ort_prev = intervalType(homoInterval_F_ort_prev_); 536 lumoInterval_F_ort_prev = intervalType(lumoInterval_F_ort_prev_); 537 } 538 get_eigs_F_ort_prev(intervalType & homoInterval_F_ort_prev_,intervalType & lumoInterval_F_ort_prev_)539 inline void get_eigs_F_ort_prev(intervalType& homoInterval_F_ort_prev_, 540 intervalType& lumoInterval_F_ort_prev_) const 541 { 542 homoInterval_F_ort_prev_ = intervalType(homoInterval_F_ort_prev); 543 lumoInterval_F_ort_prev_ = intervalType(lumoInterval_F_ort_prev); 544 } 545 get_eigvalueErrorLimit()546 inline ergo_real get_eigvalueErrorLimit() const 547 { return eigvalueErrorLimit; } get_subspaceErrorLimit()548 inline ergo_real get_subspaceErrorLimit() const 549 { return subspaceErrorLimit; } get_puri_eig_acc_factor_for_guess()550 inline ergo_real get_puri_eig_acc_factor_for_guess() const 551 { return puri_eig_acc_factor_for_guess; } 552 compute_eigenvectors(std::string eigenvectors_method_,std::string eigenvectors_iterative_method_,ergo_real eigensolver_accuracy_,int eigensolver_maxiter_,int use_prev_vector_as_initial_guess_,int try_eigv_on_next_iteration_if_fail_)553 inline void compute_eigenvectors(std::string eigenvectors_method_, 554 std::string eigenvectors_iterative_method_, 555 ergo_real eigensolver_accuracy_, 556 int eigensolver_maxiter_, 557 int use_prev_vector_as_initial_guess_, 558 int try_eigv_on_next_iteration_if_fail_) 559 { 560 assert(eigenvectors_method_ == "square" || eigenvectors_method_ == "projection"); 561 eigenvectors_method = eigenvectors_method_; 562 563 assert(eigenvectors_iterative_method_ == "power" || eigenvectors_iterative_method_ == "lanczos"); 564 eigenvectors_iterative_method = eigenvectors_iterative_method_; 565 566 eigensolver_accuracy = eigensolver_accuracy_; 567 eigensolver_maxiter = eigensolver_maxiter_; 568 569 if (use_prev_vector_as_initial_guess_ > 0) 570 { 571 set_use_prev_vector_as_initial_guess(); 572 } 573 else 574 { 575 unset_use_prev_vector_as_initial_guess(); 576 } 577 if (try_eigv_on_next_iteration_if_fail_ > 0) 578 { 579 set_try_eigv_on_next_iteration_if_fail(); 580 } 581 else 582 { 583 unset_try_eigv_on_next_iteration_if_fail(); 584 } 585 586 } 587 compute_eigenvectors_extra(int puri_compute_eigv_in_each_iteration_,int run_shift_and_square_method_on_F_)588 inline void compute_eigenvectors_extra(int puri_compute_eigv_in_each_iteration_, int run_shift_and_square_method_on_F_) 589 { 590 if (puri_compute_eigv_in_each_iteration_ > 0) 591 { 592 set_puri_compute_eigv_in_each_iteration(); 593 } 594 else 595 { 596 unset_puri_compute_eigv_in_each_iteration(); 597 } 598 599 if (run_shift_and_square_method_on_F_ > 0) 600 { 601 set_run_shift_and_square_method_on_F(); 602 } 603 else 604 { 605 unset_run_shift_and_square_method_on_F(); 606 } 607 } 608 609 private: 610 611 int SCF_step; 612 613 614 615 bool use_diagonalization; /**< Flag to turn on diagonalization. */ 616 617 bool use_purification; /**< Flag to turn on purification. */ 618 619 bool store_all_eigenvalues_to_file; /**< Store eigenvalues to the file when doing diagonalization. 620 * NOTE: works just with diagonalization */ 621 622 bool try_eigv_on_next_iteration_if_fail; /**< For square method: if eigenvector is not computed in iteration 623 * i, try to compute it in iteration i+1 */ 624 625 ergo_real electronicTemperature; /**< Electronic temperature */ 626 627 ergo_real gap_expected_lower_bound; /**< Expected lower bound for the gap to be used in early iterations. */ 628 629 ergo_real eigvalueErrorLimit; /**< Tolerated deviation of eigenvalues from 0 and 1 in the computed density matrix. */ 630 631 ergo_real subspaceErrorLimit; /**< Tolerated error in the occupied subspace as measured by the sinus of the largest canonical angle. */ 632 633 ergo_real puri_eig_acc_factor_for_guess; /**< With this number will be multiplied the tolerated deviation of 634 * eigenvalues from 0 and 1 in the computed density matrix for the 635 * initial guess density matrix */ 636 637 bool use_diag_on_error; /**< Flag to fall back on diagonalization if purification fails. */ 638 bool use_diag_on_error_guess; 639 640 bool create_m_files; /**< Flag to create m-files with information about the purification process. */ 641 642 bool output_homo_and_lumo_eigenvectors; /**< Compute homo and lumo eigenvectors and write them to the file */ 643 644 int number_of_occupied_eigenvectors; /**< Number of occupied eigenvectors to compute. */ 645 646 int number_of_unoccupied_eigenvectors; /**< Number of unoccupied eigenvectors to compute. */ 647 648 int go_back_X_iter_proj_method; /**< Parameter used in the projection method for computing eigenvectors. Defines the iteration to start computation of the eigenvectors. */ 649 650 int jump_over_X_iter_proj_method; /**< Parameter used in the projection method for computing eigenvectors. Defines how many iterations will be skipped before the next attempt if some eigenvectors are not computed. */ 651 652 bool use_prev_vector_as_initial_guess; /**< Use eigenvector from the previous SCF cycle as an initial guess in this cycle */ 653 654 bool puri_compute_eigv_in_each_iteration; /**< Compute eigenvectors in each iteration of the recursive expansion. */ 655 656 bool run_shift_and_square_method_on_F; /**< (for comparison) Run shift_and_square method to 657 * get eigenvectors of the matrix F for various shifts. */ 658 659 bool save_permuted_F_matrix_in_bin; /**< Save sparse matrix F into bin file in the current permutation of rows and columns. */ 660 661 bool ignore_purification_failure; /**< Continue even if purification fails to converge. */ 662 663 bool use_rand_perturbation_for_alleigsint; /**< Apply a random 664 * perturbation to (try 665 * to) improve the 666 * convergence speed of 667 * Lanczos calculation of 668 * extremal 669 * eigenvalues. */ 670 671 std::string stats_prefix; /**< Prefix to be added to statistics files. */ 672 673 bool plot_puri_results; /**< Plot results of the purification from this function call */ 674 675 std::string plot_puri_results_str; 676 677 bool use_acceleration; /**< Use acceleration in the purification */ 678 679 bool use_new_stopping_criterion; /**< Use new parameterless stopping criterion */ 680 681 std::string eigenvectors_method; /**< Method for computing eigenvectors: square or projection */ 682 683 std::string eigenvectors_iterative_method; /**< Iterative method for computing eigenvectors: power or lanczos */ 684 685 ergo_real eigensolver_accuracy; /**< The accuracy for the eigenvalue problem solver */ 686 687 int eigensolver_maxiter; /**< Maximum number of iterations for the eigenvalue problem solver */ 688 689 int n; /**< System size. */ 690 691 int noOfOccupiedOrbs; /**< Number of occupied orbitals. */ 692 693 ergo_real factor; /**< Factor to scale the resulting density matrix. 694 * (for restricted vs unrestricted calc) */ 695 696 symmMatrix overlapMatrix; /**< Overlap matrix (written to file) */ 697 698 symmMatrix D_ort_prev; /**< Density matrix from previous SCF cycle (written to file) */ 699 700 triangMatrix invCholFactor; /**< Inverse Cholesky factor (written to file) */ 701 702 ergo_real invCholFactor_euclnorm; /**< Euclidean norm of inverse Cholesky factor. */ 703 704 mat::normType truncationNormPurification; /**< Norm to be used for truncation. */ 705 706 mat::normType stopCriterionNormPurification; /**< Norm to be used for stopping criterion. */ 707 708 int maxMul; /**< Maximum allowed number of matrix multiplications in the purification */ 709 710 mat::SizesAndBlocks matrixSizesAndBlocks; /**< Information about HML matrix block sizes etc. */ 711 712 int leavesSizeMax; /**< Information about leavesSizeMax and blocksize for CHTMatrix */ 713 int blocksize; /**< Information about leavesSizeMax and blocksize for CHTMatrix */ 714 715 intervalType homoInterval_Finput; 716 intervalType lumoInterval_Finput; 717 718 intervalType homoInterval_F_ort_prev; 719 intervalType lumoInterval_F_ort_prev; 720 721 ergo_real resultEntropyTerm; 722 std::map<std::string, double> puri_stats; 723 724 std::vector<generalVector> eigVecOCC; 725 std::vector<generalVector> eigVecUNOCC; 726 std::vector<ergo_real> eigValOCC; 727 std::vector<ergo_real> eigValUNOCC; 728 729 // Names of files needed for checkpoints 730 const char *filenameFinput; 731 const char *filenameF_ort_prev; 732 const char *filenameeigVecLUMO; 733 const char *filenameeigVecHOMO; 734 const char *filenameOverlap; 735 const char *filenameD_ort_prev; 736 const char *filenameinvCholFactor; 737 const char *file_for_basic_types; 738 }; 739 740 741 742 #endif // GETDENSFROMFOCKHEADER 743