ergo
scf_utils.h
Go to the documentation of this file.
1 /* Ergo, version 3.4, a program for linear scaling electronic structure
2  * calculations.
3  * Copyright (C) 2014 Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek.
4  *
5  * This program is free software: you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation, either version 3 of the License, or
8  * (at your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program. If not, see <http://www.gnu.org/licenses/>.
17  *
18  * Primary academic reference:
19  * Kohn−Sham Density Functional Theory Electronic Structure Calculations
20  * with Linearly Scaling Computational Time and Memory Usage,
21  * Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek,
22  * J. Chem. Theory Comput. 7, 340 (2011),
23  * <http://dx.doi.org/10.1021/ct100611z>
24  *
25  * For further information about Ergo, see <http://www.ergoscf.org>.
26  */
27 
28 #ifndef SCF_UTILS_HEADER
29 #define SCF_UTILS_HEADER
30 
31 #include "molecule.h"
32 #include "basisinfo.h"
33 #include "integrals_2el.h"
34 #include "matrix_typedefs.h"
35 #include "densityfitting.h"
36 #include "grid_stream.h"
37 #include "SCF_statistics.h"
38 
39 
40 void output_sparsity(int n, const normalMatrix & M, const char* matrixName);
41 void output_sparsity_symm(int n, const symmMatrix & M, const char* matrixName);
42 void output_sparsity_triang(int n, const triangMatrix & M, const char* matrixName);
43 
44 int
45 compute_h_core_matrix_sparse(const IntegralInfo& integralInfo,
46  const Molecule& molecule,
47  const Molecule& extraCharges,
48  ergo_real electric_field_x,
49  ergo_real electric_field_y,
50  ergo_real electric_field_z,
51  const BasisInfoStruct& basisInfo,
52  symmMatrix & H_core_Matrix_sparse,
53  ergo_real threshold_integrals_1el,
54  int noOfThreadsForV,
55  mat::SizesAndBlocks const & matrix_size_block_info,
56  std::vector<int> const & permutationHML,
57  int const create_dipole_mtx = 0,
58  std::vector<int> const * const inversePermutationHML = 0,
59  std::string const * const calculation_identifier = 0,
60  std::string const * const method_and_basis_set = 0);
61 
62 int
64  const Molecule& molecule,
65  const BasisInfoStruct& basisInfo,
66  symmMatrix & H_core_Matrix_sparse,
67  ergo_real threshold_integrals_1el,
68  int noOfThreadsForV,
69  mat::SizesAndBlocks const & matrix_size_block_info,
70  std::vector<int> const & permutationHML);
71 
72 int
74  const Molecule& molecule,
75  const BasisInfoStruct& basisInfo,
76  const symmMatrix & D,
77  ergo_real threshold_integrals_1el,
78  mat::SizesAndBlocks const & matrix_size_block_info,
79  std::vector<int> const & permutationHML,
80  ergo_real* result_gradient_list);
81 
83  const BasisInfoStruct & basisInfo,
84  const char *name,
85  std::vector<int> const & inversePermutationHML);
86 
87 int
89  symmMatrix & A,
90  ergo_real disturbance,
91  int specificElementCount,
92  const int* elementIndexVector,
93  std::vector<int> const & permutationHML);
94 
95 int
97  int noOfElectrons,
98  symmMatrix & densityMatrix);
99 
100 int
102  const symmMatrix & M,
103  const char* fileName,
104  std::vector<int> const & permutationHML);
105 
106 int
108  symmMatrix & M,
109  const char* fileName,
110  std::vector<int> const & permutationHML);
111 
112 int
113 write_full_matrix(int n,
114  const symmMatrix & M,
115  const char* fileName,
116  std::vector<int> const & inversePermutationHML);
117 
118 int
120 
121 int
122 write_2el_integral_m_file(const BasisInfoStruct & basisInfo, const IntegralInfo & integralInfo);
123 
124 int
126  const BasisInfoStruct & basisInfoDensFit,
127  const Molecule& molecule,
128  const IntegralInfo& integralInfo,
129  symmMatrix & twoelMatrix_sparse,
130  symmMatrix & densityMatrix_sparse,
131  const JK::Params& J_K_params,
132  const JK::ExchWeights & CAM_params,
133  const Dft::GridParams& gridParams,
134  int do_xc,
135  ergo_real* energy_2el,
136  int noOfElectrons,
137  DensfitData* df_data,
138  mat::SizesAndBlocks const & matrix_size_block_info,
139 
140  std::vector<int> const & permutationHML,
141  std::vector<int> const & inversePermutationHML,
142  int get_J_K_Fxc_matrices,
143  symmMatrix & J_matrix,
144  symmMatrix & K_matrix,
145  symmMatrix & Fxc_matrix,
146  SCF_statistics & stats);
147 
148 int
150  const BasisInfoStruct & basisInfoDensFit,
151  const Molecule& molecule,
152  const IntegralInfo& integralInfo,
153  const JK::ExchWeights & CAM_params,
154  symmMatrix & twoelMatrix_sparse_alpha,
155  symmMatrix & twoelMatrix_sparse_beta,
156  symmMatrix & densityMatrix_sparse_alpha,
157  symmMatrix & densityMatrix_sparse_beta,
158  const JK::Params& J_K_params,
159  const Dft::GridParams& gridParams,
160  int do_xc,
161  ergo_real* energy_2el,
162  int noOfElectrons,
163  DensfitData* df_data,
164  mat::SizesAndBlocks const & matrix_size_block_info,
165  std::vector<int> const & permutationHML,
166  std::vector<int> const & inversePermutationHML);
167 
168 int
170  const BasisInfoStruct & basisInfoDensFit,
171  const Molecule& molecule,
172  const IntegralInfo& integralInfo,
173  const JK::ExchWeights & CAM_params,
174  symmMatrix & twoelMatrix_Fc,
175  symmMatrix & twoelMatrix_Fo,
176  symmMatrix & densityMatrix_sparse_alpha,
177  symmMatrix & densityMatrix_sparse_beta,
178  const JK::Params& J_K_params,
179  const Dft::GridParams& gridParams,
180  int do_xc,
181  ergo_real* energy_2el,
182  int noOfElectrons,
183  DensfitData* df_data,
184  mat::SizesAndBlocks const & matrix_size_block_info,
185  std::vector<int> const & permutationHML,
186  std::vector<int> const & inversePermutationHML);
187 
188 int
190  symmMatrix & F_symm,
191  symmMatrix & D_symm,
192  symmMatrix & S_symm,
193  normalMatrix & result,
194  ergo_real sparse_threshold);
195 
196 int
198  int alpha_beta_diff,
199  int* noOfElectrons_alpha,
200  int* noOfElectrons_beta);
201 
202 void
203 get_hf_weight_and_cam_params(int use_dft,
204  ergo_real* exch_param_alpha,
205  ergo_real* exch_param_beta,
206  ergo_real* exch_param_mu);
207 
208 int
210  int alpha_beta_diff,
211  int* noOfElectrons_alpha,
212  int* noOfElectrons_beta);
213 
214 void
215 get_dipole_moment(const symmMatrix & densityMatrix,
216  const BasisInfoStruct & basisInfo,
217  mat::SizesAndBlocks const & matrix_size_block_info,
218  std::vector<int> const & permutationHML,
219  const Molecule& molecule);
220 
221 void
222 do_mulliken_atomic_charges(const symmMatrix & densityMatrix,
223  const symmMatrix & S_symm,
224  const BasisInfoStruct & basisInfo,
225  mat::SizesAndBlocks const & matrix_size_block_info,
226  std::vector<int> const & permutationHML,
227  std::vector<int> const & inversePermutationHML,
228  const Molecule& molecule);
229 
230 void
231 do_mulliken_spin_densities(const symmMatrix & spinDensityMatrix,
232  const symmMatrix & S_symm,
233  const BasisInfoStruct & basisInfo,
234  mat::SizesAndBlocks const & matrix_size_block_info,
235  std::vector<int> const & permutationHML,
236  std::vector<int> const & inversePermutationHML,
237  const Molecule& molecule);
238 
239 void
240 do_density_images(const BasisInfoStruct & basisInfo,
241  const Molecule& molecule,
242  const ergo_real* densityMatrixFull_tot,
243  const ergo_real* densityMatrixFull_spin,
244  double output_density_images_boxwidth);
245 
246 void
247 do_acc_scan_J(const symmMatrix & D,
248  const IntegralInfo & integralInfo,
249  const BasisInfoStruct & basisInfo,
250  triangMatrix & invCholFactor,
251  bool doInvCholFactorTransformation,
252  const JK::Params & J_K_params,
253  mat::SizesAndBlocks const & matrix_size_block_info,
254  std::vector<int> const & permutationHML,
255  int nSteps,
256  ergo_real startThresh,
257  ergo_real stepFactor);
258 
259 void
261  const IntegralInfo & integralInfo,
262  const BasisInfoStruct & basisInfo,
263  triangMatrix & invCholFactor,
264  bool doInvCholFactorTransformation,
265  const JK::ExchWeights & CAM_params,
266  const JK::Params & J_K_params,
267  mat::SizesAndBlocks const & matrix_size_block_info,
268  std::vector<int> const & permutationHML,
269  std::vector<int> const & inversePermutationHML,
270  int nSteps,
271  ergo_real startThresh,
272  ergo_real stepFactor);
273 
274 void
276  const IntegralInfo & integralInfo,
277  const BasisInfoStruct & basisInfo,
278  const Molecule & molecule,
279  const Dft::GridParams & gridParams,
280  int noOfElectrons,
281  triangMatrix & invCholFactor,
282  bool doInvCholFactorTransformation,
283  mat::SizesAndBlocks const & matrix_size_block_info,
284  std::vector<int> const & permutationHML,
285  std::vector<int> const & inversePermutationHML,
286  int nSteps,
287  ergo_real startThresh,
288  ergo_real stepFactor);
289 
290 
291 #endif
#define A
void output_sparsity_symm(int n, const symmMatrix &M, const char *matrixName)
Definition: scf_utils.cc:373
double ergo_real
Definition: realtype.h:53
int get_simple_starting_guess_sparse(int n, int noOfElectrons, symmMatrix &densityMatrix)
Definition: scf_utils.cc:862
Normal matrix.
Definition: MatrixBase.h:47
A structure describing the grid settings.
Definition: grid_params.h:49
Streaming grid generator.
void do_density_images(const BasisInfoStruct &basisInfo, const Molecule &molecule, const ergo_real *densityMatrixFull_tot, const ergo_real *densityMatrixFull_spin, double output_density_images_boxwidth)
Definition: scf_utils.cc:2512
void do_mulliken_spin_densities(const symmMatrix &spinDensityMatrix, const symmMatrix &S_symm, const BasisInfoStruct &basisInfo, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML, const Molecule &molecule)
Definition: scf_utils.cc:2475
Definition: integrals_2el.h:36
Definition: SCF_statistics.h:47
int get_2e_matrix_and_energy_sparse(const BasisInfoStruct &basisInfo, const BasisInfoStruct &basisInfoDensFit, const Molecule &molecule, const IntegralInfo &integralInfo, symmMatrix &twoelMatrix_sparse, symmMatrix &densityMatrix_sparse, const JK::Params &J_K_params, const JK::ExchWeights &CAM_params, const Dft::GridParams &gridParams, int do_xc, ergo_real *energy_2el, int noOfElectrons, DensfitData *df_data, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML, int get_J_K_Fxc_matrices, symmMatrix &J_matrix, symmMatrix &K_matrix, symmMatrix &Fxc_matrix, SCF_statistics &stats)
General routine for computing the two-electron part of the Fock/KS matrix.
Definition: scf_utils.cc:1574
int get_diag_matrix_from_file(int n, symmMatrix &M, const char *fileName, std::vector< int > const &permutationHML)
Definition: scf_utils.cc:884
void output_sparsity_triang(int n, const triangMatrix &M, const char *matrixName)
Definition: scf_utils.cc:378
Describes dimensions of matrix and its blocks on all levels.
Definition: SizesAndBlocks.h:37
Representation of a molecule as a set of nuclei and total charge.
Definition: molecule.h:76
Definition: densityfitting.h:37
void output_sparsity(int n, const normalMatrix &M, const char *matrixName)
Definition: scf_utils.cc:368
Definition: integral_info.h:130
int write_diag_elements_to_file(int n, const symmMatrix &M, const char *fileName, std::vector< int > const &permutationHML)
Definition: scf_utils.cc:961
Contains coefficients needed for quick integral evaluation.
Definition: integral_info.h:81
int get_2e_matrices_and_energy_restricted_open(const BasisInfoStruct &basisInfo, const BasisInfoStruct &basisInfoDensFit, const Molecule &molecule, const IntegralInfo &integralInfo, const JK::ExchWeights &CAM_params, symmMatrix &twoelMatrix_Fc, symmMatrix &twoelMatrix_Fo, symmMatrix &densityMatrix_sparse_alpha, symmMatrix &densityMatrix_sparse_beta, const JK::Params &J_K_params, const Dft::GridParams &gridParams, int do_xc, ergo_real *energy_2el, int noOfElectrons, DensfitData *df_data, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML)
Computes G_c and G_o.
Definition: scf_utils.cc:2022
int write_basis_func_coord_file(const BasisInfoStruct &basisInfo)
Definition: scf_utils.cc:1020
Definition: basisinfo.h:111
Header file with typedefs for matrix and vector types.
int compute_h_core_matrix_simple_dense(const IntegralInfo &integralInfo, const Molecule &molecule, const BasisInfoStruct &basisInfo, symmMatrix &H_core_Matrix_sparse, ergo_real threshold_integrals_1el, int noOfThreadsForV, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML)
Definition: scf_utils.cc:459
void do_acc_scan_Vxc(symmMatrix &D, const IntegralInfo &integralInfo, const BasisInfoStruct &basisInfo, const Molecule &molecule, const Dft::GridParams &gridParams, int noOfElectrons, triangMatrix &invCholFactor, bool doInvCholFactorTransformation, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML, int nSteps, ergo_real startThresh, ergo_real stepFactor)
Definition: scf_utils.cc:317
int save_symmetric_matrix(symmMatrix &A, const BasisInfoStruct &basisInfo, const char *name, std::vector< int > const &inversePermutationHML)
Saves specified symmetic matrix to a file of specified name.
Definition: scf_utils.cc:779
int determine_number_of_electrons_unrestricted(int noOfElectrons, int alpha_beta_diff, int *noOfElectrons_alpha, int *noOfElectrons_beta)
Definition: scf_utils.cc:2303
int get_2e_matrices_and_energy_sparse_unrestricted(const BasisInfoStruct &basisInfo, const BasisInfoStruct &basisInfoDensFit, const Molecule &molecule, const IntegralInfo &integralInfo, const JK::ExchWeights &CAM_params, symmMatrix &twoelMatrix_sparse_alpha, symmMatrix &twoelMatrix_sparse_beta, symmMatrix &densityMatrix_sparse_alpha, symmMatrix &densityMatrix_sparse_beta, const JK::Params &J_K_params, const Dft::GridParams &gridParams, int do_xc, ergo_real *energy_2el, int noOfElectrons, DensfitData *df_data, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML)
Definition: scf_utils.cc:1803
void get_hf_weight_and_cam_params(int use_dft, ergo_real *exch_param_alpha, ergo_real *exch_param_beta, ergo_real *exch_param_mu)
Definition: scf_utils.cc:2708
void get_dipole_moment(const symmMatrix &densityMatrix, const BasisInfoStruct &basisInfo, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, const Molecule &molecule)
Definition: scf_utils.cc:2362
void do_acc_scan_J(const symmMatrix &D, const IntegralInfo &integralInfo, const BasisInfoStruct &basisInfo, triangMatrix &invCholFactor, bool doInvCholFactorTransformation, const JK::Params &J_K_params, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, int nSteps, ergo_real startThresh, ergo_real stepFactor)
Definition: scf_utils.cc:155
void do_acc_scan_K(symmMatrix &D, const IntegralInfo &integralInfo, const BasisInfoStruct &basisInfo, triangMatrix &invCholFactor, bool doInvCholFactorTransformation, const JK::ExchWeights &CAM_params, const JK::Params &J_K_params, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML, int nSteps, ergo_real startThresh, ergo_real stepFactor)
Definition: scf_utils.cc:232
int add_disturbance_to_matrix(int n, symmMatrix &A, ergo_real disturbance, int specificElementCount, const int *elementIndexVector, std::vector< int > const &permutationHML)
Definition: scf_utils.cc:805
void do_mulliken_atomic_charges(const symmMatrix &densityMatrix, const symmMatrix &S_symm, const BasisInfoStruct &basisInfo, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, std::vector< int > const &inversePermutationHML, const Molecule &molecule)
Definition: scf_utils.cc:2446
int compute_FDSminusSDF_sparse(int n, symmMatrix &F_symm, symmMatrix &D_symm, symmMatrix &S_symm, normalMatrix &result, ergo_real sparse_threshold)
Definition: scf_utils.cc:2186
Symmetric matrix.
Definition: MatrixBase.h:49
int get_gradient_for_given_mol_and_dens(const IntegralInfo &integralInfo, const Molecule &molecule, const BasisInfoStruct &basisInfo, const symmMatrix &D, ergo_real threshold_integrals_1el, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, ergo_real *result_gradient_list)
Definition: scf_utils.cc:745
int write_2el_integral_m_file(const BasisInfoStruct &basisInfo, const IntegralInfo &integralInfo)
Definition: scf_utils.cc:1039
int write_full_matrix(int n, const symmMatrix &M, const char *fileName, std::vector< int > const &inversePermutationHML)
Definition: scf_utils.cc:1003
int compute_h_core_matrix_sparse(const IntegralInfo &integralInfo, const Molecule &molecule, const Molecule &extraCharges, ergo_real electric_field_x, ergo_real electric_field_y, ergo_real electric_field_z, const BasisInfoStruct &basisInfo, symmMatrix &H_core_Matrix_sparse, ergo_real threshold_integrals_1el, int noOfThreadsForV, mat::SizesAndBlocks const &matrix_size_block_info, std::vector< int > const &permutationHML, int const create_dipole_mtx=0, std::vector< int > const *const inversePermutationHML=0, std::string const *const calculation_identifier=0, std::string const *const method_and_basis_set=0)
Definition: scf_utils.cc:505