Psi4
dfhelper.h
Go to the documentation of this file.
1 /*
2  * @BEGIN LICENSE
3  *
4  * Psi4: an open-source quantum chemistry software package
5  *
6  * Copyright (c) 2007-2018 The Psi4 Developers.
7  *
8  * The copyrights for code used from other parties are included in
9  * the corresponding files.
10  *
11  * This file is part of Psi4.
12  *
13  * Psi4 is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU Lesser General Public License as published by
15  * the Free Software Foundation, version 3.
16  *
17  * Psi4 is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU Lesser General Public License for more details.
21  *
22  * You should have received a copy of the GNU Lesser General Public License along
23  * with Psi4; if not, write to the Free Software Foundation, Inc.,
24  * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
25  *
26  * @END LICENSE
27  */
28 
29 #ifndef three_index_dfhelper
30 #define three_index_dfhelper
31 
32 #include "psi4/psi4-dec.h"
33 #include <psi4/libmints/typedefs.h>
34 
35 #include <map>
36 #include <list>
37 #include <vector>
38 #include <tuple>
39 #include <string>
40 
41 namespace psi {
42 
43 class BasisSet;
44 class Options;
45 class Matrix;
46 class ERISieve;
47 class TwoBodyAOInt;
48 
50  public:
51  DFHelper(std::shared_ptr<BasisSet> primary, std::shared_ptr<BasisSet> aux);
52  ~DFHelper();
53 
62  void set_method(std::string method) { method_ = method; }
63  std::string get_method() { return method_; }
64 
70  void set_nthreads(size_t nthreads) { nthreads_ = nthreads; }
71  size_t get_nthreads() { return nthreads_; }
72 
78  void set_memory(size_t doubles) { memory_ = doubles; }
79  size_t get_memory() { return memory_; }
80 
82  size_t get_AO_size() { return big_skips_[nao_]; }
83 
85  double ao_sparsity() { return (1.0 - (double)small_skips_[nao_] / (double)(nao_ * nao_)); }
86 
95  void set_AO_core(bool core) { AO_core_ = core; }
96  bool get_AO_core() { return AO_core_; }
97 
104  void set_MO_core(bool core) { MO_core_ = core; }
105  bool get_MO_core() { return MO_core_; }
106 
108  void set_schwarz_cutoff(double cutoff) { cutoff_ = cutoff; }
109  double get_schwarz_cutoff() { return cutoff_; }
110 
112  void set_metric_pow(double pow) { mpower_ = pow; }
113  double get_metric_pow() { return mpower_; }
114 
120  void hold_met(bool hold) { hold_met_ = hold; }
121  bool get_hold_met() { return hold_met_; }
122 
127  void set_fitting_condition(double condition) { condition_ = condition; }
128  bool get_fitting_condition() { return condition_; }
129 
134  void set_do_wK(bool do_wK) { do_wK_ = do_wK; }
135  size_t get_do_wK() { return do_wK_; }
136 
141  void set_omega(double omega) { omega_ = omega; }
142  size_t get_omega() { return omega_; }
143 
148  void set_print_lvl(int print_lvl) { print_lvl_ = print_lvl; }
149 
151  void initialize();
152 
154  void print_header();
155 
161  void add_space(std::string key, SharedMatrix space);
162 
170  void add_transformation(std::string name, std::string key1, std::string key2, std::string order = "Qpq");
171 
173  void transform();
174 
175  // => Tensor IO <=
176  // many ways to access the 3-index tensors.
177 
187  void fill_tensor(std::string name, SharedMatrix M);
188  void fill_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1);
189  void fill_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2);
190  void fill_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2,
191  std::vector<size_t> a3);
192 
199  void fill_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2,
200  std::vector<size_t> a3);
201  void fill_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2);
202  void fill_tensor(std::string name, double* b, std::vector<size_t> a1);
203  void fill_tensor(std::string name, double* b);
204 
212  SharedMatrix get_tensor(std::string name);
213  SharedMatrix get_tensor(std::string name, std::vector<size_t> a1);
214  SharedMatrix get_tensor(std::string name, std::vector<size_t> a1, std::vector<size_t> a2);
215  SharedMatrix get_tensor(std::string name, std::vector<size_t> a1, std::vector<size_t> a2, std::vector<size_t> a3);
216 
222  void add_disk_tensor(std::string key, std::tuple<size_t, size_t, size_t> dimensions);
223 
231  void write_disk_tensor(std::string name, SharedMatrix M);
232  void write_disk_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1);
233  void write_disk_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2);
234  void write_disk_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2,
235  std::vector<size_t> a3);
236 
244  void write_disk_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2,
245  std::vector<size_t> a3);
246  void write_disk_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2);
247  void write_disk_tensor(std::string name, double* b, std::vector<size_t> a1);
248  void write_disk_tensor(std::string name, double* b);
249 
251  void transpose(std::string name, std::tuple<size_t, size_t, size_t> order);
252 
254  void clear_spaces();
255 
257  void clear_all();
258 
260  size_t get_space_size(std::string key);
261  size_t get_tensor_size(std::string key);
262  std::tuple<size_t, size_t, size_t> get_tensor_shape(std::string key);
263  size_t get_naux() { return naux_; }
264 
266  void build_JK(std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright,
267  std::vector<SharedMatrix> D, std::vector<SharedMatrix> J,
268  std::vector<SharedMatrix> K, size_t max_nocc,
269  bool do_J, bool do_K, bool do_wK, bool lr_symmetric);
270 
271  protected:
272  // => basis sets <=
273  std::shared_ptr<BasisSet> primary_;
274  std::shared_ptr<BasisSet> aux_;
275  size_t nao_;
276  size_t naux_;
277 
278  // => memory in doubles <=
279  size_t memory_ = 256000000;
280 
281  // => internal holders <=
282  std::string method_ = "STORE";
283  bool direct_;
286  bool AO_core_ = true;
287  bool MO_core_ = false;
288  size_t nthreads_ = 1;
289  double cutoff_ = 1e-12;
290  double condition_ = 1e-12;
291  double mpower_ = -0.5;
292  bool hold_met_ = false;
293  bool built_ = false;
294  bool transformed_ = false;
295  std::pair<size_t, size_t> info_;
296  bool ordered_ = false;
297  bool do_wK_ = false;
298  double omega_;
299  bool debug_ = false;
300  int print_lvl_ = 1;
301 
302  // => in-core machinery <=
303  void AO_core();
304  std::unique_ptr<double[]> Ppq_;
305  std::map<double, SharedMatrix> metrics_;
306 
307  // => AO building machinery <=
308  void prepare_AO();
309  void prepare_AO_core();
310  void compute_dense_Qpq_blocking_Q(const size_t start, const size_t stop, double* Mp,
311  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
312  void compute_sparse_pQq_blocking_Q(const size_t start, const size_t stop, double* Mp,
313  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
314  void compute_sparse_pQq_blocking_p(const size_t start, const size_t stop, double* Mp,
315  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
316  void compute_sparse_pQq_blocking_p_symm(const size_t start, const size_t stop, double* Mp,
317  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
318  void contract_metric_AO_core_symm(double* Qpq, double* metp, size_t begin, size_t end);
319  void grab_AO(const size_t start, const size_t stop, double* Mp);
320 
321  // first integral transforms
322  void first_transform_pQq(size_t nao, size_t naux, size_t bsize, size_t bcount, size_t block_size,
323  double* Mp, double* Tp, double* Bp, std::vector<std::vector<double>>& C_buffers);
324 
325 
326  // => index vectors for screened AOs <=
327  std::vector<size_t> small_skips_;
328  std::vector<size_t> big_skips_;
329  std::vector<size_t> symm_ignored_columns_;
330  std::vector<size_t> symm_small_skips_;
331  std::vector<size_t> symm_big_skips_;
332 
333  // => shell info and blocking <=
334  size_t pshells_;
335  size_t Qshells_;
336  double Qshell_max_;
337  std::vector<size_t> pshell_aggs_;
338  std::vector<size_t> Qshell_aggs_;
339  void prepare_blocking();
340 
341  // => generalized blocking <=
342  std::pair<size_t, size_t> pshell_blocks_for_AO_build(const size_t mem, size_t symm,
343  std::vector<std::pair<size_t, size_t>>& b);
344  std::pair<size_t, size_t> Qshell_blocks_for_transform(const size_t mem, size_t wtmp, size_t wfinal,
345  std::vector<std::pair<size_t, size_t>>& b);
346  void metric_contraction_blocking(std::vector<std::pair<size_t, size_t>>& steps,
347  size_t blocking_index, size_t block_sizes, size_t total_mem, size_t memory_factor, size_t memory_bump);
348 
349  // => Schwarz Screening <=
350  std::vector<size_t> schwarz_fun_mask_;
351  std::vector<size_t> schwarz_shell_mask_;
352  std::vector<size_t> schwarz_fun_count_;
353  void prepare_sparsity();
354 
355  // => Coulomb metric handling <=
356  std::vector<std::pair<double, std::string>> metric_keys_;
357  void prepare_metric();
358  void prepare_metric_core();
359  double* metric_prep_core(double pow);
360  std::string return_metfile(double pow);
361  std::string compute_metric(double pow);
362 
363  // => metric operations <=
364  void contract_metric_Qpq(std::string file, double* metp, double* Mp, double* Fp, const size_t tots);
365  void contract_metric(std::string file, double* metp, double* Mp, double* Fp, const size_t tots);
366  void contract_metric_core(std::string file);
367  void contract_metric_AO(double* Mp);
368  void contract_metric_AO_core(double* Qpq, double* metp);
369 
370  // => spaces and transformation maps <=
371  std::map<std::string, std::tuple<SharedMatrix, size_t>> spaces_;
372  std::map<std::string, std::tuple<std::string, std::string, size_t>> transf_;
373  std::map<std::string, std::unique_ptr<double[]>> transf_core_;
374 
375  // => transformation machinery <=
376  std::pair<size_t, size_t> identify_order();
377  void print_order();
378  void put_transformations_Qpq(int naux, int begin, int end,
379  int wsize, int bsize, double* Fp, int ind, bool bleft);
380  void put_transformations_pQq(int naux, int begin, int end, int block_size, int bcount,
381  int wsize, int bsize, double* Np, double* Fp, int ind, bool bleft);
382  std::vector<std::pair<std::string, size_t>> sorted_spaces_;
383  std::vector<std::string> order_;
384  std::vector<std::string> bspace_;
385  std::vector<size_t> strides_;
386 
387  // => FILE IO maintenence <=
388  typedef struct StreamStruct {
389 
390  StreamStruct();
391  StreamStruct(std::string filename, std::string op, bool activate=true);
392  ~StreamStruct();
393 
394  FILE* get_stream(std::string op);
395  void change_stream(std::string op);
396  void close_stream();
397 
398  FILE* fp_;
399  std::string op_;
400  bool open_ = false;
401  std::string filename_;
402 
403  } Stream;
404 
405  std::map<std::string, std::shared_ptr<Stream>> file_streams_;
406  FILE* stream_check(std::string filename, std::string op);
407 
408  // => FILE IO machinery <=
409  void put_tensor(std::string file, double* b, std::pair<size_t, size_t> a1, std::pair<size_t, size_t> a2,
410  std::pair<size_t, size_t> a3, std::string op);
411 
412  void put_tensor(std::string file, double* b, const size_t start1, const size_t stop1, const size_t start2,
413  const size_t stop2, std::string op);
414  void get_tensor_(std::string file, double* b, std::pair<size_t, size_t> a1, std::pair<size_t, size_t> a2,
415  std::pair<size_t, size_t> a3);
416  void get_tensor_(std::string file, double* b, const size_t start1, const size_t stop1, const size_t start2,
417  const size_t stop2);
418  void put_tensor_AO(std::string file, double* Mp, size_t size, size_t start, std::string op);
419  void get_tensor_AO(std::string file, double* Mp, size_t size, size_t start);
420 
421  // => internal handlers for FILE IO <=
422  std::map<std::string, std::tuple<std::string, std::string>> files_;
423  std::map<std::string, std::tuple<size_t, size_t, size_t>> sizes_;
424  std::map<std::string, std::tuple<size_t, size_t, size_t>> tsizes_;
425  std::map<std::string, std::string> AO_files_;
426  std::vector<size_t> AO_file_sizes_;
427  std::vector<std::string> AO_names_;
428  std::string start_filename(std::string start);
429  void filename_maker(std::string name, size_t a0, size_t a1, size_t a2, size_t op=0);
430  void AO_filename_maker(size_t i);
431  void check_file_key(std::string);
432  void check_file_tuple(std::string name, std::pair<size_t, size_t> t0, std::pair<size_t, size_t> t1,
433  std::pair<size_t, size_t> t2);
434  void check_matrix_size(std::string name, SharedMatrix M, std::pair<size_t, size_t> t0, std::pair<size_t, size_t> t1,
435  std::pair<size_t, size_t> t2);
436 
437  // => transpose a tensor <=
438  void transpose_core(std::string name, std::tuple<size_t, size_t, size_t> order);
439  void transpose_disk(std::string name, std::tuple<size_t, size_t, size_t> order);
440 
441  // => JK <=
442  void compute_JK(std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright,
443  std::vector<SharedMatrix> D, std::vector<SharedMatrix> J,
444  std::vector<SharedMatrix> K, size_t max_nocc,
445  bool do_J, bool do_K, bool do_wK, bool lr_symmetric);
446  void compute_D(std::vector<SharedMatrix> D, std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright);
447  void compute_J(std::vector<SharedMatrix> D, std::vector<SharedMatrix> J, double* Mp, double* T1p, double* T2p,
448  std::vector<std::vector<double>>& D_buffers, size_t bcount, size_t block_size);
449  void compute_J_symm(std::vector<SharedMatrix> D, std::vector<SharedMatrix> J, double* Mp, double* T1p, double* T2p,
450  std::vector<std::vector<double>>& D_buffers, size_t bcount, size_t block_size);
451  void compute_K(std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright, std::vector<SharedMatrix> K,
452  double* Tp, double* Jtmp, double* Mp, size_t bcount, size_t block_size,
453  std::vector<std::vector<double>>& C_buffers, bool lr_symmetric);
454  std::tuple<size_t, size_t> Qshell_blocks_for_JK_build(
455  std::vector<std::pair<size_t, size_t>>& b, size_t max_nocc, bool lr_symmetric);
456 
457  // => misc <=
458  void fill(double* b, size_t count, double value);
459 
460 }; // End DF Helper class
461 } // psi4 namespace
462 #endif
std::map< std::string, std::tuple< SharedMatrix, size_t > > spaces_
Definition: dfhelper.h:371
void hold_met(bool hold)
Definition: dfhelper.h:120
std::vector< size_t > small_skips_
Definition: dfhelper.h:327
double ao_sparsity()
Returns the amount of sparsity in the AO integrals.
Definition: dfhelper.h:85
std::vector< size_t > AO_file_sizes_
Definition: dfhelper.h:426
size_t pshells_
Definition: dfhelper.h:334
void set_metric_pow(double pow)
fitting metric power (defaults to -0.5)
Definition: dfhelper.h:112
double get_metric_pow()
Definition: dfhelper.h:113
std::map< double, SharedMatrix > metrics_
Definition: dfhelper.h:305
std::vector< std::string > order_
Definition: dfhelper.h:383
void set_nthreads(size_t nthreads)
Definition: dfhelper.h:70
size_t get_do_wK()
Definition: dfhelper.h:135
void set_AO_core(bool core)
Definition: dfhelper.h:95
std::string get_method()
Definition: dfhelper.h:63
std::vector< size_t > pshell_aggs_
Definition: dfhelper.h:337
bool direct_
Definition: dfhelper.h:283
double Qshell_max_
Definition: dfhelper.h:336
std::map< std::string, std::tuple< size_t, size_t, size_t > > tsizes_
Definition: dfhelper.h:424
double get_schwarz_cutoff()
Definition: dfhelper.h:109
void set_method(std::string method)
Definition: dfhelper.h:62
Definition: pointgrp.h:104
void set_do_wK(bool do_wK)
Definition: dfhelper.h:134
bool get_hold_met()
Definition: dfhelper.h:121
std::map< std::string, std::tuple< std::string, std::string > > files_
Definition: dfhelper.h:422
int nao
Definition: dx_write.cc:56
FILE * fp_
Definition: dfhelper.h:398
size_t Qshells_
Definition: dfhelper.h:335
void set_fitting_condition(double condition)
Definition: dfhelper.h:127
size_t get_AO_size()
Returns the number of doubles in the screened AO integrals.
Definition: dfhelper.h:82
std::vector< std::pair< std::string, size_t > > sorted_spaces_
Definition: dfhelper.h:382
void set_MO_core(bool core)
Definition: dfhelper.h:104
size_t nao_
Definition: dfhelper.h:275
std::vector< std::string > bspace_
Definition: dfhelper.h:384
std::map< std::string, std::shared_ptr< Stream > > file_streams_
Definition: dfhelper.h:405
std::vector< size_t > Qshell_aggs_
Definition: dfhelper.h:338
size_t naux_
Definition: dfhelper.h:276
bool direct_iaQ_
Definition: dfhelper.h:284
std::vector< size_t > schwarz_fun_count_
Definition: dfhelper.h:352
std::string filename_
Definition: dfhelper.h:401
std::unique_ptr< double[]> Ppq_
Definition: dfhelper.h:304
std::vector< size_t > big_skips_
Definition: dfhelper.h:328
std::map< std::string, std::tuple< size_t, size_t, size_t > > sizes_
Definition: dfhelper.h:423
std::string op_
Definition: dfhelper.h:399
Definition: dfhelper.h:49
Definition: dfhelper.h:388
size_t get_nthreads()
Definition: dfhelper.h:71
bool get_MO_core()
Definition: dfhelper.h:105
std::vector< size_t > strides_
Definition: dfhelper.h:385
std::map< std::string, std::unique_ptr< double[]> > transf_core_
Definition: dfhelper.h:373
std::vector< size_t > symm_ignored_columns_
Definition: dfhelper.h:329
std::vector< std::string > AO_names_
Definition: dfhelper.h:427
std::vector< size_t > symm_small_skips_
Definition: dfhelper.h:330
std::shared_ptr< BasisSet > primary_
Definition: dfhelper.h:273
void set_schwarz_cutoff(double cutoff)
schwarz screening cutoff (defaults to 1e-12)
Definition: dfhelper.h:108
void set_memory(size_t doubles)
Definition: dfhelper.h:78
std::map< std::string, std::tuple< std::string, std::string, size_t > > transf_
Definition: dfhelper.h:372
std::shared_ptr< BasisSet > aux_
Definition: dfhelper.h:274
void set_print_lvl(int print_lvl)
Definition: dfhelper.h:148
void set_omega(double omega)
Definition: dfhelper.h:141
std::shared_ptr< Matrix > SharedMatrix
Definition: adc.h:49
#define PSI_API
Definition: pragma.h:155
double omega_
Definition: dfhelper.h:298
bool get_AO_core()
Definition: dfhelper.h:96
std::vector< size_t > schwarz_shell_mask_
Definition: dfhelper.h:351
std::pair< size_t, size_t > info_
Definition: dfhelper.h:295
size_t get_omega()
Definition: dfhelper.h:142
std::vector< size_t > symm_big_skips_
Definition: dfhelper.h:331
size_t get_memory()
Definition: dfhelper.h:79
std::vector< std::pair< double, std::string > > metric_keys_
Definition: dfhelper.h:356
std::vector< size_t > schwarz_fun_mask_
Definition: dfhelper.h:350
size_t get_naux()
Definition: dfhelper.h:263
bool space(char c)
Definition: stl_string.cc:109
bool get_fitting_condition()
Definition: dfhelper.h:128
std::map< std::string, std::string > AO_files_
Definition: dfhelper.h:425
bool symm_compute_
Definition: dfhelper.h:285