diff --git a/arbdmodel/__init__.py b/arbdmodel/__init__.py
index 03ba1d2ef781b16fb27d53a7f2a2d6c37ce1b6b0..981b3ccc210034d48ac4b760eb9e9a01d0daf0d3 100644
--- a/arbdmodel/__init__.py
+++ b/arbdmodel/__init__.py
@@ -417,6 +417,10 @@ class PointParticle(Transformable, Child):
segname = p.segname
except:
segname = "A"
+ try:
+ chain = p.chain
+ except:
+ chain = "A"
try:
resname = p.resname
except:
@@ -438,7 +442,7 @@ class PointParticle(Transformable, Child):
data = dict(segname = segname,
resname = resname,
name = str(p.name)[:4],
- chain = "A",
+ chain = chain[0],
resid = int(resid),
idx = p.idx+1,
type = p.type_.name[:4],
diff --git a/arbdmodel/abstract_polymer.py b/arbdmodel/abstract_polymer.py
index 3a25e1169435359859ca455af7c7ade53a51439d..732a9562645669297d7b4ea6b1218fe302815faa 100644
--- a/arbdmodel/abstract_polymer.py
+++ b/arbdmodel/abstract_polymer.py
@@ -174,6 +174,9 @@ class PolymerSection(ConnectableElement):
if 'segname' not in kwargs:
self.segname = name
+ for key,val in kwargs.items():
+ self.__dict__[key] = val
+
self.num_monomers = int(num_monomers)
self.monomer_length = monomer_length
diff --git a/arbdmodel/hps_polymer_model.py b/arbdmodel/hps_polymer_model.py
index f0cd7237c88bb9703072dcc8d7020c269192679e..c8a63274937b74b0986e610cfbe5919d6ccb87c3 100644
--- a/arbdmodel/hps_polymer_model.py
+++ b/arbdmodel/hps_polymer_model.py
@@ -17,7 +17,7 @@ _types = dict(
mass = 71.08,
charge = 0,
sigma = 5.04,
- lambda_ = 0.730,
+ lambda_ = 0.72973,
),
R = ParticleType("ARG",
mass = 156.2,
@@ -29,67 +29,67 @@ _types = dict(
mass = 114.1,
charge = 0,
sigma = 5.68,
- lambda_ = 0.432,
+ lambda_ = 0.432432,
),
D = ParticleType("ASP",
mass = 115.1,
charge = -1,
sigma = 5.58,
- lambda_ = 0.378,
+ lambda_ = 0.378378,
),
C = ParticleType("CYS",
mass = 103.1,
charge = 0,
sigma = 5.48,
- lambda_ = 0.595,
+ lambda_ = 0.594595,
),
Q = ParticleType("GLN",
mass = 128.1,
charge = 0,
sigma = 6.02,
- lambda_ = 0.514,
+ lambda_ = 0.513514,
),
E = ParticleType("GLU",
mass = 129.1,
charge = -1,
sigma = 5.92,
- lambda_ = 0.459,
+ lambda_ = 0.459459,
),
G = ParticleType("GLY",
mass = 57.05,
charge = 0,
sigma = 4.5,
- lambda_ = 0.649,
+ lambda_ = 0.648649,
),
H = ParticleType("HIS",
mass = 137.1,
charge = 0.5,
sigma = 6.08,
- lambda_ = 0.514,
+ lambda_ = 0.513514,
),
I = ParticleType("ILE",
mass = 113.2,
charge = 0,
sigma = 6.18,
- lambda_ = 0.973,
+ lambda_ = 0.972973,
),
L = ParticleType("LUE",
mass = 113.2,
charge = 0,
sigma = 6.18,
- lambda_ = 0.973,
+ lambda_ = 0.972973,
),
K = ParticleType("LYS",
mass = 128.2,
charge = 1,
sigma = 6.36,
- lambda_ = 0.514,
+ lambda_ = 0.513514,
),
M = ParticleType("MET",
mass = 131.2,
charge = 0,
sigma = 6.18,
- lambda_ = 0.838,
+ lambda_ = 0.837838,
),
F = ParticleType("PHE",
mass = 147.2,
@@ -107,33 +107,35 @@ _types = dict(
mass = 87.08,
charge = 0,
sigma = 5.18,
- lambda_ = 0.595,
+ lambda_ = 0.594595,
),
T = ParticleType("THR",
mass = 101.1,
charge = 0,
sigma = 5.62,
- lambda_ = 0.676,
+ lambda_ = 0.675676,
),
W = ParticleType("TRP",
mass = 186.2,
charge = 0,
sigma = 6.78,
- lambda_ = 0.946,
+ lambda_ = 0.945946,
),
Y = ParticleType("TYR",
mass = 163.2,
charge = 0,
sigma = 6.46,
- lambda_ = 0.865,
+ lambda_ = 0.864865,
),
V = ParticleType("VAL",
mass = 99.07,
charge = 0,
sigma = 5.86,
- lambda_ = 0.892,
+ lambda_ = 0.891892,
)
)
+for k,t in _types.items():
+ t.resname = t.name
class HpsNonbonded(NonbondedScheme):
def __init__(self, debye_length=10, resolution=0.1, rMin=0):
@@ -193,6 +195,16 @@ class HpsBeadsFromPolymer(Group):
self.polymer = polymer
self.sequence = sequence
+
+ for prop in ('segname','chain'):
+ if prop not in kwargs:
+ # import pdb
+ # pdb.set_trace()
+ try:
+ self.__dict__[prop] = polymer.__dict__[prop]
+ except:
+ pass
+
if len(sequence) != polymer.num_monomers:
raise ValueError("Length of sequence does not match length of polymer")
Group.__init__(self, **kwargs)
@@ -208,8 +220,13 @@ class HpsBeadsFromPolymer(Group):
r = self.polymer.contour_to_position(c)
s = self.sequence[i]
- bead = PointParticle(_types[s], r )
+ bead = PointParticle(_types[s], r,
+ name = s,
+ resid = i+1)
self.add(bead)
+ # import pdb
+ # pdb.set_trace()
+ # continue
## Two consecutive nts
for i in range(len(self.children)-1):
@@ -236,25 +253,23 @@ class HpsModel(ArbdModel):
[debye_length]: angstroms
[damping_coefficient]: ns
"""
-
kwargs['timestep'] = 10e-6
kwargs['cutoff'] = max(4*debye_length,20)
+ if 'decompPeriod' not in kwargs:
+ kwargs['decompPeriod'] = 1000
+
""" Assign sequences """
if sequences is None:
raise NotImplementedError("HpsModel must be provided a sequences argument")
self.polymer_group = AbstractPolymerGroup(polymers)
- self.peptides = [HpsBeadsFromPolymer(p,s)
- for p,s in zip(self.polymer_group.polymers,sequences)]
-
- ArbdModel.__init__(self, self.peptides, **kwargs)
+ self.sequences = sequences
+ ArbdModel.__init__(self, [], **kwargs)
""" Update type diffusion coefficients """
- all_types = [t for key,t in _types.items()]
- for t in all_types:
- t.damping_coefficient = damping_coefficient
- # t.diffusivity = 831447.2 * temperature / (t.mass * damping_coefficient)
+ self.types = all_types = [t for key,t in _types.items()]
+ self.set_damping_coefficient( damping_coefficient )
""" Set up nonbonded interactions """
nonbonded = HpsNonbonded(debye_length)
@@ -268,12 +283,28 @@ class HpsModel(ArbdModel):
self.generate_beads()
def update_splines(self, coords):
- ...
-
+ i = 0
+ for p in self.polymer_group.polymers:
+ n = p.num_monomers
+ p.set_splines(np.linspace(0,1,n), coords[i:i+n])
+ i += n
+
+ self.clear_all()
+ self.generate_beads()
+ ## TODO Apply restraints, etc
+
def generate_beads(self):
+ self.peptides = [HpsBeadsFromPolymer(p,s)
+ for p,s in zip(self.polymer_group.polymers,self.sequences)]
+
for s in self.peptides:
+ self.add(s)
s._generate_beads()
+ def set_damping_coefficient(self, damping_coefficient):
+ for t in self.types:
+ t.damping_coefficient = damping_coefficient
+ # t.diffusivity = 831447.2 * temperature / (t.mass * damping_coefficient)
if __name__ == "__main__":