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()