Psi4
libmints/molecule.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-2019 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_libmints_molecule_h_
30 #define _psi_src_lib_libmints_molecule_h_
31 
32 #include <vector>
33 #include <string>
34 #include <cstdio>
35 #include <map>
36 #include <memory>
37 
38 #define LINEAR_A_TOL 1.0E-2 // When sin(a) is below this, we consider the angle to be linear
39 #define DEFAULT_SYM_TOL 1.0E-8
40 #define FULL_PG_TOL 1.0e-8 // default
41 
42 #include "typedefs.h"
43 #include "coordentry.h"
44 
45 namespace psi {
46 class PointGroup;
47 class BasisSet;
66 };
67 
68 const std::string RotorTypeList[] = {"ASYMMETRIC_TOP", "SYMMETRIC_TOP", "SPHERICAL_TOP", "LINEAR", "ATOM"};
69 
70 const std::string FullPointGroupList[] = {"ATOM", "C_inf_v", "D_inf_h", "C1", "Cs", "Ci", "Cn", "Cnv",
71  "Cnh", "Sn", "Dn", "Dnd", "Dnh", "Td", "Oh", "Ih"};
72 
78  public:
85  };
89  enum GeometryUnits { Angstrom, Bohr };
93  enum FragmentType {
95  Real,
96  Ghost
97  };
98 
99  typedef std::vector<std::shared_ptr<CoordEntry>> EntryVector;
100  typedef EntryVector::iterator EntryVectorIter;
101  typedef std::map<std::string, std::string> Provenance;
102  typedef std::vector<std::tuple<int, int, double>> Connectivity;
103 
104  protected:
106  std::string name_;
108  std::string comment_;
118  std::vector<int> fragment_charges_;
120  std::vector<int> fragment_multiplicities_;
121 
126 
136  std::vector<std::string> all_variables_;
138  void clear();
139 
146  CoordValue* get_coord_value(const std::string& str);
147 
155  int get_anchor_atom(const std::string& str, const std::string& line);
156 
158  std::shared_ptr<PointGroup> pg_;
163 
165  int nunique_;
167  int* nequiv_;
169  int** equiv_;
172 
174  std::map<std::string, double> geometry_variables_;
176  std::vector<FragmentType> fragment_types_;
178  std::vector<std::pair<int, int>> fragments_;
180  std::string symmetry_from_input_;
188  bool zmat_;
190  bool cart_;
191 
192  public:
193  Molecule();
195  Molecule(const Molecule& other);
196  virtual ~Molecule();
197 
199  Molecule new_obj(*this);
200  return new_obj;
201  }
202 
206  Molecule& operator=(const Molecule& other);
208 
223  void add_atom(double Z, double x, double y, double z, std::string sym = "", double mass = 0.0, double charge = 0.0,
224  std::string lbl = "", int A = -1);
225  void add_unsettled_atom(double Z, std::vector<std::string> anchor, std::string sym = "", double mass = 0.0,
226  double charge = 0.0, std::string lbl = "", int A = -1);
227 
229  int nfragments() const { return fragments_.size(); }
231  int nactive_fragments();
233  void set_lock_frame(bool tf) { lock_frame_ = tf; }
234 
236  const std::string name() const { return name_; }
238  const std::string& basis_on_atom(int atom) const;
240  void set_name(const std::string& _name) { name_ = _name; }
242  const std::string comment() const { return comment_; }
244  void set_comment(const std::string &_comment) { comment_ = _comment; }
246  const Provenance provenance() const { return provenance_; }
248  void set_provenance(const Provenance &_provenance) { provenance_ = _provenance; }
250  const Connectivity connectivity() const { return connectivity_; }
252  void set_connectivity(const Connectivity &_connectivity) { connectivity_ = _connectivity; }
254  int natom() const;
256  int nallatom() const { return full_atoms_.size(); }
258  const double& Z(int atom) const;
260  double fZ(int atom) const;
262  double x(int atom) const;
264  double y(int atom) const;
266  double z(int atom) const;
268  double fx(int atom) const;
270  double fy(int atom) const;
272  double fz(int atom) const;
274  Vector3 xyz(int atom) const;
275  Vector3 fxyz(int atom) const;
277  double xyz(int atom, int _xyz) const;
279  double mass(int atom) const;
280 
282  void set_mass(int atom, double mass);
284  void set_nuclear_charge(int atom, double newZ);
285 
287  std::string symbol(int atom) const;
289  std::string fsymbol(int atom) const;
291  std::string label(int atom) const;
293  double charge(int atom) const;
295  int mass_number(int atom) const;
297  int true_atomic_number(int atom) const;
298  int ftrue_atomic_number(int atom) const;
300  double fmass(int atom) const;
302  std::string flabel(int atom) const;
304  double fcharge(int atom) const;
306  int fmass_number(int atom) const;
308  const std::shared_ptr<CoordEntry>& atom_entry(int atom) const;
309 
310  void set_basis_all_atoms(const std::string& name, const std::string& type = "BASIS");
311  void set_basis_by_symbol(const std::string& symbol, const std::string& name, const std::string& type = "BASIS");
312  void set_basis_by_number(int number, const std::string& name, const std::string& type = "BASIS");
313  void set_basis_by_label(const std::string& label, const std::string& name, const std::string& type = "BASIS");
314  void set_shell_by_label(const std::string& label, const std::string& name, const std::string& type = "BASIS");
315 
318  int atom_at_position1(double*, double tol = 0.05) const;
319  int atom_at_position2(Vector3&, double tol = 0.05) const;
320  int atom_at_position3(const std::array<double, 3>&, const double tol = 0.05) const;
322 
324  void set_reinterpret_coordentry(bool rc);
325 
327  Matrix geometry() const;
329  Matrix full_geometry() const;
330 
335  void set_geometry(double** geom);
336 
340  void set_geometry(const Matrix& geom);
341 
345  void set_full_geometry(double** geom);
346 
350  void set_full_geometry(const Matrix& geom);
351 
355  void rotate(const Matrix& R);
356  void rotate_full(const Matrix& R);
357 
361  void reinterpret_coordentries();
362 
367  void symmetrize_to_abelian_group(double tol);
368 
370  Vector3 center_of_mass() const;
373  double nuclear_repulsion_energy(const std::array<double, 3>& dipole_field) const;
375  Vector3 nuclear_dipole() const;
377  Vector3 nuclear_dipole(const Vector3& origin) const;
380  Matrix nuclear_repulsion_energy_deriv1(const std::array<double, 3>& dipole_field) const;
382  Matrix nuclear_repulsion_energy_deriv2() const;
383 
385  double pairwise_nuclear_repulsion_energy(std::shared_ptr<Molecule> other) const;
386 
388  void translate(const Vector3& r);
390  void move_to_com();
396  // void reorient();
397 
399  Matrix distance_matrix() const;
400 
402  Matrix* inertia_tensor() const;
403 
405  Vector rotational_constants(double tol = FULL_PG_TOL) const;
406 
408  void print_rotational_constants() const;
410  RotorType rotor_type(double tol = FULL_PG_TOL) const;
411 
413  void print() const;
414 
416  void print_cluster() const;
417 
419  void print_full() const;
420 
422  void print_in_input_format() const;
423 
425  void print_in_bohr() const;
426 
428  void print_in_angstrom() const;
429 
431  void print_distances() const;
432  void print_bond_angles() const;
433  void print_dihedrals() const;
434  void print_out_of_planes() const;
435 
437  void save_xyz_file(const std::string& filename, bool save_ghosts = true) const;
439  std::string save_string_xyz_file() const;
440 
442  std::string save_string_xyz() const;
443 
449  int nunique() const { return nunique_; }
451  int unique(int iuniq) const { return equiv_[iuniq][0]; }
453  int nequivalent(int iuniq) const { return nequiv_[iuniq]; }
455  int equivalent(int iuniq, int j) const { return equiv_[iuniq][j]; }
458  int atom_to_unique(int iatom) const { return atom_to_unique_[iatom]; }
461  int atom_to_unique_offset(int iatom) const;
463  int max_nequivalent() const;
465 
469  bool has_symmetry_element(Vector3& op, double tol = DEFAULT_SYM_TOL) const;
470  std::shared_ptr<PointGroup> point_group() const;
471  void set_point_group(std::shared_ptr<PointGroup> pg);
473  void set_full_point_group(double tol = FULL_PG_TOL);
475  bool has_inversion(Vector3& origin, double tol = DEFAULT_SYM_TOL) const;
477  bool is_plane(Vector3& origin, Vector3& uperp, double tol = DEFAULT_SYM_TOL) const;
479  bool is_axis(Vector3& origin, Vector3& axis, int order, double tol = DEFAULT_SYM_TOL) const;
481  void is_linear_planar(bool& linear, bool& planar, double tol = DEFAULT_SYM_TOL) const;
483  std::shared_ptr<PointGroup> find_point_group(double tol = DEFAULT_SYM_TOL) const;
485  void reset_point_group(const std::string& pgname);
487  std::shared_ptr<PointGroup> find_highest_point_group(double tol = DEFAULT_SYM_TOL) const;
490  std::shared_ptr<Matrix> symmetry_frame(double tol = DEFAULT_SYM_TOL);
492  void release_symmetry_information();
495  void form_symmetry_information(double tol = DEFAULT_SYM_TOL);
497  std::string sym_label();
499  std::vector<std::string> irrep_labels();
500  const std::string& symmetry_from_input() const { return symmetry_from_input_; }
501 
506  void symmetrize(double tol = 0.05, bool suppress_mol_print_in_exc = false);
508 
515  std::string create_psi4_string_from_molecule() const;
516 
520  void activate_all_fragments();
521 
525  void deactivate_all_fragments();
526 
531  void set_active_fragments(std::vector<int> reals);
532 
537  void set_active_fragment(int fragment);
538 
543  void set_ghost_fragments(std::vector<int> ghosts);
544 
549  void set_ghost_fragment(int fragment);
550 
558  std::shared_ptr<Molecule> extract_subsets(const std::vector<int>& real_list,
559  const std::vector<int>& ghost_list) const;
560 
567  std::shared_ptr<Molecule> py_extract_subsets_1(std::vector<int> reals, std::vector<int> ghost);
568 
575  std::shared_ptr<Molecule> py_extract_subsets_2(std::vector<int> reals, int ghost = -1);
576 
583  std::shared_ptr<Molecule> py_extract_subsets_3(int reals, std::vector<int> ghost);
584 
591  std::shared_ptr<Molecule> py_extract_subsets_4(int reals, int ghost = -1);
592 
598  std::shared_ptr<Molecule> py_extract_subsets_5(std::vector<int> reals);
599 
605  std::shared_ptr<Molecule> py_extract_subsets_6(int reals);
606 
607  // => Fragment Composition <= //
608 
611  const std::vector<std::pair<int, int>>& get_fragments() const { return fragments_; }
613  const std::vector<FragmentType>& get_fragment_types() const { return fragment_types_; }
615  const std::vector<int>& get_fragment_charges() const { return fragment_charges_; }
617  const std::vector<int>& get_fragment_multiplicities() const { return fragment_multiplicities_; }
619  void set_fragment_pattern(const std::vector<std::pair<int, int>>, const std::vector<FragmentType>,
620  const std::vector<int>, const std::vector<int>);
621 
623  void set_has_cartesian(bool tf) { cart_ = tf; }
625  void set_has_zmatrix(bool tf) { zmat_ = tf; }
627  bool has_zmatrix() const { return zmat_; }
630  void set_variable(const std::string& str, double val);
632  void set_geometry_variable(const std::string& str, double val) { geometry_variables_[str] = val; }
635  double get_variable(const std::string& str);
638  bool is_variable(const std::string& str) const;
639 
641  void set_molecular_charge(int charge) { molecular_charge_ = charge; }
643  int molecular_charge() const { return molecular_charge_; }
645  void set_multiplicity(int mult) { multiplicity_ = mult; }
647  int multiplicity() const { return multiplicity_; }
648 
650  void set_units(GeometryUnits units);
652  GeometryUnits units() const { return units_; }
655  void set_input_units_to_au(double conv);
657  double input_units_to_au() const { return input_units_to_au_; }
659  bool orientation_fixed() const { return fix_orientation_; }
661  void set_orientation_fixed(bool fix = true) { fix_orientation_ = fix; }
663  bool com_fixed() const { return !move_to_com_; }
665  void set_com_fixed(bool fix = true) { move_to_com_ = !fix; }
667  std::string schoenflies_symbol() const;
669  bool valid_atom_map(double tol = 0.05) const;
671  std::string full_point_group() const;
673  std::string full_point_group_with_n() const { return FullPointGroupList[full_pg_]; }
675  int full_pg_n() { return full_pg_n_; }
678  int rotational_symmetry_number() const;
679 
684  void update_geometry();
685 };
686 
687 } // namespace psi
688 
689 #endif
int multiplicity() const
Get the multiplicity (defined as 2Ms + 1)
Definition: libmints/molecule.h:647
Definition: libmints/molecule.h:59
int ** equiv_
Equivalent atom mapping array.
Definition: libmints/molecule.h:169
bool reinterpret_coordentries_
Definition: libmints/molecule.h:183
Definition: libmints/molecule.h:53
Definition: libmints/molecule.h:55
const std::string & symmetry_from_input() const
Definition: libmints/molecule.h:500
Definition: libmints/molecule.h:64
Definition: libmints/molecule.h:63
const std::vector< std::pair< int, int > > & get_fragments() const
Definition: libmints/molecule.h:611
Definition: libmints/molecule.h:48
std::string symmetry_from_input_
Symmetry string from geometry specification.
Definition: libmints/molecule.h:180
int molecular_charge() const
Gets the molecular charge.
Definition: libmints/molecule.h:643
void set_connectivity(const Connectivity &_connectivity)
Set molecule connectivity.
Definition: libmints/molecule.h:252
const std::string comment() const
Get molecule comment.
Definition: libmints/molecule.h:242
Definition: manybody.h:66
bool fix_orientation_
Reorient or not?
Definition: libmints/molecule.h:123
Provenance provenance_
Molecule origin.
Definition: libmints/molecule.h:110
std::map< std::string, std::string > Provenance
Definition: libmints/molecule.h:101
const std::vector< FragmentType > & get_fragment_types() const
A list describing how to handle each fragment.
Definition: libmints/molecule.h:613
Definition: libmints/molecule.h:83
#define FULL_PG_TOL
Definition: libmints/molecule.h:40
FragmentType
Definition: libmints/molecule.h:93
void set_has_cartesian(bool tf)
Sets whether this molecule contains at least one cartesian entry.
Definition: libmints/molecule.h:623
void set_name(const std::string &_name)
Set molecule name.
Definition: libmints/molecule.h:240
std::vector< std::shared_ptr< CoordEntry > > EntryVector
Definition: libmints/molecule.h:99
Definition: libmints/molecule.h:60
Definition: libmints/molecule.h:54
EntryVector::iterator EntryVectorIter
Definition: libmints/molecule.h:100
void set_comment(const std::string &_comment)
Set molecule comment.
Definition: libmints/molecule.h:244
double input_units_to_au() const
Gets the geometry unit conversion.
Definition: libmints/molecule.h:657
Definition: libmints/molecule.h:48
const std::vector< int > & get_fragment_charges() const
The charge of each fragment.
Definition: libmints/molecule.h:615
int molecular_charge_
The molecular charge.
Definition: libmints/molecule.h:128
int nunique_
Number of unique atoms.
Definition: libmints/molecule.h:165
Definition: vector3.h:40
RotorType
Definition: libmints/molecule.h:48
const std::string RotorTypeList[]
Definition: libmints/molecule.h:68
GeometryUnits units_
The units used to define the geometry.
Definition: libmints/molecule.h:132
std::vector< FragmentType > fragment_types_
A list describing how to handle each fragment.
Definition: libmints/molecule.h:176
Definition: libmints/molecule.h:51
std::string comment_
Molecule comment.
Definition: libmints/molecule.h:108
std::vector< std::pair< int, int > > fragments_
The list of atom ranges defining each fragment from parent molecule.
Definition: libmints/molecule.h:178
std::map< std::string, double > geometry_variables_
A listing of the variables used to define the geometries.
Definition: libmints/molecule.h:174
Definition: libmints/molecule.h:57
int nequivalent(int iuniq) const
Returns the number of atoms equivalent to iuniq.
Definition: libmints/molecule.h:453
int equivalent(int iuniq, int j) const
Returns the j&#39;th atom equivalent to iuniq.
Definition: libmints/molecule.h:455
std::shared_ptr< PointGroup > pg_
Point group to use with this molecule.
Definition: libmints/molecule.h:158
void set_lock_frame(bool tf)
Set whether to leave the geometry alone upon update_geometry()
Definition: libmints/molecule.h:233
Connectivity connectivity_
Molecule connectivity.
Definition: libmints/molecule.h:112
int full_pg_n()
Return n in Cnv, etc.; If there is no n (e.g. Td) it&#39;s the highest-order rotation axis...
Definition: libmints/molecule.h:675
std::vector< int > fragment_multiplicities_
The multiplicity of each fragment.
Definition: libmints/molecule.h:120
Definition: libmints/molecule.h:50
Definition: gshell.h:42
GeometryFormat
Definition: libmints/molecule.h:82
const std::string name() const
Get molecule name.
Definition: libmints/molecule.h:236
int multiplicity_
The multiplicity (defined as 2Ms + 1)
Definition: libmints/molecule.h:130
const std::string FullPointGroupList[]
Definition: libmints/molecule.h:70
std::string full_point_group_with_n() const
Return point group name such as Cnv or Sn.
Definition: libmints/molecule.h:673
bool lock_frame_
Definition: libmints/molecule.h:186
Definition: libmints/molecule.h:62
int nfragments() const
The number of fragments in the molecule.
Definition: libmints/molecule.h:229
bool has_zmatrix() const
Whether this molecule has at least one zmatrix entry.
Definition: libmints/molecule.h:627
bool com_fixed() const
Get whether or not COM is fixed.
Definition: libmints/molecule.h:663
std::vector< std::tuple< int, int, double > > Connectivity
Definition: libmints/molecule.h:102
int full_pg_n_
n of the highest rotational axis Cn
Definition: libmints/molecule.h:162
Definition: libmints/molecule.h:58
EntryVector full_atoms_
Atom info vector (includes dummy atoms)
Definition: libmints/molecule.h:116
Definition: libmints/molecule.h:48
int nunique() const
Definition: libmints/molecule.h:449
Makes using matrices just a little easier.
Definition: libmints/matrix.h:112
FullPointGroup full_pg_
Full point group.
Definition: libmints/molecule.h:160
int atom_to_unique(int iatom) const
Definition: libmints/molecule.h:458
bool orientation_fixed() const
Get whether or not orientation is fixed.
Definition: libmints/molecule.h:659
void set_geometry_variable(const std::string &str, double val)
Plain assigns the vlue val to the variable labeled string in the list of geometry variables...
Definition: libmints/molecule.h:632
Definition: coordentry.h:51
Definition: libmints/molecule.h:48
int * nequiv_
Number of equivalent atoms per unique atom (length nunique_)
Definition: libmints/molecule.h:167
#define DEFAULT_SYM_TOL
Definition: libmints/molecule.h:39
const Provenance provenance() const
Get molecule provenance.
Definition: libmints/molecule.h:246
Definition: libmints/molecule.h:48
void set_com_fixed(bool fix=true)
Fix the center of mass at its current frame.
Definition: libmints/molecule.h:665
Definition: libmints/molecule.h:94
FullPointGroup
Definition: libmints/molecule.h:49
void set_orientation_fixed(bool fix=true)
Fix the orientation at its current frame.
Definition: libmints/molecule.h:661
#define PSI_API
Definition: pragma.h:155
bool zmat_
Whether this molecule has at least one zmatrix entry.
Definition: libmints/molecule.h:188
Molecule clone()
Definition: libmints/molecule.h:198
const std::vector< int > & get_fragment_multiplicities() const
The multiplicity of each fragment.
Definition: libmints/molecule.h:617
bool cart_
Whether this molecule has at least one cartesian entry.
Definition: libmints/molecule.h:190
std::vector< int > fragment_charges_
The charge of each fragment.
Definition: libmints/molecule.h:118
std::string name_
Molecule (or fragment) name.
Definition: libmints/molecule.h:106
void set_molecular_charge(int charge)
Sets the molecular charge.
Definition: libmints/molecule.h:641
EntryVector atoms_
Atom info vector (no knowledge of dummy atoms)
Definition: libmints/molecule.h:114
const Connectivity connectivity() const
Get molecule connectivity.
Definition: libmints/molecule.h:250
Definition: libmints/molecule.h:52
int unique(int iuniq) const
Returns the overall number of the iuniq&#39;th unique atom.
Definition: libmints/molecule.h:451
void set_provenance(const Provenance &_provenance)
Set molecule provenance.
Definition: libmints/molecule.h:248
Molecule information class.
Definition: libmints/molecule.h:77
double input_units_to_au_
The conversion factor to take input units to Bohr.
Definition: libmints/molecule.h:134
Definition: libmints/molecule.h:65
Definition: libmints/molecule.h:95
void set_multiplicity(int mult)
Sets the multiplicity (defined as 2Ms + 1)
Definition: libmints/molecule.h:645
GeometryUnits units() const
Gets the geometry units.
Definition: libmints/molecule.h:652
Definition: libmints/molecule.h:61
int nallatom() const
Number of all atoms (includes dummies)
Definition: libmints/molecule.h:256
std::vector< std::string > all_variables_
A list of all variables known, whether they have been set or not.
Definition: libmints/molecule.h:136
Definition: vector.h:44
void set_has_zmatrix(bool tf)
Sets whether this molecule contains at least one zmatrix entry.
Definition: libmints/molecule.h:625
GeometryUnits
Definition: libmints/molecule.h:89
int * atom_to_unique_
Atom to unique atom mapping array (length natom)
Definition: libmints/molecule.h:171
Definition: libmints/molecule.h:56
bool move_to_com_
Move to center of mass or not?
Definition: libmints/molecule.h:125