1 //
2 // obint.h
3 //
4 // Copyright (C) 1996 Limit Point Systems, Inc.
5 //
6 // Author: Curtis Janssen <cljanss@limitpt.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_obint_h
29 #define _chemistry_qc_basis_obint_h
30 
31 #ifdef __GNUC__
32 #pragma interface
33 #endif
34 
35 #include <util/ref/ref.h>
36 #include <util/state/state.h>
37 #include <math/scmat/matrix.h>
38 #include <math/scmat/elemop.h>
39 
40 #include <chemistry/qc/basis/gaussbas.h>
41 #include <chemistry/qc/basis/dercent.h>
42 
43 namespace sc {
44 
45 class Integral;
46 
47 // //////////////////////////////////////////////////////////////////////////
48 
49 class EfieldDotVectorData: public RefCount
50 {
51   public:
EfieldDotVectorData()52     EfieldDotVectorData() {};
53     ~EfieldDotVectorData();
54 
55     double position[3];
56     double vector[3];
57 
58     void set_position(double*);
59     void set_vector(double*);
60 };
61 
62 
63 class DipoleData: public RefCount
64 {
65   public:
66     double origin[3];
67 
DipoleData(double * d)68     DipoleData(double *d) {origin[0]=d[0]; origin[1]=d[1]; origin[2]=d[2];}
DipoleData()69     DipoleData() {origin[0]=origin[1]=origin[2]=0.0;}
70     ~DipoleData();
71     void set_origin(double*);
72 };
73 
74 
75 class PointChargeData: public RefCount
76 {
77   private:
78     int ncharges_;
79     const double *charges_;
80     const double *const*positions_;
81     double *alloced_charges_;
82     double **alloced_positions_;
83 
84   public:
85     // If copy_data is 0, the passed positions and charges will
86     // be stored (but not freed).
87     PointChargeData(int ncharge,
88                     const double *const*positions, const double *charges,
89                     int copy_data = 0);
90     ~PointChargeData();
91 
ncharges()92     int ncharges() const { return ncharges_; }
charges()93     const double *charges() const { return charges_; }
positions()94     const double *const*positions() const { return positions_; }
95 };
96 
97 
98 /** OneBodyInt is an abstract base class for objects that
99     compute integrals between two basis functions. */
100 class OneBodyInt : public RefCount {
101   protected:
102     // this is who created me
103     Integral *integral_;
104 
105     Ref<GaussianBasisSet> bs1_;
106     Ref<GaussianBasisSet> bs2_;
107 
108     double *buffer_;
109 
110     OneBodyInt(Integral *integral,
111                const Ref<GaussianBasisSet>&b1,
112                const Ref<GaussianBasisSet>&b2 = 0);
113 
114   public:
115     virtual ~OneBodyInt();
116 
117     /// Returns the number of basis functions on center one.
118     int nbasis() const;
119 
120     /// Returns the number of basis functions on the center one.
121     int nbasis1() const;
122     /// Returns the number of basis functions on the center two.
123     int nbasis2() const;
124 
125     /// Return the number of shells on center one.
126     int nshell() const;
127 
128     /// Return the number of shells on the center one.
129     int nshell1() const;
130     /// Return the number of shells on the center two.
131     int nshell2() const;
132 
133     /// Return the basis set on center one.
134     Ref<GaussianBasisSet> basis();
135 
136     /// Return the basis set on the center one.
137     Ref<GaussianBasisSet> basis1();
138     /// Return the basis set on the center two.
139     Ref<GaussianBasisSet> basis2();
140 
141     /// Returns the buffer where the integrals are placed.
142     const double * buffer() const;
143 
144     /** Computes the integrals between basis functions in the
145         given shell pair. */
146     virtual void compute_shell(int,int) = 0;
147 
148     /** This is called for one body integrals that take data to let
149         them know that the data they reference has changed. */
150     virtual void reinitialize();
151 
152     /** Return true if the clone member can be called.  The default
153      * implementation returns false. */
154     virtual bool cloneable();
155 
156     /** Returns a clone of this.  The default implementation throws an
157      * exception. */
158     virtual Ref<OneBodyInt> clone();
159 
integral()160     Integral *integral() const { return integral_; }
161 };
162 
163 // //////////////////////////////////////////////////////////////////////////
164 
165 /** OneBodyOneCenterInt is an abstract base class for objects that
166     compute integrals between two basis functions. */
167 class OneBodyOneCenterInt : public RefCount {
168   protected:
169     // this is who created me
170     Integral *integral_;
171 
172     Ref<GaussianBasisSet> bs1_;
173 
174     double *buffer_;
175 
176     OneBodyOneCenterInt(Integral *integral,
177                const Ref<GaussianBasisSet>&b1);
178 
179   public:
180     virtual ~OneBodyOneCenterInt();
181 
182     /// Returns the number of basis functions on center one.
183     int nbasis() const;
184 
185     /// Returns the number of basis functions on the center one.
186     int nbasis1() const;
187 
188     /// Return the number of shells on center one.
189     int nshell() const;
190 
191     /// Return the number of shells on the center one.
192     int nshell1() const;
193 
194     /// Return the basis set on center one.
195     Ref<GaussianBasisSet> basis();
196 
197     /// Return the basis set on the center one.
198     Ref<GaussianBasisSet> basis1();
199 
200     /// Returns the buffer where the integrals are placed.
201     const double * buffer() const;
202 
203     /** Computes the integrals for basis functions on the
204         given shell. */
205     virtual void compute_shell(int) = 0;
206 
207     /** This is called for one body integrals that take data to let
208         them know that the data they reference has changed. */
209     virtual void reinitialize();
210 
211     /** Return true if the clone member can be called.  The default
212      * implementation returns false. */
213     virtual bool cloneable();
214 
215     /** Returns a clone of this.  The default implementation throws an
216      * exception. */
217     virtual Ref<OneBodyOneCenterInt> clone();
218 
integral()219     Integral *integral() const { return integral_; }
220 };
221 
222 // //////////////////////////////////////////////////////////////////////////
223 
224 class OneBodyOneCenterWrapper : public OneBodyOneCenterInt {
225     Ref<OneBodyInt> ob_;
226     int jsh_;
227   public:
228     OneBodyOneCenterWrapper(const Ref<OneBodyInt>& ob,
229                             int sh2 = 0);
230     void compute_shell(int);
231 };
232 
233 // //////////////////////////////////////////////////////////////////////////
234 
235 class ShellPairIter {
236   private:
237     const double * buf;
238     double scale_;
239 
240     int e12;
241 
242     int index;
243 
244     int ioffset;
245     int joffset;
246 
247     int iend;
248     int jend;
249 
250     int icur;
251     int jcur;
252 
253   public:
254     ShellPairIter();
255     ~ShellPairIter();
256 
257     void init(const double * buffer, int ishell, int jshell,
258               int ioff, int joff, int nfunci, int nfuncj, int redund=0,
259               double scale=1.0);
260 
start()261     void start() { icur=jcur=index=0; }
ready()262     int ready() const { return (icur < iend); }
263 
next()264     void next() {
265       if (jcur < ((e12)?(icur):((jend)-1))) {
266         index++;
267         jcur++;
268         return;
269       }
270 
271       jcur=0;
272       icur++;
273 
274       index = icur*jend;
275     }
276 
current_i()277     int current_i() const { return icur; }
current_j()278     int current_j() const { return jcur; }
279 
i()280     int i() const { return icur+ioffset; }
j()281     int j() const { return jcur+joffset; }
282 
nint()283     int nint() const { return iend*jend; }
284 
val()285     double val() const { return buf[index]*scale_; }
286 };
287 
288 // //////////////////////////////////////////////////////////////////////////
289 
290 class OneBodyIntIter : public RefCount {
291   protected:
292     Ref<OneBodyInt> obi; // help me obi wan
293     ShellPairIter spi;
294 
295     int redund;
296 
297     int istart;
298     int jstart;
299 
300     int iend;
301     int jend;
302 
303     int icur;
304     int jcur;
305 
306     int ij;
307 
308   public:
309     OneBodyIntIter();
310     OneBodyIntIter(const Ref<OneBodyInt>&);
311     virtual ~OneBodyIntIter();
312 
313     virtual void start(int ist=0, int jst=0, int ien=0, int jen=0);
314     virtual void next();
315 
ready()316     int ready() const { return (icur < iend); }
317 
ishell()318     int ishell() const { return icur; }
jshell()319     int jshell() const { return jcur; }
320 
ijshell()321     int ijshell() const { return ij; }
322 
redundant()323     int redundant() const { return redund; }
set_redundant(int i)324     void set_redundant(int i) { redund=i; }
325 
326     virtual double scale() const;
327 
one_body_int()328     Ref<OneBodyInt> one_body_int() { return obi; }
329 
330     ShellPairIter& current_pair();
331 
332     virtual bool cloneable();
333     virtual Ref<OneBodyIntIter> clone();
334 };
335 
336 
337 
338 // //////////////////////////////////////////////////////////////////////////
339 
340 class OneBodyIntOp: public SCElementOp {
341   protected:
342     Ref<OneBodyIntIter> iter;
343 
344   public:
345     OneBodyIntOp(const Ref<OneBodyInt>&);
346     OneBodyIntOp(const Ref<OneBodyIntIter>&);
347     virtual ~OneBodyIntOp();
348 
349     void process(SCMatrixBlockIter&);
350     void process_spec_rect(SCMatrixRectBlock*);
351     void process_spec_ltri(SCMatrixLTriBlock*);
352     void process_spec_rectsub(SCMatrixRectSubBlock*);
353     void process_spec_ltrisub(SCMatrixLTriSubBlock*);
354 
355     bool cloneable();
356     Ref<SCElementOp> clone();
357 
358     int has_side_effects();
359 };
360 
361 class OneBody3IntOp: public SCElementOp3 {
362   private:
363     Ref<OneBodyIntIter> iter;
364 
365   public:
366     OneBody3IntOp(const Ref<OneBodyInt>&b);
367     OneBody3IntOp(const Ref<OneBodyIntIter>&);
368     virtual ~OneBody3IntOp();
369 
370     void process(SCMatrixBlockIter&,
371                  SCMatrixBlockIter&,
372                  SCMatrixBlockIter&);
373     void process_spec_rect(SCMatrixRectBlock*,
374                            SCMatrixRectBlock*,
375                            SCMatrixRectBlock*);
376     void process_spec_ltri(SCMatrixLTriBlock*,
377                            SCMatrixLTriBlock*,
378                            SCMatrixLTriBlock*);
379 
380     int has_side_effects();
381     int has_side_effects_in_arg1();
382     int has_side_effects_in_arg2();
383 
384 };
385 
386 // //////////////////////////////////////////////////////////////////////////
387 
388 /** OneBodyDerivInt is an abstract base class for objects that
389     compute one body derivative integrals. */
390 class OneBodyDerivInt : public RefCount {
391   protected:
392     // this is who created me
393     Integral *integral_;
394 
395     Ref<GaussianBasisSet> bs1;
396     Ref<GaussianBasisSet> bs2;
397 
398     double *buffer_;
399 
400   public:
401     OneBodyDerivInt(Integral *, const Ref<GaussianBasisSet>&b);
402     OneBodyDerivInt(Integral *,
403                     const Ref<GaussianBasisSet>&b1,
404                     const Ref<GaussianBasisSet>&b2);
405     virtual ~OneBodyDerivInt();
406 
407     /// Return the number of basis functions on center one.
408     int nbasis() const;
409     /// Return the number of basis functions on the center one.
410     int nbasis1() const;
411     /// Return the number of basis functions on the center two.
412     int nbasis2() const;
413 
414     /// Return the number of shells on center one.
415     int nshell() const;
416     /// Return the number of shells on center one.
417     int nshell1() const;
418     /// Return the number of shells on center two.
419     int nshell2() const;
420 
421     /// Return the basis set on center one.
422     Ref<GaussianBasisSet> basis();
423     /// Return the basis set on center one.
424     Ref<GaussianBasisSet> basis1();
425     /// Return the basis set on center two.
426     Ref<GaussianBasisSet> basis2();
427 
428     /** The computed shell integrals will be put in the buffer returned by
429         this member.  */
430     const double * buffer() const;
431 
432     /** Compute the derivative integrals and place the result in the buffer
433         returned by buffer(). */
434     virtual void compute_shell(int ish, int jsh, DerivCenters&) = 0;
435     /** Compute the derivative integrals with respect to the given center
436         and place the result in the buffer returned by buffer(). */
437     virtual void compute_shell(int ish, int jsh, int center) = 0;
438 };
439 
440 // //////////////////////////////////////////////////////////////////////////
441 
442 /** OneBodyOneCenterDerivInt is an abstract base class for objects that
443     compute one body derivative integrals on a single center. */
444 class OneBodyOneCenterDerivInt : public RefCount {
445   protected:
446     // this is who created me
447     Integral *integral_;
448 
449     Ref<GaussianBasisSet> bs1;
450 
451     double *buffer_;
452 
453   public:
454     OneBodyOneCenterDerivInt(Integral *, const Ref<GaussianBasisSet>&b);
455     virtual ~OneBodyOneCenterDerivInt();
456 
457     /// Return the number of basis functions on center one.
458     int nbasis() const;
459     /// Return the number of basis functions on center one.
460     int nbasis1() const;
461 
462     /// Return the number of shells on center one.
463     int nshell() const;
464     /// Return the number of shells on center one.
465     int nshell1() const;
466 
467     /// Return the basis set on center one.
468     Ref<GaussianBasisSet> basis();
469     /// Return the basis set on center one.
470     Ref<GaussianBasisSet> basis1();
471 
472     /** The computed shell integrals will be put in the buffer returned by
473         this member.  */
474     const double * buffer() const;
475 
476     /** Compute the derivative integrals and place the result in the buffer
477         returned by buffer(). */
478     virtual void compute_shell(int ish, DerivCenters&) = 0;
479     /** Compute the derivative integrals with respect to the given center
480         and place the result in the buffer returned by buffer(). */
481     virtual void compute_shell(int ish, int center) = 0;
482 };
483 
484 }
485 
486 #endif
487 
488 // Local Variables:
489 // mode: c++
490 // c-file-style: "ETS"
491 // End:
492