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
|
from __future__ import print_function
class LAMMPSPairPotential(object):
def __init__(self):
self.pmap=dict()
self.units='lj'
def map_coeff(self,name,ltype):
self.pmap[ltype]=name
def check_units(self,units):
if (units != self.units):
raise Exception("Conflicting units: %s vs. %s" % (self.units,units))
class LJCutFourMol(LAMMPSPairPotential):
def __init__(self):
super(LJCutFourMol,self).__init__()
self.units='real'
self.coeff={'1': {'1': (57220.458984375, 117.1875, 4768.37158203125, 19.53125), '2': (396.00240328788755, 6.89349609375, 33.00020027399063, 1.148916015625), '3': (275682.72513361875, 257.2233543675, 22973.560427801563, 42.870559061250006), '4': (193062.6401184066, 200.31803764867476, 16088.553343200549, 33.38633960811246), '5': (219401.24669848208, 213.54552859539564, 18283.437224873505, 35.590921432565935)}, '2': {'1': (396.00240328788755, 6.89349609375, 33.00020027399063, 1.148916015625), '2': (0.24, 0.12, 0.02, 0.02), '3': (3530.7972054655893, 20.583869040000007, 294.2331004554658, 3.430644840000001), '4': (5.859375e-05, 0.001875, 4.8828125e-06, 0.0003125), '5': (2726.0378943788915, 16.831463637644212, 227.1698245315743, 2.8052439396073683)}, '3': {'1': (275682.72513361875, 257.2233543675, 22973.560427801563, 42.870559061250006), '2': (3530.7972054655893, 20.583869040000007, 294.2331004554658, 3.430644840000001), '3': (1106804.6444225737, 515.3960755200002, 92233.72036854782, 85.89934592000003), '4': (793466.3439000695, 406.10205934924875, 66122.19532500578, 67.68367655820813), '5': (889052.2891953762, 429.8674775516488, 74087.69076628136, 71.64457959194146)}, '4': {'1': (193062.6401184066, 200.31803764867476, 16088.553343200549, 33.38633960811246), '2': (5.859375e-05, 0.001875, 4.8828125e-06, 0.0003125), '3': (793466.3439000695, 406.10205934924875, 66122.19532500578, 67.68367655820813), '4': (567117.2043277561, 319.50132516, 47259.767027313, 53.250220860000006), '5': (636596.9403948045, 338.50767978151515, 53049.74503290038, 56.41794663025252)}, '5': {'1': (219401.24669848208, 213.54552859539564, 18283.437224873505, 35.590921432565935), '2': (2726.0378943788915, 16.831463637644212, 227.1698245315743, 2.8052439396073683), '3': (889052.2891953762, 429.8674775516488, 74087.69076628136, 71.64457959194146), '4': (636596.9403948045, 338.50767978151515, 53049.74503290038, 56.41794663025252), '5': (713801.551928242, 358.44703841304585, 59483.46266068683, 59.74117306884097)}}
def compute_force(self,rsq,itype,jtype):
coeff = self.coeff[self.pmap[itype]][self.pmap[jtype]]
r2inv = 1.0/rsq
r6inv = r2inv*r2inv*r2inv
lj1 = coeff[0]
lj2 = coeff[1]
return (r6inv * (lj1*r6inv - lj2))*r2inv
def compute_energy(self,rsq,itype,jtype):
coeff = self.coeff[self.pmap[itype]][self.pmap[jtype]]
r2inv = 1.0/rsq
r6inv = r2inv*r2inv*r2inv
lj3 = coeff[2]
lj4 = coeff[3]
return (r6inv * (lj3*r6inv - lj4))
|
