1# -*- coding: utf-8 -*- 

2# 

3# Copyright (c) 2018, the cclib development team 

4# 

5# This file is part of cclib (http://cclib.github.io) and is distributed under 

6# the terms of the BSD 3-Clause License. 

7 

8"""Calculation of overlap population analysis based on cclib data.""" 

9 

10import random 

11 

12import numpy 

13 

14from cclib.method.calculationmethod import Method 

15from cclib.method.population import Population 

16 

17 

18def func(x): 

19 if x==1: 

20 return 1 

21 else: 

22 return x+func(x-1) 

23 

24 

25class OPA(Population): 

26 """Overlap population analysis.""" 

27 

28 def __init__(self, *args): 

29 

30 # Call the __init__ method of the superclass. 

31 super(OPA, self).__init__(logname="OPA", *args) 

32 

33 def __str__(self): 

34 """Return a string representation of the object.""" 

35 return "OPA of %s" % (self.data) 

36 

37 def __repr__(self): 

38 """Return a representation of the object.""" 

39 return 'OPA("%s")' % (self.data) 

40 

41 def calculate(self, indices=None, fupdate=0.05): 

42 """Perform an overlap population analysis given the results of a parser""" 

43 if not indices: 

44 # Build list of groups of orbitals in each atom for atomresults. 

45 if hasattr(self.data, "aonames"): 

46 names = self.data.aonames 

47 elif hasattr(self.data, "foonames"): 

48 names = self.data.fonames 

49 

50 atoms = [] 

51 indices = [] 

52 

53 name = names[0].split('_')[0] 

54 atoms.append(name) 

55 indices.append([0]) 

56 

57 for i in range(1, len(names)): 

58 name = names[i].split('_')[0] 

59 try: 

60 index = atoms.index(name) 

61 except ValueError: #not found in atom list 

62 atoms.append(name) 

63 indices.append([i]) 

64 else: 

65 indices[index].append(i) 

66 

67 # Determine number of steps, and whether process involves beta orbitals. 

68 nfrag = len(indices) #nfrag 

69 nstep = func(nfrag - 1) 

70 unrestricted = (len(self.data.mocoeffs) == 2) 

71 alpha = len(self.data.mocoeffs[0]) 

72 nbasis = self.data.nbasis 

73 

74 self.logger.info("Creating attribute results: array[4]") 

75 results= [ numpy.zeros([nfrag, nfrag, alpha], "d") ] 

76 if unrestricted: 

77 beta = len(self.data.mocoeffs[1]) 

78 results.append(numpy.zeros([nfrag, nfrag, beta], "d")) 

79 nstep *= 2 

80 

81 if hasattr(self.data, "aooverlaps"): 

82 overlap = self.data.aooverlaps 

83 elif hasattr(self.data,"fooverlaps"): 

84 overlap = self.data.fooverlaps 

85 

86 #intialize progress if available 

87 if self.progress: 

88 self.progress.initialize(nstep) 

89 

90 size = len(self.data.mocoeffs[0]) 

91 step = 0 

92 

93 preresults = [] 

94 for spin in range(len(self.data.mocoeffs)): 

95 two = numpy.array([2.0]*len(self.data.mocoeffs[spin]),"d") 

96 

97 

98 # OP_{AB,i} = \sum_{a in A} \sum_{b in B} 2 c_{ai} c_{bi} S_{ab} 

99 

100 for A in range(len(indices)-1): 

101 

102 for B in range(A+1, len(indices)): 

103 

104 if self.progress: #usually only a handful of updates, so remove random part 

105 self.progress.update(step, "Overlap Population Analysis") 

106 

107 for a in indices[A]: 

108 

109 ca = self.data.mocoeffs[spin][:,a] 

110 

111 for b in indices[B]: 

112 

113 cb = self.data.mocoeffs[spin][:,b] 

114 temp = ca * cb * two *overlap[a,b] 

115 results[spin][A,B] = numpy.add(results[spin][A,B],temp) 

116 results[spin][B,A] = numpy.add(results[spin][B,A],temp) 

117 

118 step += 1 

119 

120 temparray2 = numpy.swapaxes(results[0],1,2) 

121 self.results = [ numpy.swapaxes(temparray2,0,1) ] 

122 if unrestricted: 

123 temparray2 = numpy.swapaxes(results[1],1,2) 

124 self.results.append(numpy.swapaxes(temparray2, 0, 1)) 

125 

126 if self.progress: 

127 self.progress.update(nstep, "Done") 

128 

129 return True