MPQC 2.3.1
petite.h
1//
2// petite.h --- definition of the PetiteList class
3//
4// Copyright (C) 1996 Limit Point Systems, Inc.
5//
6// Author: Edward Seidl <seidl@janed.com>
7// Maintainer: LPS
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 _chemistry_qc_basis_petite_h
29#define _chemistry_qc_basis_petite_h
30
31#ifdef __GNUC__
32#pragma interface
33#endif
34
35#include <scconfig.h>
36#include <iostream>
37#include <scconfig.h>
38
39#include <util/misc/scint.h>
40#include <util/ref/ref.h>
41#include <math/scmat/blocked.h>
42#include <math/scmat/offset.h>
43#include <chemistry/molecule/molecule.h>
44#include <chemistry/qc/basis/gaussbas.h>
45#include <chemistry/qc/basis/integral.h>
46
47// //////////////////////////////////////////////////////////////////////////
48
49namespace sc {
50
51inline sc_int_least64_t
52ij_offset64(sc_int_least64_t i, sc_int_least64_t j)
53{
54 return (i>j) ? (((i*(i+1)) >> 1) + j) : (((j*(j+1)) >> 1) + i);
55}
56
57inline sc_int_least64_t
58i_offset64(sc_int_least64_t i)
59{
60 return ((i*(i+1)) >> 1);
61}
62
63// //////////////////////////////////////////////////////////////////////////
64// These are helper functions for PetiteList and GPetite4
65
66int **compute_atom_map(const Ref<GaussianBasisSet> &);
67void delete_atom_map(int **atom_map, const Ref<GaussianBasisSet> &);
68
69int **compute_shell_map(int **atom_map, const Ref<GaussianBasisSet> &);
70void delete_shell_map(int **shell_map, const Ref<GaussianBasisSet> &);
71
72// //////////////////////////////////////////////////////////////////////////
73
75 int bfn;
76 double coef;
77
79 contribution(int b, double c);
81};
82
83struct SO {
84 int len;
85 int length;
86 contribution *cont;
87
88 SO();
89 SO(int);
90 ~SO();
91
92 SO& operator=(const SO&);
93
94 void set_length(int);
95 void reset_length(int);
96
97 // is this equal to so to within a sign
98 int equiv(const SO& so);
99};
100
101struct SO_block {
102 int len;
103 SO *so;
104
105 SO_block();
106 SO_block(int);
107 ~SO_block();
108
109 void set_length(int);
110 void reset_length(int);
111
112 int add(SO& s, int i);
113 void print(const char *title);
114};
115
116// //////////////////////////////////////////////////////////////////////////
117// this should only be used from within a SymmGaussianBasisSet
118
119class PetiteList : public RefCount {
120 private:
121 int natom_;
122 int nshell_;
123 int ng_;
124 int nirrep_;
125 int nblocks_;
126 int c1_;
127
129 Ref<Integral> ints_;
130
131 char *p1_; // p1[n] is 1 if shell n is in the group P1
132 int **atom_map_; // atom_map[n][g] is the atom that symop g maps atom n
133 // into
134 int **shell_map_; // shell_map[n][g] is the shell that symop g maps shell n
135 // into
136 char *lamij_; // see Dupuis & King, IJQC 11,613,(1977)
137
138 int *nbf_in_ir_;
139
140 void init();
141
142 public:
144 ~PetiteList();
145
146 Ref<GaussianBasisSet> basis() { return gbs_; }
147 Ref<Integral> integral() { return ints_; }
148 Ref<PetiteList> clone() { return new PetiteList(gbs_, ints_); }
149
150 int nirrep() const { return nirrep_; }
151 int order() const { return ng_; }
152 int atom_map(int n, int g) const { return (c1_) ? n : atom_map_[n][g]; }
153 int shell_map(int n, int g) const { return (c1_) ? n : shell_map_[n][g]; }
154 int lambda(int ij) const { return (c1_) ? 1 : (int) lamij_[ij]; }
155 int lambda(int i, int j) const
156 { return (c1_) ? 1 : (int) lamij_[ij_offset(i,j)]; }
157
158 int in_p1(int n) const { return (c1_) ? 1 : (int) p1_[n]; }
159 int in_p2(int ij) const { return (c1_) ? 1 : (int) lamij_[ij]; }
161 int in_p2(int i, int j) const { return (c1_) ? 1 : (int) lamij_[ij_offset(i,j)]; }
162 int in_p4(int ij, int kl, int i, int j, int k, int l) const;
164 int in_p4(int i, int j, int k, int l) const;
165
166 int nfunction(int i) const
167 { return (c1_) ? gbs_->nbasis() : nbf_in_ir_[i]; }
168
169 int nblocks() const { return nblocks_; }
170
171 void print(std::ostream& =ExEnv::out0(), int verbose=1);
172
173 // these return blocked dimensions
174 RefSCDimension AO_basisdim();
175 RefSCDimension SO_basisdim();
176
177 // return the basis function rotation matrix R(g)
178 RefSCMatrix r(int g);
179
180 // return information about the transformation from AOs to SOs
181 SO_block * aotoso_info();
182 RefSCMatrix aotoso();
183 RefSCMatrix sotoao();
184
185 // given a skeleton matrix, form the symmetrized matrix in the SO basis
186 void symmetrize(const RefSymmSCMatrix& skel, const RefSymmSCMatrix& sym);
187
188 // transform a matrix from AO->SO or SO->AO.
189 // this can take either a blocked or non-blocked AO basis matrix.
190 RefSymmSCMatrix to_SO_basis(const RefSymmSCMatrix&);
191
192 // this returns a non-blocked AO basis matrix.
193 RefSymmSCMatrix to_AO_basis(const RefSymmSCMatrix&);
194
195 // these two are just for eigenvectors
196 // returns non-blocked AO basis eigenvectors
197 RefSCMatrix evecs_to_AO_basis(const RefSCMatrix&);
198 // returns blocked SO basis eigenvectors
199 RefSCMatrix evecs_to_SO_basis(const RefSCMatrix&);
200};
201
202inline int
203PetiteList::in_p4(int ij, int kl, int i, int j, int k, int l) const
204{
205 if (c1_)
206 return 1;
207
208 sc_int_least64_t ijkl = i_offset64(ij)+kl;
209 int nijkl=1;
210
211 for (int g=1; g < ng_; g++) {
212 int gij = ij_offset(shell_map_[i][g],shell_map_[j][g]);
213 int gkl = ij_offset(shell_map_[k][g],shell_map_[l][g]);
214 sc_int_least64_t gijkl = ij_offset64(gij,gkl);
215
216 if (gijkl > ijkl)
217 return 0;
218 else if (gijkl == ijkl)
219 nijkl++;
220 }
221
222 return ng_/nijkl;
223}
224
225inline int
226PetiteList::in_p4(int i, int j, int k, int l) const
227{
228 if (c1_)
229 return 1;
230
231 int ij = ij_offset(i,j);
232 int kl = ij_offset(k,l);
233 sc_int_least64_t ijkl = ij_offset64(ij,kl);
234 int nijkl=1;
235
236 for (int g=1; g < ng_; g++) {
237 int gij = ij_offset(shell_map_[i][g],shell_map_[j][g]);
238 int gkl = ij_offset(shell_map_[k][g],shell_map_[l][g]);
239 sc_int_least64_t gijkl = ij_offset64(gij,gkl);
240
241 if (gijkl > ijkl)
242 return 0;
243 else if (gijkl == ijkl)
244 nijkl++;
245 }
246
247 return ng_/nijkl;
248}
249
250}
251
252
253
254#endif
255
256// Local Variables:
257// mode: c++
258// c-file-style: "ETS"
259// End:
static std::ostream & out0()
Return an ostream that writes from node 0.
Definition petite.h:119
int in_p2(int i, int j) const
Same as previous, except for it takes i and j separately.
Definition petite.h:161
The base class for all reference counted objects.
Definition ref.h:194
A template class that maintains references counts.
Definition ref.h:332
Definition petite.h:101
Definition petite.h:83
Definition petite.h:74

Generated at Thu Jul 18 2024 00:00:00 for MPQC 2.3.1 using the documentation package Doxygen 1.11.0.