-
Notifications
You must be signed in to change notification settings - Fork 1
/
ex2b.py
84 lines (62 loc) · 1.7 KB
/
ex2b.py
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
from numpy import *
def energy(state):
if(1<=state<=2):
return 1.0
elif (3<=state<=6):
return 2.0
elif (7<=state<=12):
return 3.0
def delta(i,j):
if(i == j):
return 1
else:
return 0
def mapping(state):
if((state == 1) | (state==2)):
return "1" #(n=0, m=0)
elif((state == 3) | (state == 4)):
return "2" #(n=0,m=-1)
elif ((state == 5) | (state == 6)):
return "3" #(n=0,m=1)
elif ((state == 7) | (state == 8)):
return "4" #(n=0,m=-2)
elif ((state==9) | (state == 10)):
return "5" #(n=1, m=0)
elif ((state == 11) | (state == 12)):
return "6" #(n=0,m=2)
def spin(state):
if((state%2 == 0)):
return -2 #Spin down if the state is an even number
else:
return 2 #Spin up if the state is an odd number
def expref(N,w):
#<ref|H|ref>
E0 = 0.0
for i in range(1,N+1):
E0+=energy(i)
for j in range(1,N+1):
key1 = mapping(i) + mapping(j) + mapping(i) + mapping(j)
key2 = mapping(i) + mapping(j) + mapping(j) + mapping(i)
E0+= 0.5*(delta(spin(i),spin(i))*delta(spin(j),spin(j))*w[key1] - delta(spin(i),spin(j))*delta(spin(j),spin(i))*w[key2])
return E0
def exp_ref_refia(i,a,N,w):
#<ref|H|ref_i^a>
val = 0.0
for j in range(1,N+1):
key1 = mapping(i) + mapping(j) + mapping(a) + mapping(j)
key2 = mapping(i) + mapping(j) + mapping(j) + mapping(a)
val+= delta(spin(i),spin(a))*delta(spin(j),spin(j))*w[key1] - delta(spin(i),spin(j))*delta(spin(j),spin(a))*w[key2]
return val
inFile = open('coulomb.dat','r')
w = {}
for line in inFile:
tmp = line.split()
key = tmp[0] + tmp[1] + tmp[2] + tmp[3]
val = float(tmp[4])
w[key] = val
inFile.close()
N = 12
exp_ref = expref(N,w)
exp_ref_19 = exp_ref_refia(1,9,N,w)
exp_ref_210 = exp_ref_refia(2,10,N,w)
print exp_ref