blob: 272c7154791b70148bbbb3f1606378deeaf10d66 (
plain)
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
|
#! /usr/bin/python
from math import *
r0=1
De=1
ks=100
beta = sqrt(ks/(2*De))
def morse(r):
return De*(1-exp(-sqrt(ks/(2*De))*(r-r0)))**2
def lipp(r):
return De*(1-exp(-ks*r0*(r-r0)**2/(2*De*r)))
def lippmor(r):
if r<r0: return morse(r) - De
else: return lipp(r) - De
def flp(x):
k1 = 4.55328
k2 = -8.24688
k3 = 4.7872
c = 4.00777
return c*(1/(k1*x**2+k2*x+k3)**16-1/(k1*x**2+k2*x+k3)**8)
def fmo(x):
k1 = 2303.27
k2 = 1.12845
k3 = 1.45109
c = 59.1548
k1 = 95.1566
k2 = 1.3563
k3 = 0.416429
c = 46234.9
return k1*(c/(k2+x**2)**16-1/(k3+x**2)**8)
def flipmor(r):
if r<1.03946: return fmo(r)
else: return flp(r)
for i in range(600):
x=i*0.001 + .9
print x, flp(x), fmo(x)
# print x, lippmor(x)
# print x, flipmor(x),(flipmor(x+0.001)-flipmor(x-0.001))*1e3
## r<0, quartic
## a = 2689.48 +/- 12.49 (0.4645%)
## b = -11027.6 +/- 48.73 (0.4419%)
## c = 16997.7 +/- 71.21 (0.419%)
## d = -11670.6 +/- 46.22 (0.396%)
## e = 3009.98 +/- 11.24 (0.3734%)
## r>0
## a = 274.58 +/- 6.785 (2.471%)
## b = -1349 +/- 29.58 (2.193%)
## c = 2459.39 +/- 48.31 (1.964%)
## d = -1970.22 +/- 35.02 (1.777%)
## e = 584.236 +/- 9.506 (1.627%)
## for r>r0, using +/ a[i]/x^i
## a = 214.69 +/- 1.915 (0.8918%)
## b = -1092.98 +/- 9.376 (0.8578%)
## c = 2068.38 +/- 17.13 (0.8281%)
## d = -1721.46 +/- 13.83 (0.8036%)
## e = 530.38 +/- 4.168 (0.7859%)
|
