/usr/share/pyshared/cclib/method/cda.py is in python-cclib 1.1-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 | # This file is part of cclib (http://cclib.sf.net), a library for parsing
# and interpreting the results of computational chemistry packages.
#
# Copyright (C) 2007, the cclib development team
#
# The library is free software, distributed under the terms of
# the GNU Lesser General Public version 2.1 or later. You should have
# received a copy of the license along with cclib. You can also access
# the full license online at http://www.gnu.org/copyleft/lgpl.html.
__revision__ = "$Revision: 961 $"
import random # For sometimes running the progress updater
import numpy
from fragments import FragmentAnalysis
class CDA(FragmentAnalysis):
"""Charge Decomposition Analysis (CDA)"""
def __init__(self, *args):
# Call the __init__ method of the superclass.
super(FragmentAnalysis, self).__init__(logname="CDA", *args)
def __str__(self):
"""Return a string representation of the object."""
return "CDA of" % (self.data)
def __repr__(self):
"""Return a representation of the object."""
return 'CDA("%s")' % (self.data)
def calculate(self, fragments, cupdate=0.05):
"""Perform the charge decomposition analysis.
Inputs:
fragments - list of ccData data objects
"""
retval = super(CDA, self).calculate(fragments, cupdate)
if not retval:
return False
# At this point, there should be a mocoeffs and fooverlaps
# in analogy to a ccData object.
donations = []
bdonations = []
repulsions = []
residuals = []
if len(self.mocoeffs) == 2:
occs = 1
else:
occs = 2
# Intialize progress if available.
nstep = self.data.homos[0]
if len(self.data.homos) == 2:
nstep += self.data.homos[1]
if self.progress:
self.progress.initialize(nstep)
# Begin the actual method.
step = 0
for spin in range(len(self.mocoeffs)):
size = len(self.mocoeffs[spin])
homo = self.data.homos[spin]
if len(fragments[0].homos) == 2:
homoa = fragments[0].homos[spin]
else:
homoa = fragments[0].homos[0]
if len(fragments[1].homos) == 2:
homob = fragments[1].homos[spin]
else:
homob = fragments[1].homos[0]
print "handling spin unrestricted"
if spin == 0:
fooverlaps = self.fooverlaps
elif spin == 1 and hasattr(self, "fooverlaps2"):
fooverlaps = self.fooverlaps2
offset = fragments[0].nbasis
self.logger.info("Creating donations, bdonations, and repulsions: array[]")
donations.append(numpy.zeros(size, "d"))
bdonations.append(numpy.zeros(size, "d"))
repulsions.append(numpy.zeros(size, "d"))
residuals.append(numpy.zeros(size, "d"))
for i in range(self.data.homos[spin] + 1):
# Calculate donation for each MO.
for k in range(0, homoa + 1):
for n in range(offset + homob + 1, self.data.nbasis):
donations[spin][i] += 2 * occs * self.mocoeffs[spin][i,k] \
* self.mocoeffs[spin][i,n] * fooverlaps[k][n]
for l in range(offset, offset + homob + 1):
for m in range(homoa + 1, offset):
bdonations[spin][i] += 2 * occs * self.mocoeffs[spin][i,l] \
* self.mocoeffs[spin][i,m] * fooverlaps[l][m]
for k in range(0, homoa + 1):
for m in range(offset, offset+homob + 1):
repulsions[spin][i] += 2 * occs * self.mocoeffs[spin][i,k] \
* self.mocoeffs[spin][i, m] * fooverlaps[k][m]
for m in range(homoa + 1, offset):
for n in range(offset + homob + 1, self.data.nbasis):
residuals[spin][i] += 2 * occs * self.mocoeffs[spin][i,m] \
* self.mocoeffs[spin][i, n] * fooverlaps[m][n]
step += 1
if self.progress and random.random() < cupdate:
self.progress.update(step, "Charge Decomposition Analysis...")
if self.progress:
self.progress.update(nstep, "Done.")
self.donations = donations
self.bdonations = bdonations
self.repulsions = repulsions
self.residuals = residuals
return True
|