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
|
from math import pi as M_PI
theta0 = [0.0, 110.1, 111.0, 120.0, 108.5]
k = [0.0, 75.0, 45.0, 50.0, 100.0]
def compute(itype,theta):
dtheta = theta-theta0[itype]
tk = k[itype]*dtheta*M_PI/180.0*M_PI/180.0
phi = tk*dtheta
force = -2.0*tk
return (phi,force)
fp = open("angle_table.txt","w")
fp.write("# generated by gen_angle.py\n")
num = 72
fac = 180.0/num
for i in [1,2,3,4]:
fp.write("\nharmonic_%d\n" % i)
if i > 2:
fp.write("N %d EQ %g\n\n" % (num+1,theta0[i]))
else:
grad=2.0*k[i]*M_PI/180.0*M_PI/180.0
fp.write("N %d FP %g %g EQ %g\n\n" % (num+1,grad, grad, theta0[i]))
for j in range(num+1):
theta = fac*j
(phi,force) = compute(i,theta)
fp.write("%d %.15g %.15g %.15g\n" % (j+1,theta,phi,force))
fp.close()
|
