1 //
2 // natorbval.cc
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 #include <math/scmat/matrix.h>
29 #include <util/keyval/keyval.h>
30 #include <chemistry/qc/basis/gaussshell.h>
31 #include <chemistry/qc/basis/gaussbas.h>
32 #include <chemistry/qc/wfn/wfn.h>
33
34 using namespace sc;
35
36 // Function for returning an orbital value at a point
natural_orbital(const SCVector3 & r,int iorb)37 double Wavefunction::natural_orbital(const SCVector3& r, int iorb)
38 {
39 return orbital(r,iorb,natural_orbitals());
40 }
41
42 // Function for returning an orbital value at a point
natural_orbital_density(const SCVector3 & r,int iorb,double * orbval)43 double Wavefunction::natural_orbital_density(const SCVector3& r,
44 int iorb,
45 double* orbval)
46 {
47 return orbital_density(r,iorb,natural_orbitals(),orbval);
48 }
49
50 // Function for returning an orbital value at a point
orbital(const SCVector3 & r,int iorb,const RefSCMatrix & orbs)51 double Wavefunction::orbital(const SCVector3& r,
52 int iorb,
53 const RefSCMatrix& orbs)
54 {
55 int nbasis = basis()->nbasis();
56 if (!bs_values) bs_values=new double[nbasis];
57
58 // compute the basis set values
59 GaussianBasisSet::ValueData *valdat
60 = new GaussianBasisSet::ValueData(basis(), integral_);
61 basis()->values(r,valdat,bs_values);
62 delete valdat;
63
64 // loop over basis functions
65 double orb_value = 0;
66 for (int i=0; i<nbasis; i++)
67 orb_value += orbs.get_element(i,iorb)*bs_values[i];
68
69 return orb_value;
70 }
71
orbital_density(const SCVector3 & r,int iorb,const RefSCMatrix & orbs,double * orbvalue)72 double Wavefunction::orbital_density(const SCVector3& r,
73 int iorb,
74 const RefSCMatrix& orbs,
75 double* orbvalue)
76 {
77 double tmp = orbital(r,iorb,orbs);
78 if (orbvalue) *orbvalue = tmp;
79 return 2.0 * tmp * tmp;
80 }
81
82 /////////////////////////////////////////////////////////////////////////////
83
84 // Local Variables:
85 // mode: c++
86 // c-file-style: "CLJ"
87 // End:
88