Psi4
csg.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-2017 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 _psi_src_lib_libcubeprop_csg_h_
30 #define _psi_src_lib_libcubeprop_csg_h_
31 
32 #include <map>
33 #include <set>
34 
35 #include "psi4/libmints/typedefs.h"
36 
37 namespace psi {
38 
39 class Options;
40 class BasisExtents;
41 class RKSFunctions;
42 class BlockOPoints;
43 
45 
46 protected:
47 
48  // => Input Specification <= //
49 
53  std::shared_ptr<Molecule> mol_;
55  std::shared_ptr<BasisSet> primary_;
57  std::shared_ptr<BasisSet> auxiliary_;
59  std::string filepath_;
60 
61  // => Physical Grid <= //
62 
64  int* N_;
66  double* D_;
68  double* O_;
69 
71  size_t npoints_;
73  size_t nxyz_;
74 
76  double* x_;
78  double* y_;
80  double* z_;
82  double* w_;
83 
84  // => Grid Computers <= //
85 
87  std::vector<std::shared_ptr<BlockOPoints> > blocks_;
89  std::shared_ptr<BasisExtents> extents_;
91  std::shared_ptr<RKSFunctions> points_;
92 
93  // => Helper Routines <= //
94 
96  void populate_grid();
97 
98 public:
99  // => Constructors <= //
100 
101  CubicScalarGrid(std::shared_ptr<BasisSet> primary, Options& options);
102  virtual ~CubicScalarGrid();
103 
104  // => High-Level Setup Routines <= //
105 
107  void build_grid();
109  void build_grid(const std::string filepath, int* N, double* D, double* O);
111  void build_grid(std::shared_ptr<CubicScalarGrid> other);
113  void print_header();
114 
116  void set_auxiliary_basis(std::shared_ptr<BasisSet> aux) { auxiliary_ = aux; }
117 
118  // => High-Level Set Routines <= //
119 
121  void set_filepath(const std::string filepath) { filepath_ = filepath; }
122 
123  // => High-Level Accessor Routines <= //
124 
126  int* N() const { return N_; }
128  double* D() const { return D_; }
130  double* O() const { return O_; }
132  std::string filepath() const { return filepath_; }
133 
134  // => Low-Level Accessor Routines (Use only if you know what you are doing) <= //
135 
137  size_t npoints() const { return npoints_; }
139  size_t nxyz() const { return nxyz_; }
140 
142  double* x() const { return x_; }
144  double* y() const { return y_; }
146  double* z() const { return z_; }
148  double* w() const { return w_; }
149 
150  // => Low-Level Write Routines (Use only if you know what you are doing) <= //
151 
153  void write_gen_file(double* v, const std::string& name, const std::string& type);
155  void write_cube_file(double* v, const std::string& name);
156 
157  // => Low-Level Scalar Field Computation (Use only if you know what you are doing) <= //
158 
160  void add_density(double* v, std::shared_ptr<Matrix> D);
162  void add_esp(double* v, std::shared_ptr<Matrix> D, const std::vector<double>& nuc_weights = std::vector<double>());
164  void add_basis_functions(double** v, const std::vector<int>& indices);
166  void add_orbitals(double** v, std::shared_ptr<Matrix> C);
168  void add_LOL(double* v, std::shared_ptr<Matrix> D);
170  void add_ELF(double* v, std::shared_ptr<Matrix> D);
171 
172  // => High-Level Scalar Field Computation <= //
173 
175  void compute_density(std::shared_ptr<Matrix> D, const std::string& name, const std::string& type = "CUBE");
177  void compute_esp(std::shared_ptr<Matrix> D, const std::vector<double>& nuc_weights, const std::string& name, const std::string& type = "CUBE");
179  void compute_basis_functions(const std::vector<int>& indices, const std::string& name, const std::string& type = "CUBE");
181  void compute_orbitals(std::shared_ptr<Matrix> C, const std::vector<int>& indices, const std::vector<std::string>& labels, const std::string& name, const std::string& type = "CUBE");
183  void compute_LOL(std::shared_ptr<Matrix> D, const std::string& name, const std::string& type = "CUBE");
185  void compute_ELF(std::shared_ptr<Matrix> D, const std::string& name, const std::string& type = "CUBE");
186 
187 };
188 
189 } // End namespace
190 
191 #endif
double * y_
y coordinates of grid
Definition: csg.h:78
size_t npoints_
number of points of grid
Definition: csg.h:71
std::string filepath_
File path for grid storage.
Definition: csg.h:59
double * w() const
w weights (rectangular) in fast ordering
Definition: csg.h:148
void set_filepath(const std::string filepath)
Set the directory for cube file storage (defaults to &quot;./&quot;)
Definition: csg.h:121
int * N() const
Number of voxels in [x,y,z]. Number of points along each dimensions in N_i + 1.
Definition: csg.h:126
int * N_
Voxel quanta in x, y, z [(N_x+1) x (N_y + 1) x (N_z + 1) points].
Definition: csg.h:64
void build_grid()
Build grid with options overages.
Definition: csg.cc:103
Definition: csg.h:44
void write_gen_file(double *v, const std::string &name, const std::string &type)
Write a general file of the scalar field v (in fast ordering) to filepath/name.ext.
Definition: csg.cc:237
double * z() const
z points in fast ordering
Definition: csg.h:146
std::shared_ptr< BasisSet > auxiliary_
The auxiliary basisset for ESP contractions.
Definition: csg.h:57
std::string filepath() const
Filepath where grid output will be stored.
Definition: csg.h:132
void compute_orbitals(std::shared_ptr< Matrix > C, const std::vector< int > &indices, const std::vector< std::string > &labels, const std::string &name, const std::string &type="CUBE")
Compute a set of orbital-type properties and drop files corresponding to name, index, symmetry label, and type.
Definition: csg.cc:646
void compute_density(std::shared_ptr< Matrix > D, const std::string &name, const std::string &type="CUBE")
Compute a density-type property and drop a file corresponding to name and type.
Definition: csg.cc:618
virtual ~CubicScalarGrid()
Definition: csg.cc:69
size_t npoints() const
Total number of points in grid.
Definition: csg.h:137
void compute_ELF(std::shared_ptr< Matrix > D, const std::string &name, const std::string &type="CUBE")
Compute an ELF-type property and drop a file corresponding to name and type (TODO: this seems very un...
Definition: csg.cc:672
void add_LOL(double *v, std::shared_ptr< Matrix > D)
Add a LOL-type property to the scalar field.
Definition: csg.cc:558
void add_esp(double *v, std::shared_ptr< Matrix > D, const std::vector< double > &nuc_weights=std::vector< double >())
Add an ESP-type property to the scalar field (total density matrix, must set DF_BASIS_SCF option) ...
Definition: csg.cc:323
double * x() const
x points in fast ordering
Definition: csg.h:142
void add_basis_functions(double **v, const std::vector< int > &indices)
Add a basis function property for desired indices to the scalar fields in v (rows are basis functions...
Definition: csg.cc:513
size_t nxyz_
Sparsity blocking in all cardinal directions.
Definition: csg.h:73
std::shared_ptr< BasisSet > primary_
Basis set this grid is built around.
Definition: csg.h:55
std::vector< std::shared_ptr< BlockOPoints > > blocks_
Vector of blocks.
Definition: csg.h:87
void compute_esp(std::shared_ptr< Matrix > D, const std::vector< double > &nuc_weights, const std::string &name, const std::string &type="CUBE")
Compute an ESP-type property and drop a file corresponding to name and type.
Definition: csg.cc:626
void add_ELF(double *v, std::shared_ptr< Matrix > D)
Add an ELF-type property to the scalar field.
Definition: csg.cc:586
double * D() const
Voxel width in [x,y,z], in bohr.
Definition: csg.h:128
double * D_
Voxel spacing in x, y, z.
Definition: csg.h:66
void compute_LOL(std::shared_ptr< Matrix > D, const std::string &name, const std::string &type="CUBE")
Compute a LOL-type property and drop a file corresponding to name and type.
Definition: csg.cc:664
void write_cube_file(double *v, const std::string &name)
Write a Gaussian cube file of the scalar field v (in fast ordering) to filepath/name.cube.
Definition: csg.cc:245
size_t nxyz() const
Number of points.
Definition: csg.h:139
void add_density(double *v, std::shared_ptr< Matrix > D)
Add a density-type property to the scalar field.
Definition: csg.cc:309
void populate_grid()
Setup grid from info in N_, D_, O_.
Definition: csg.cc:152
void set_auxiliary_basis(std::shared_ptr< BasisSet > aux)
Set the auxiliary for ESP if desired.
Definition: csg.h:116
std::shared_ptr< BasisExtents > extents_
Points to basis extents, built internally.
Definition: csg.h:89
Definition: liboptions.h:360
const char * labels[]
Definition: petitelist.cc:1076
std::shared_ptr< Molecule > mol_
Molecule this grid is built around.
Definition: csg.h:53
void compute_basis_functions(const std::vector< int > &indices, const std::string &name, const std::string &type="CUBE")
Compute a set of basis function-type properties and drop files corresponding to name, index, and type.
Definition: csg.cc:634
double * O_
Voxel origin in x, y, z.
Definition: csg.h:68
Options & options_
Options object for overages and voxel spacing.
Definition: csg.h:51
void print_header()
Header info.
Definition: csg.cc:213
double * x_
x coordinates of grid
Definition: csg.h:76
double * y() const
y points in fast ordering
Definition: csg.h:144
indices
Definition: libdpd/dpd.h:246
void add_orbitals(double **v, std::shared_ptr< Matrix > C)
Add orbital property for desired indices to the scalar fields in v (rows are orbitals) ...
Definition: csg.cc:538
std::shared_ptr< RKSFunctions > points_
RKS points object.
Definition: csg.h:91
double * w_
w quadrature weights of grid (rectangular)
Definition: csg.h:82
CubicScalarGrid(std::shared_ptr< BasisSet > primary, Options &options)
Definition: csg.cc:50
double * z_
z coordinates of grid
Definition: csg.h:80
double * O() const
Lower-left corner of th grid, in bohr.
Definition: csg.h:130