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 
145  void initialize();
146 
148  void print_header();
149 
155  void add_space(std::string key, SharedMatrix space);
156 
164  void add_transformation(std::string name, std::string key1, std::string key2, std::string order = "Qpq");
165 
167  void transform();
168 
169  // => Tensor IO <=
170  // many ways to access the 3-index tensors.
171 
181  void fill_tensor(std::string name, SharedMatrix M);
182  void fill_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1);
183  void fill_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2);
184  void fill_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2,
185  std::vector<size_t> a3);
186 
193  void fill_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2,
194  std::vector<size_t> a3);
195  void fill_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2);
196  void fill_tensor(std::string name, double* b, std::vector<size_t> a1);
197  void fill_tensor(std::string name, double* b);
198 
206  SharedMatrix get_tensor(std::string name);
207  SharedMatrix get_tensor(std::string name, std::vector<size_t> a1);
208  SharedMatrix get_tensor(std::string name, std::vector<size_t> a1, std::vector<size_t> a2);
209  SharedMatrix get_tensor(std::string name, std::vector<size_t> a1, std::vector<size_t> a2, std::vector<size_t> a3);
210 
216  void add_disk_tensor(std::string key, std::tuple<size_t, size_t, size_t> dimensions);
217 
225  void write_disk_tensor(std::string name, SharedMatrix M);
226  void write_disk_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1);
227  void write_disk_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2);
228  void write_disk_tensor(std::string name, SharedMatrix M, std::vector<size_t> a1, std::vector<size_t> a2,
229  std::vector<size_t> a3);
230 
238  void write_disk_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2,
239  std::vector<size_t> a3);
240  void write_disk_tensor(std::string name, double* b, std::vector<size_t> a1, std::vector<size_t> a2);
241  void write_disk_tensor(std::string name, double* b, std::vector<size_t> a1);
242  void write_disk_tensor(std::string name, double* b);
243 
245  void transpose(std::string name, std::tuple<size_t, size_t, size_t> order);
246 
248  void clear_spaces();
249 
251  void clear_all();
252 
254  size_t get_space_size(std::string key);
255  size_t get_tensor_size(std::string key);
256  std::tuple<size_t, size_t, size_t> get_tensor_shape(std::string key);
257  size_t get_naux() { return naux_; }
258 
260  void build_JK(std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright,
261  std::vector<SharedMatrix> D, std::vector<SharedMatrix> J,
262  std::vector<SharedMatrix> K, size_t max_nocc,
263  bool do_J, bool do_K, bool do_wK, bool lr_symmetric);
264 
265  protected:
266  // => basis sets <=
267  std::shared_ptr<BasisSet> primary_;
268  std::shared_ptr<BasisSet> aux_;
269  size_t nao_;
270  size_t naux_;
271 
272  // => memory in doubles <=
273  size_t memory_ = 256000000;
274 
275  // => internal holders <=
276  std::string method_ = "STORE";
277  bool direct_;
280  bool AO_core_ = true;
281  bool MO_core_ = false;
282  size_t nthreads_ = 1;
283  double cutoff_ = 1e-12;
284  double condition_ = 1e-12;
285  double mpower_ = -0.5;
286  bool hold_met_ = false;
287  bool built_ = false;
288  bool transformed_ = false;
289  std::pair<size_t, size_t> info_;
290  bool ordered_ = false;
291  bool do_wK_ = false;
292  double omega_;
293  bool debug_ = false;
294 
295  // => in-core machinery <=
296  void AO_core();
297  std::vector<double> Ppq_;
298  std::map<double, SharedMatrix> metrics_;
299 
300  // => AO building machinery <=
301  void prepare_AO();
302  void prepare_AO_core();
303  void compute_dense_Qpq_blocking_Q(const size_t start, const size_t stop, double* Mp,
304  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
305  void compute_sparse_pQq_blocking_Q(const size_t start, const size_t stop, double* Mp,
306  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
307  void compute_sparse_pQq_blocking_p(const size_t start, const size_t stop, double* Mp,
308  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
309  void compute_sparse_pQq_blocking_p_symm(const size_t start, const size_t stop, double* Mp,
310  std::vector<std::shared_ptr<TwoBodyAOInt>> eri);
311  void contract_metric_AO_core_symm(double* Qpq, double* metp, size_t begin, size_t end);
312  void grab_AO(const size_t start, const size_t stop, double* Mp);
313 
314  // first integral transforms
315  void first_transform_pQq(size_t nao, size_t naux, size_t bsize, size_t bcount, size_t block_size,
316  double* Mp, double* Tp, double* Bp, std::vector<std::vector<double>>& C_buffers);
317 
318 
319  // => index vectors for screened AOs <=
320  std::vector<size_t> small_skips_;
321  std::vector<size_t> big_skips_;
322  std::vector<size_t> symm_ignored_columns_;
323  std::vector<size_t> symm_small_skips_;
324  std::vector<size_t> symm_big_skips_;
325 
326  // => shell info and blocking <=
327  size_t pshells_;
328  size_t Qshells_;
329  double Qshell_max_;
330  std::vector<size_t> pshell_aggs_;
331  std::vector<size_t> Qshell_aggs_;
332  void prepare_blocking();
333 
334  // => generalized blocking <=
335  std::pair<size_t, size_t> pshell_blocks_for_AO_build(const size_t mem, size_t symm,
336  std::vector<std::pair<size_t, size_t>>& b);
337  std::pair<size_t, size_t> Qshell_blocks_for_transform(const size_t mem, size_t wtmp, size_t wfinal,
338  std::vector<std::pair<size_t, size_t>>& b);
339  void metric_contraction_blocking(std::vector<std::pair<size_t, size_t>>& steps,
340  size_t blocking_index, size_t block_sizes, size_t total_mem, size_t memory_factor, size_t memory_bump);
341 
342  // => Schwarz Screening <=
343  std::vector<size_t> schwarz_fun_mask_;
344  std::vector<size_t> schwarz_shell_mask_;
345  std::vector<size_t> schwarz_fun_count_;
346  void prepare_sparsity();
347 
348  // => Coulomb metric handling <=
349  std::vector<std::pair<double, std::string>> metric_keys_;
350  void prepare_metric();
351  void prepare_metric_core();
352  double* metric_prep_core(double pow);
353  std::string return_metfile(double pow);
354  std::string compute_metric(double pow);
355 
356  // => metric operations <=
357  void contract_metric_Qpq(std::string file, double* metp, double* Mp, double* Fp, const size_t tots);
358  void contract_metric(std::string file, double* metp, double* Mp, double* Fp, const size_t tots);
359  void contract_metric_core(std::string file);
360  void contract_metric_AO(double* Mp);
361  void contract_metric_AO_core(double* Qpq, double* metp);
362 
363  // => spaces and transformation maps <=
364  std::map<std::string, std::tuple<SharedMatrix, size_t>> spaces_;
365  std::map<std::string, std::tuple<std::string, std::string, size_t>> transf_;
366  std::map<std::string, std::vector<double>> transf_core_;
367 
368  // => transformation machinery <=
369  std::pair<size_t, size_t> identify_order();
370  void print_order();
371  void put_transformations_Qpq(int naux, int begin, int end,
372  int wsize, int bsize, double* Fp, int ind, bool bleft);
373  void put_transformations_pQq(int naux, int begin, int end, int block_size, int bcount,
374  int wsize, int bsize, double* Np, double* Fp, int ind, bool bleft);
375  std::vector<std::pair<std::string, size_t>> sorted_spaces_;
376  std::vector<std::string> order_;
377  std::vector<std::string> bspace_;
378  std::vector<size_t> strides_;
379 
380  // => FILE IO maintenence <=
381  typedef struct StreamStruct {
382 
383  StreamStruct();
384  StreamStruct(std::string filename, std::string op, bool activate=true);
385  ~StreamStruct();
386 
387  FILE* get_stream(std::string op);
388  void change_stream(std::string op);
389  void close_stream();
390 
391  FILE* fp_;
392  std::string op_;
393  bool open_ = false;
394  std::string filename_;
395 
396  } Stream;
397 
398  std::map<std::string, std::shared_ptr<Stream>> file_streams_;
399  FILE* stream_check(std::string filename, std::string op);
400 
401  // => FILE IO machinery <=
402  void put_tensor(std::string file, double* b, std::pair<size_t, size_t> a1, std::pair<size_t, size_t> a2,
403  std::pair<size_t, size_t> a3, std::string op);
404 
405  void put_tensor(std::string file, double* b, const size_t start1, const size_t stop1, const size_t start2,
406  const size_t stop2, std::string op);
407  void get_tensor_(std::string file, double* b, std::pair<size_t, size_t> a1, std::pair<size_t, size_t> a2,
408  std::pair<size_t, size_t> a3);
409  void get_tensor_(std::string file, double* b, const size_t start1, const size_t stop1, const size_t start2,
410  const size_t stop2);
411  void put_tensor_AO(std::string file, double* Mp, size_t size, size_t start, std::string op);
412  void get_tensor_AO(std::string file, double* Mp, size_t size, size_t start);
413 
414  // => internal handlers for FILE IO <=
415  std::map<std::string, std::tuple<std::string, std::string>> files_;
416  std::map<std::string, std::tuple<size_t, size_t, size_t>> sizes_;
417  std::map<std::string, std::tuple<size_t, size_t, size_t>> tsizes_;
418  std::map<std::string, std::string> AO_files_;
419  std::vector<size_t> AO_file_sizes_;
420  std::vector<std::string> AO_names_;
421  std::string start_filename(std::string start);
422  void filename_maker(std::string name, size_t a0, size_t a1, size_t a2, size_t op=0);
423  void AO_filename_maker(size_t i);
424  void check_file_key(std::string);
425  void check_file_tuple(std::string name, std::pair<size_t, size_t> t0, std::pair<size_t, size_t> t1,
426  std::pair<size_t, size_t> t2);
427  void check_matrix_size(std::string name, SharedMatrix M, std::pair<size_t, size_t> t0, std::pair<size_t, size_t> t1,
428  std::pair<size_t, size_t> t2);
429 
430  // => transpose a tensor <=
431  void transpose_core(std::string name, std::tuple<size_t, size_t, size_t> order);
432  void transpose_disk(std::string name, std::tuple<size_t, size_t, size_t> order);
433 
434  // => JK <=
435  void compute_JK(std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright,
436  std::vector<SharedMatrix> D, std::vector<SharedMatrix> J,
437  std::vector<SharedMatrix> K, size_t max_nocc,
438  bool do_J, bool do_K, bool do_wK, bool lr_symmetric);
439  void compute_D(std::vector<SharedMatrix> D, std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright);
440  void compute_J(std::vector<SharedMatrix> D, std::vector<SharedMatrix> J, double* Mp, double* T1p, double* T2p,
441  std::vector<std::vector<double>>& D_buffers, size_t bcount, size_t block_size);
442  void compute_J_symm(std::vector<SharedMatrix> D, std::vector<SharedMatrix> J, double* Mp, double* T1p, double* T2p,
443  std::vector<std::vector<double>>& D_buffers, size_t bcount, size_t block_size);
444  void compute_K(std::vector<SharedMatrix> Cleft, std::vector<SharedMatrix> Cright, std::vector<SharedMatrix> K,
445  double* Tp, double* Jtmp, double* Mp, size_t bcount, size_t block_size,
446  std::vector<std::vector<double>>& C_buffers, bool lr_symmetric);
447  std::tuple<size_t, size_t> Qshell_blocks_for_JK_build(
448  std::vector<std::pair<size_t, size_t>>& b, size_t max_nocc, bool lr_symmetric);
449 
450  // => misc <=
451  void fill(double* b, size_t count, double value);
452 
453 }; // End DF Helper class
454 } // psi4 namespace
455 #endif
std::map< std::string, std::tuple< SharedMatrix, size_t > > spaces_
Definition: dfhelper.h:364
void hold_met(bool hold)
Definition: dfhelper.h:120
std::vector< size_t > small_skips_
Definition: dfhelper.h:320
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:419
size_t pshells_
Definition: dfhelper.h:327
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:298
std::vector< std::string > order_
Definition: dfhelper.h:376
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:330
bool direct_
Definition: dfhelper.h:277
double Qshell_max_
Definition: dfhelper.h:329
std::map< std::string, std::tuple< size_t, size_t, size_t > > tsizes_
Definition: dfhelper.h:417
double get_schwarz_cutoff()
Definition: dfhelper.h:109
void set_method(std::string method)
Definition: dfhelper.h:62
Definition: pointgrp.h:106
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:415
int nao
Definition: dx_write.cc:56
FILE * fp_
Definition: dfhelper.h:391
size_t Qshells_
Definition: dfhelper.h:328
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:375
void set_MO_core(bool core)
Definition: dfhelper.h:104
size_t nao_
Definition: dfhelper.h:269
std::vector< double > Ppq_
Definition: dfhelper.h:297
std::vector< std::string > bspace_
Definition: dfhelper.h:377
std::map< std::string, std::shared_ptr< Stream > > file_streams_
Definition: dfhelper.h:398
std::vector< size_t > Qshell_aggs_
Definition: dfhelper.h:331
size_t naux_
Definition: dfhelper.h:270
bool direct_iaQ_
Definition: dfhelper.h:278
std::vector< size_t > schwarz_fun_count_
Definition: dfhelper.h:345
std::string filename_
Definition: dfhelper.h:394
std::vector< size_t > big_skips_
Definition: dfhelper.h:321
std::map< std::string, std::tuple< size_t, size_t, size_t > > sizes_
Definition: dfhelper.h:416
std::string op_
Definition: dfhelper.h:392
Definition: dfhelper.h:49
Definition: dfhelper.h:381
size_t get_nthreads()
Definition: dfhelper.h:71
bool get_MO_core()
Definition: dfhelper.h:105
std::vector< size_t > strides_
Definition: dfhelper.h:378
std::vector< size_t > symm_ignored_columns_
Definition: dfhelper.h:322
std::vector< std::string > AO_names_
Definition: dfhelper.h:420
std::vector< size_t > symm_small_skips_
Definition: dfhelper.h:323
std::shared_ptr< BasisSet > primary_
Definition: dfhelper.h:267
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:365
std::shared_ptr< BasisSet > aux_
Definition: dfhelper.h:268
void set_omega(double omega)
Definition: dfhelper.h:141
std::shared_ptr< Matrix > SharedMatrix
Definition: adc.h:49
std::map< std::string, std::vector< double > > transf_core_
Definition: dfhelper.h:366
#define PSI_API
Definition: pragma.h:128
double omega_
Definition: dfhelper.h:292
bool get_AO_core()
Definition: dfhelper.h:96
std::vector< size_t > schwarz_shell_mask_
Definition: dfhelper.h:344
std::pair< size_t, size_t > info_
Definition: dfhelper.h:289
size_t get_omega()
Definition: dfhelper.h:142
std::vector< size_t > symm_big_skips_
Definition: dfhelper.h:324
size_t get_memory()
Definition: dfhelper.h:79
std::vector< std::pair< double, std::string > > metric_keys_
Definition: dfhelper.h:349
std::vector< size_t > schwarz_fun_mask_
Definition: dfhelper.h:343
size_t get_naux()
Definition: dfhelper.h:257
bool space(char c)
Definition: stl_string.cc:119
bool get_fitting_condition()
Definition: dfhelper.h:128
std::map< std::string, std::string > AO_files_
Definition: dfhelper.h:418
bool symm_compute_
Definition: dfhelper.h:279