MPQC  3.0.0-alpha
df_runtime.h
1 //
2 // df_runtime.h
3 //
4 // Copyright (C) 2009 Edward Valeev
5 //
6 // Author: Edward Valeev <evaleev@vt.edu>
7 // Maintainer: EV
8 //
9 // This file is part of the SC Toolkit.
10 //
11 // The SC Toolkit is free software; you can redistribute it and/or modify
12 // it under the terms of the GNU Library General Public License as published by
13 // the Free Software Foundation; either version 2, or (at your option)
14 // any later version.
15 //
16 // The SC Toolkit is distributed in the hope that it will be useful,
17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 // GNU Library General Public License for more details.
20 //
21 // You should have received a copy of the GNU Library General Public License
22 // along with the SC Toolkit; see the file COPYING.LIB. If not, write to
23 // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
24 //
25 // The U.S. Government is granted a limited license as per AL 91-7.
26 //
27 
28 #ifndef _mpqc_src_lib_chemistry_qc_df_df_runtime_h
29 #define _mpqc_src_lib_chemistry_qc_df_df_runtime_h
30 
31 #include <chemistry/qc/lcao/df.h>
32 #include <chemistry/qc/lcao/tbint_runtime.h>
33 #include <Eigen/Dense>
34 #include <util/misc/sharedptr.h>
35 
36 
37 namespace sc {
38 
53  public:
54  ParsedDensityFittingKey(const std::string& key);
55 
56  const std::string& key() const { return key_; }
57  const std::string& space1() const { return space1_; }
58  const std::string& space2() const { return space2_; }
59  const std::string& fspace() const { return fspace_; }
60  const std::string& kernel() const { return kernel_pkey_.key(); }
61  const std::string& kernel_oper() const { return kernel_pkey_.oper(); }
62  const std::string& kernel_param() const { return kernel_pkey_.params(); }
63  bool ri() const { return ri_; }
64 
66  static std::string key(const std::string& space1,
67  const std::string& space2,
68  const std::string& fspace,
69  const std::string& kernel,
70  bool ri = false);
71 
72  private:
73  std::string key_;
74  std::string space1_, space2_, fspace_;
75  ParsedTwoBodyOperSetKey kernel_pkey_;
76  bool ri_;
77  };
78 
80 
84  class DensityFittingRuntime : virtual public SavableState {
85  public:
87  typedef DistArray4 Result;
88  typedef Ref<Result> ResultRef;
91 #ifdef MPQC_NEW_FEATURES
92  typedef Eigen::VectorXd CoefResult;
93  typedef std::shared_ptr<Eigen::VectorXd> CoefResultRef;
94  typedef std::pair<int, int> IntPair;
95  typedef std::pair<std::string, IntPair> CoefKey;
96 #endif // MPQC_NEW_FEATURES
97 
98  // uses MOIntsRuntime to evaluate integrals
102  void save_data_state(StateOut& so);
103 
105  void obsolete();
106 
108  const DensityFittingParams* dfparams() const { return dfparams_; }
109 
112  bool exists(const std::string& key) const;
113 
114 #ifdef MPQC_NEW_FEATURES
115 
117  bool exists(const CoefKey& key) const;
118 #endif // MPQC_NEW_FEATURES
119 
126  ResultRef get(const std::string& key); // non-const: can compute something
127 
128 #ifdef MPQC_NEW_FEATURES
129 
132  const CoefResultRef get(const CoefKey& key);
133  const CoefResultRef get(const std::string& dfkey, int bf1, int bf2){ return get(CoefKey(dfkey, IntPair(bf1, bf2))); }
134 #endif // MPQC_NEW_FEATURES
135 
137  const Ref<MOIntsRuntime>& moints_runtime() const { return moints_runtime_; }
138 
140  void remove_if(const std::string& space_key);
141 
150  static std::string default_dfbs_name(const std::string& obs_name,
151  int incX = 0,
152  bool force_aug = false);
153 
154  private:
155  Ref<MOIntsRuntime> moints_runtime_;
156  const DensityFittingParams* dfparams_;
157 
159  Ref<ResultRegistry> results_;
160 
161 #ifdef MPQC_NEW_FEATURES
163  Ref<CoefRegistry> coef_results_;
164 
165  CoefResultRef get_coefficients(const CoefKey& key);
166 
167  typedef Eigen::HouseholderQR<Eigen::MatrixXd> Decomposition;
168  typedef std::map<IntPair, std::shared_ptr<Decomposition> > DecompositionMap;
169  DecompositionMap decomps_;
170 #endif // MPQC_NEW_FEATURES
171 
172  // creates the result for a given key
173  const ResultRef& create_result(const std::string& key);
174 
175  static ClassDesc class_desc_;
176 
177  };
178 
182  class DensityFittingParams : virtual public SavableState {
183  public:
194  const std::string& kernel = std::string(),
195  const std::string& solver = std::string("cholesky_inv"));
198  void save_data_state(StateOut&);
199 
200  const Ref<GaussianBasisSet>& basis() const { return basis_; }
201  const std::string& kernel_key() const { return kernel_; }
202  DensityFitting::SolveMethod solver() const { return solver_; }
203  bool local_coulomb() const { return local_coulomb_; }
204  void local_coulomb(bool val) { local_coulomb_ = val; }
205  bool local_exchange() const { return local_exchange_; }
206  void local_exchange(bool val) { local_exchange_ = val; }
207  bool exact_diag_J() const { return exact_diag_J_; }
208  void exact_diag_J(bool val) { exact_diag_J_ = val; }
209  bool exact_diag_K() const { return exact_diag_K_; }
210  void exact_diag_K(bool val) { exact_diag_K_ = val; }
215  std::string intparams_key() const;
216 
217  void print(std::ostream& o) const;
218  private:
219  static ClassDesc class_desc_;
220 
221  bool local_coulomb_;
222  bool local_exchange_;
223  bool exact_diag_J_;
224  bool exact_diag_K_;
225  Ref<GaussianBasisSet> basis_;
226  std::string kernel_;
227  DensityFitting::SolveMethod solver_;
228  };
229 
230  inline bool operator==(const DensityFittingParams& A, const DensityFittingParams& B) {
231  return A.basis()->equiv(B.basis()) && A.kernel_key() == B.kernel_key() && A.solver() == B.solver();
232  }
233 
235  struct DensityFittingInfo : virtual public SavableState {
236  public:
238 
240  const Ref<DensityFittingRuntime>& r) :
241  params_(p), runtime_(r) {}
243  params_ << SavableState::restore_state(si);
244  runtime_ << SavableState::restore_state(si);
245  }
247  SavableState::save_state(params_.pointer(),so);
248  SavableState::save_state(runtime_.pointer(),so);
249  }
250 
251  const Ref<DensityFittingParams>& params() const { return params_; }
252  const Ref<DensityFittingRuntime>& runtime() const { return runtime_; }
253 
255  void obsolete();
256 
257  private:
258  Ref<DensityFittingParams> params_;
259  Ref<DensityFittingRuntime> runtime_;
260 
261  static ClassDesc class_desc_;
262  };
263 
264  inline bool operator==(const DensityFittingInfo& A, const DensityFittingInfo& B) {
265  // dfinfo objects are equivalent if they use equivalent params and same runtime object
266  return (*A.params() == *B.params()) && A.runtime() == B.runtime();
267  }
268 
269 
270 } // end of namespace sc
271 
272 #endif // end of header guard
273 
274 
275 // Local Variables:
276 // mode: c++
277 // c-file-style: "CLJ-CONDENSED"
278 // End:
sc::TwoBodyMOIntsRuntimeUnion23
TwoBodyMOIntsRuntimeUnion23 packages 2-center and 3-center runtimes; it also keeps track of 2-center ...
Definition: tbint_runtime.h:523
sc::DensityFittingParams
DensityFittingParams defines parameters used by DensityFittingRuntime and other runtime components to...
Definition: df_runtime.h:182
sc::DistArray4
DistArray4 contains a set of one or more distributed dense 4-index arrays.
Definition: distarray4.h:94
sc::ParsedTwoBodyOperSetKey
Parsed representation of a string key that represents a two-body operator set (TwoBodyOperSet + assoc...
Definition: tbint_runtime.h:77
sc::Ref
A template class that maintains references counts.
Definition: ref.h:361
sc::ParsedDensityFittingKey
Parsed representation of a string key that represents fitting of a product of space1 and space2 into ...
Definition: df_runtime.h:52
sc::DensityFittingRuntime
Smart runtime support for managing DensityFitting objects.
Definition: df_runtime.h:84
sc::DensityFittingParams::print
void print(std::ostream &o) const
Print the object.
sc::DensityFittingRuntime::default_dfbs_name
static std::string default_dfbs_name(const std::string &obs_name, int incX=0, bool force_aug=false)
tries to translate a library basis set label to the corresponding default value for the DF basis
sc::DensityFittingRuntime::obsolete
void obsolete()
obsoletes this object
sc::DensityFittingParams::intparams_key
std::string intparams_key() const
returns the IntParams object corresponding to the fitting kernel_key
sc::DensityFittingRuntime::save_data_state
void save_data_state(StateOut &so)
Save the base classes (with save_data_state) and the members in the same order that the StateIn CTOR ...
sc::DensityFittingRuntime::moints_runtime
const Ref< MOIntsRuntime > & moints_runtime() const
returns the runtime used to compute results
Definition: df_runtime.h:137
sc::DensityFittingRuntime::get
ResultRef get(const std::string &key)
Returns the DistArray4 object corresponding to this key.
sc::DensityFittingInfo
this class encapsulates objects needed to perform density fitting of a 4-center integral
Definition: df_runtime.h:235
sc::StateIn
Definition: statein.h:79
sc::ClassDesc
This class is used to contain information about classes.
Definition: class.h:147
sc::DensityFittingInfo::save_data_state
void save_data_state(StateOut &so)
Save the base classes (with save_data_state) and the members in the same order that the StateIn CTOR ...
Definition: df_runtime.h:246
sc::SavableState::restore_state
static SavableState * restore_state(StateIn &si)
Restores objects saved with save_state.
sc::DensityFittingRuntime::remove_if
void remove_if(const std::string &space_key)
removes all entries that contain this space
sc::DensityFittingRuntime::dfparams
const DensityFittingParams * dfparams() const
which density fitting solver will be used to compute?
Definition: df_runtime.h:108
sc::StateOut
Definition: stateout.h:71
sc::DensityFittingInfo::obsolete
void obsolete()
obsoletes all dependent runtimes
sc::SavableState::save_state
void save_state(StateOut &)
Save the state of the object as specified by the StateOut object.
sc::DensityFittingParams::DensityFittingParams
DensityFittingParams(const Ref< GaussianBasisSet > &basis, const std::string &kernel=std::string(), const std::string &solver=std::string("cholesky_inv"))
sc::DensityFittingParams::kernel_otype
TwoBodyOper::type kernel_otype() const
returns the TwoBodyInt::oper_type object that specifies the type of the operator kernel_key used for ...
sc::TwoBodyOper::type
type
types of known two-body operators
Definition: operator.h:318
sc::DensityFittingParams::save_data_state
void save_data_state(StateOut &)
Save the base classes (with save_data_state) and the members in the same order that the StateIn CTOR ...
sc::SavableState
Base class for objects that can save/restore state.
Definition: state.h:45
sc::operator==
bool operator==(const Atom &a, const Atom &b)
sc::DensityFittingRuntime::exists
bool exists(const std::string &key) const
Returns true if the given DensityFitting is available.
sc
Contains all MPQC code up to version 3.
Definition: mpqcin.h:14
sc::Registry
Registry wraps std::map and can be policy-configured to act as a Singleton or a regular object.
Definition: registry.h:112

Generated at Sun Jan 26 2020 23:23:58 for MPQC 3.0.0-alpha using the documentation package Doxygen 1.8.16.