--- /home/sjplimp/oldpizza/src/pair.py 2005-09-07 10:13:09.000000000 -0600
+++ /home/sjplimp/pizza/src/pair.py 2005-09-13 12:46:56.000000000 -0600
@@ -13,16 +13,29 @@
docstr = """
p = pair("lj/charmm/coul/charmm") create pair object for specific pair style
+ available styles: lj/cut, lj/cut/coul/cut, lj/charmm/coul/charmm
+
p.coeff(d) extract pairwise coeffs from data object
p.init(cut1,cut2,...) setup based on coeffs and cutoffs
-
+
+ init args are specific to pair style:
+ lj/cut = cutlj
+ lj/cut/coul/cut = cutlj,cut_coul (cut_coul optional)
+ lj/charmm/coul/charmm = cutlj_inner,cutlj,cutcoul_inner,cut_coul
+ (last 2 optional)
+
e_vdwl,e_coul = p.single(rsq,itype,jtype,q1,q2,...) compute LJ/Coul energy
pairwise energy between 2 atoms at distance rsq with their attributes
+ args are specific to pair style:
+ lj/cut = rsq,itype,jtype
+ lj/cut/coul/cut = rsq,itype,jtype,q1,q2
+ lj/charmm/coul/charmm = rsq,itype,jtype,q1,q2
"""
# History
# 8/05, Steve Plimpton and Paul Crozier (SNL): original version
+# 9/05, Paul Crozier (SNL): added lj/cut and lj/cut/coul/cut
# ToDo list
@@ -42,7 +55,15 @@
# map generic methods to style-specific methods
- if style == "lj/charmm/coul/charmm":
+ if style == "lj/cut":
+ self.coeff_func = self.coeff_lj_cut
+ self.init_func = self.init_lj_cut
+ self.single_func = self.single_lj_cut
+ elif style == "lj/cut/coul/cut":
+ self.coeff_func = self.coeff_lj_cut_coul_cut
+ self.init_func = self.init_lj_cut_coul_cut
+ self.single_func = self.single_lj_cut_coul_cut
+ elif style == "lj/charmm/coul/charmm":
self.coeff_func = self.coeff_lj_charmm_coul_charmm
self.init_func = self.init_lj_charmm_coul_charmm
self.single_func = self.single_lj_charmm_coul_charmm
@@ -68,6 +89,105 @@
return self.single_func(list)
# --------------------------------------------------------------------
+ # --------------------------------------------------------------------
+ # lj/cut methods
+
+ def coeff_lj_cut(self,data):
+ epsilon = data.get("Pair Coeffs",2)
+ sigma = data.get("Pair Coeffs",3)
+ ntypes = len(epsilon)
+
+ self.lj3 = []
+ self.lj4 = []
+ for i in xrange(ntypes):
+ self.lj3.append(ntypes * [0])
+ self.lj4.append(ntypes * [0])
+ for j in xrange(ntypes):
+ epsilon_ij = sqrt(epsilon[i]*epsilon[j])
+ sigma_ij = 0.5 * (sigma[i] + sigma[j])
+ self.lj3[i][j] = 4.0 * epsilon_ij * pow(sigma_ij,12.0);
+ self.lj4[i][j] = 4.0 * epsilon_ij * pow(sigma_ij,6.0);
+
+ # --------------------------------------------------------------------
+ # args = cutlj
+
+ def init_lj_cut(self,list):
+ cut_lj = list[0]
+ self.cut_ljsq = cut_lj*cut_lj
+
+ # --------------------------------------------------------------------
+ # args = rsq,itype,jtype
+
+ def single_lj_cut(self,list):
+ rsq = list[0]
+ itype = list[1]
+ jtype = list[2]
+
+ r2inv = 1.0/rsq
+
+ if rsq < self.cut_ljsq:
+ r6inv = r2inv*r2inv*r2inv
+ eng_vdwl = r6inv*(self.lj3[itype][jtype]*r6inv-self.lj4[itype][jtype])
+ else: eng_vdwl = 0.0
+
+ return eng_vdwl
+
+ # --------------------------------------------------------------------
+ # --------------------------------------------------------------------
+ # lj/cut/coul/cut methods
+
+ def coeff_lj_cut_coul_cut(self,data):
+ epsilon = data.get("Pair Coeffs",2)
+ sigma = data.get("Pair Coeffs",3)
+ ntypes = len(epsilon)
+
+ self.lj3 = []
+ self.lj4 = []
+ for i in xrange(ntypes):
+ self.lj3.append(ntypes * [0])
+ self.lj4.append(ntypes * [0])
+ for j in xrange(ntypes):
+ epsilon_ij = sqrt(epsilon[i]*epsilon[j])
+ sigma_ij = 0.5 * (sigma[i] + sigma[j])
+ self.lj3[i][j] = 4.0 * epsilon_ij * pow(sigma_ij,12.0);
+ self.lj4[i][j] = 4.0 * epsilon_ij * pow(sigma_ij,6.0);
+
+ # --------------------------------------------------------------------
+ # args = cutlj, cut_coul (cut_coul optional)
+
+ def init_lj_cut_coul_cut(self,list):
+ self.qqr2e = 332.0636 # convert energy to kcal/mol
+ cut_lj = list[0]
+ self.cut_ljsq = cut_lj*cut_lj
+
+ if len(list) == 1: cut_coul = cut_lj
+ else: cut_coul = list[1]
+ self.cut_coulsq = cut_coul*cut_coul
+
+ # --------------------------------------------------------------------
+ # args = rsq,itype,jtype,q1,q2
+
+ def single_lj_cut_coul_cut(self,list):
+ rsq = list[0]
+ itype = list[1]
+ jtype = list[2]
+ q1 = list[3]
+ q2 = list[4]
+
+ r2inv = 1.0/rsq
+
+ if rsq < self.cut_coulsq: eng_coul = self.qqr2e * q1*q2*sqrt(r2inv)
+ else: eng_coul = 0.0
+
+ if rsq < self.cut_ljsq:
+ r6inv = r2inv*r2inv*r2inv
+ eng_vdwl = r6inv*(self.lj3[itype][jtype]*r6inv-self.lj4[itype][jtype])
+ else: eng_vdwl = 0.0
+
+ return eng_coul,eng_vdwl
+
+ # --------------------------------------------------------------------
+ # --------------------------------------------------------------------
# lj/charmm/coul/charmm methods
def coeff_lj_charmm_coul_charmm(self,data):
--- /home/sjplimp/oldpizza/doc/pair.txt 2005-09-07 10:13:10.000000000 -0600
+++ /home/sjplimp/pizza/doc/pair.txt 2005-09-13 12:49:29.000000000 -0600
@@ -19,10 +19,10 @@
analysis script to compute energies between groups of atoms from a
LAMMPS snapshot file.
-The pair constructor specifies the force field style. Currently only
-the lj/charmm/coul/charmm style is included in pair.py, but new styles
-can easily be added. Code from the LAMMPS pair*.cpp file needs to be
-re-coded in Python to make this work.
+The pair constructor specifies the force field style. Only some of
+the LAMMPS pair styles are currently included in this tool, but new
+styles can easily be added. Code from the LAMMPS pair*.cpp file needs
+to be re-coded in Python to make this work.
The coeff() method reads the pairwise coefficients for the force field
from a data file object (see the "data"_data.html tool). The init()
@@ -37,12 +37,24 @@
p = pair("lj/charmm/coul/charmm") create pair object for specific pair style :pre
+ available styles: lj/cut, lj/cut/coul/cut, lj/charmm/coul/charmm :pre
+
p.coeff(d) extract pairwise coeffs from data object
p.init(cut1,cut2,...) setup based on coeffs and cutoffs :pre
+ init args are specific to pair style:
+ lj/cut = cutlj
+ lj/cut/coul/cut = cutlj,cut_coul (cut_coul optional)
+ lj/charmm/coul/charmm = cutlj_inner,cutlj,cutcoul_inner,cut_coul
+ (last 2 optional) :pre
+
e_vdwl,e_coul = p.single(rsq,itype,jtype,q1,q2,...) compute LJ/Coul energy :pre
- pairwise energy between 2 atoms at distance rsq with their attributes :pre
+ pairwise energy between 2 atoms at distance rsq with their attributes
+ args are specific to pair style:
+ lj/cut = rsq,itype,jtype
+ lj/cut/coul/cut = rsq,itype,jtype,q1,q2
+ lj/charmm/coul/charmm = rsq,itype,jtype,q1,q2 :pre
[Related tools:]
--- /home/sjplimp/oldpizza/doc/pair.html 2005-09-07 10:13:10.000000000 -0600
+++ /home/sjplimp/pizza/doc/pair.html 2005-09-13 12:49:35.000000000 -0600
@@ -22,10 +22,10 @@
analysis script to compute energies between groups of atoms from a
LAMMPS snapshot file.
-The pair constructor specifies the force field style. Currently only
-the lj/charmm/coul/charmm style is included in pair.py, but new styles
-can easily be added. Code from the LAMMPS pair*.cpp file needs to be
-re-coded in Python to make this work.
+
The pair constructor specifies the force field style. Only some of
+the LAMMPS pair styles are currently included in this tool, but new
+styles can easily be added. Code from the LAMMPS pair*.cpp file needs
+to be re-coded in Python to make this work.
The coeff() method reads the pairwise coefficients for the force field
from a data file object (see the data tool). The init()
@@ -40,12 +40,24 @@
p = pair("lj/charmm/coul/charmm") create pair object for specific pair style
+ available styles: lj/cut, lj/cut/coul/cut, lj/charmm/coul/charmm
+
p.coeff(d) extract pairwise coeffs from data object
p.init(cut1,cut2,...) setup based on coeffs and cutoffs
+ init args are specific to pair style:
+ lj/cut = cutlj
+ lj/cut/coul/cut = cutlj,cut_coul (cut_coul optional)
+ lj/charmm/coul/charmm = cutlj_inner,cutlj,cutcoul_inner,cut_coul
+ (last 2 optional)
+
e_vdwl,e_coul = p.single(rsq,itype,jtype,q1,q2,...) compute LJ/Coul energy
- pairwise energy between 2 atoms at distance rsq with their attributes
+ pairwise energy between 2 atoms at distance rsq with their attributes
+ args are specific to pair style:
+ lj/cut = rsq,itype,jtype
+ lj/cut/coul/cut = rsq,itype,jtype,q1,q2
+ lj/charmm/coul/charmm = rsq,itype,jtype,q1,q2
Related tools: