diff --git a/arbdmodel.py b/arbdmodel.py
index cbc165e8bff01d6f93a56bfa0ecb3354fe56222b..4405da74d2d8b9c051e7763e131f5cb854606e70 100644
--- a/arbdmodel.py
+++ b/arbdmodel.py
@@ -57,12 +57,13 @@ class Transformable():
return obj
class Parent():
- def __init__(self, children=[]):
+ def __init__(self, children=[], remove_duplicate_bonded_terms=False):
self.children = []
for x in children:
self.add(x)
# self.children = children
-
+
+ self.remove_duplicate_bonded_terms = remove_duplicate_bonded_terms
self.bonds = []
self.angles = []
self.dihedrals = []
@@ -130,31 +131,63 @@ class Parent():
ret = self.bonds
for c in self.children:
if isinstance(c,Parent): ret.extend( c.get_bonds() )
- return ret
+ if self.remove_duplicate_bonded_terms:
+ return list(set(ret))
+ else:
+ return ret
+
def get_angles(self):
ret = self.angles
for c in self.children:
if isinstance(c,Parent): ret.extend( c.get_angles() )
- return ret
+ if self.remove_duplicate_bonded_terms:
+ return list(set(ret))
+ else:
+ return ret
def get_dihedrals(self):
ret = self.dihedrals
for c in self.children:
if isinstance(c,Parent): ret.extend( c.get_dihedrals() )
- return ret
+ if self.remove_duplicate_bonded_terms:
+ return list(set(ret))
+ else:
+ return ret
def get_impropers(self):
ret = self.impropers
for c in self.children:
if isinstance(c,Parent): ret.extend( c.get_impropers() )
- return ret
+ if self.remove_duplicate_bonded_terms:
+ return list(set(ret))
+ else:
+ return ret
def get_exclusions(self):
ret = self.exclusions
for c in self.children:
if isinstance(c,Parent): ret.extend( c.get_exclusions() )
- return ret
+ if self.remove_duplicate_bonded_terms:
+ return list(set(ret))
+ else:
+ return ret
+
+ ## Removed because prohibitively slow
+ # def remove_duplicate_terms(self):
+ # for key in "bonds angles dihedrals impropers exclusions".split():
+ # self.remove_duplicate_item(key)
+
+ # def remove_duplicate_item(self, dict_key, existing=None):
+ # if existing is None: existing = []
+ # ret = [i for i in list(set(self.__dict__[dict_key])) if i not in existing]
+ # self.__dict__[dict_key] = ret
+ # existing.extend(ret)
+ # for c in self.children:
+ # if isinstance(c,Parent):
+ # ret = ret + c.remove_duplicate_item(dict_key, existing)
+ # return ret
+
def __iter__(self):
## TODO: decide if this is the nicest way to do it!
@@ -354,10 +387,11 @@ class Group(Transformable, Parent, Child):
def __init__(self, name=None, children = [], parent=None,
position = np.array((0,0,0)),
orientation = np.array(((1,0,0),(0,1,0),(0,0,1))),
+ remove_duplicate_bonded_terms = False,
):
Transformable.__init__(self, position, orientation)
Child.__init__(self, parent) # Initialize Child first
- Parent.__init__(self, children)
+ Parent.__init__(self, children, remove_duplicate_bonded_terms)
self.name = name
self.isClone = False
@@ -393,9 +427,9 @@ class Group(Transformable, Parent, Child):
class PdbModel(Transformable, Parent):
- def __init__(self, children=[], dimensions=None):
+ def __init__(self, children=[], dimensions=None, remove_duplicate_bonded_terms=False):
Transformable.__init__(self,(0,0,0))
- Parent.__init__(self, children)
+ Parent.__init__(self, children, remove_duplicate_bonded_terms)
self.dimensions = dimensions
self.particles = [p for p in self]
self.cacheInvalid = True
@@ -554,8 +588,8 @@ class PdbModel(Transformable, Parent):
class ArbdModel(PdbModel):
- def __init__(self, children, dimensions=(1000,1000,1000), temperature=291, timestep=50e-6, cutoff=50, decompPeriod=10000, pairlistDistance=None, nonbondedResolution=0.1):
- PdbModel.__init__(self, children, dimensions)
+ def __init__(self, children, dimensions=(1000,1000,1000), temperature=291, timestep=50e-6, cutoff=50, decompPeriod=10000, pairlistDistance=None, nonbondedResolution=0.1, remove_duplicate_bonded_terms=True):
+ PdbModel.__init__(self, children, dimensions, remove_duplicate_bonded_terms)
self.temperature = temperature
self.timestep = timestep
diff --git a/cadnano_segments.py b/cadnano_segments.py
index 1170f1243edbf18fa871a3ef7997ef61f6f7cfa2..ba2199e9f75203e8d0ebabf4d91beaa824715f07 100644
--- a/cadnano_segments.py
+++ b/cadnano_segments.py
@@ -531,8 +531,8 @@ def run_simulation_protocol( output_name, job_id, json_data,
remove_long_bonds=False,
gpu = 0,
directory=None,
- coarse_steps = 1e5+1,
- fine_steps = 1e5+1
+ coarse_steps = 1e7+1,
+ fine_steps = 1e7+1
):
ret = None
diff --git a/run.py b/run.py
index 1e2ff54f38dbdbab4c00055e2c44ed9d4aaaa560..71f2f4f118a0a15fe6a182b9a9c3022ce51dd7ab 100644
--- a/run.py
+++ b/run.py
@@ -19,8 +19,8 @@ if __name__ == '__main__':
print("WARNING: skipping unreadable json file {}".format(f))
continue
- run_simulation_protocol( out, "job-"+out+"-", data, gpu=1 )
- # try:
- # run_simulation_protocol( out, "job-"+out+"-", data, gpu=0 )
- # except:
- # print("ERROR: failed to simulate {}".format(f))
+ # run_simulation_protocol( out, "job-"+out+"-", data, gpu=0 )
+ try:
+ run_simulation_protocol( out, "job-"+out+"-", data, gpu=0 )
+ except:
+ print("WARNING: failed to simulate {}".format(f))
diff --git a/segmentmodel.py b/segmentmodel.py
index 48886590a8a336a83cc8ff923717a45c5acc8cc9..c6b97bbb20c552ea9672435b79ad350ee5cdc459 100644
--- a/segmentmodel.py
+++ b/segmentmodel.py
@@ -735,7 +735,7 @@ class Segment(ConnectableElement, Group):
## Stop if we found the 3prime end
if l.on_fwd_strand == is_fwd and l.type_ == "3prime":
- print(" found end at",l)
+ # print(" found end at",l)
return pos, None, None, None, None
## Check location connections
@@ -745,15 +745,15 @@ class Segment(ConnectableElement, Group):
## Found a location on the same strand?
if l.on_fwd_strand == is_fwd:
- print(" passing through",l)
- print("from {}, connection {} to {}".format(nt_pos,l,B))
+ # print(" passing through",l)
+ # print("from {}, connection {} to {}".format(nt_pos,l,B))
Bpos = B.get_nt_pos()
return pos, B.container, Bpos, B.on_fwd_strand, 0.5
## Stop at other strand crossovers so basepairs line up
elif c.type_ == "crossover":
if nt_pos == pos: continue
- print(" pausing at",l)
+ # print(" pausing at",l)
return pos, l.container, pos+(2*is_fwd-1), is_fwd, 0
import pdb
@@ -805,14 +805,15 @@ class Segment(ConnectableElement, Group):
self.beads = []
if True:
- print("WARNING: DEBUG")
+ ## TODO: remove this if duplicates are never found
+ # print("Searching for duplicate particles...")
## Remove duplicates, preserving order
tmp = []
for c in new_children:
if c not in tmp:
tmp.append(c)
else:
- print(" duplicate particle found!")
+ print(" DUPLICATE PARTICLE FOUND!")
new_children = tmp
for b in new_children:
@@ -878,7 +879,7 @@ class Segment(ConnectableElement, Group):
# import pdb
# pdb.set_trace()
- print(" %s working on %d to %d" % (self.name, eb1.position[2], eb2.position[2]))
+ # print(" %s working on %d to %d" % (self.name, eb1.position[2], eb2.position[2]))
e_ds = eb2.get_contour_position(self) - eb1.get_contour_position(self)
num_beads = self._get_num_beads( e_ds, max_basepairs_per_bead, max_nucleotides_per_bead )
ds = e_ds / (num_beads+1)
@@ -1636,7 +1637,7 @@ class SegmentModel(ArbdModel):
# if s1.name in ("49-3","49-4") and s2.name in ("49-3","49-4"):
# if s1.name in ("51-2","51-3") and s2.name in ("51-2","51-3"):
# pdb.set_trace()
- print("Working on {}".format(c))
+ # print("Working on {}".format(c))
## TODO be more elegant!
# if isinstance(s1, DoubleStrandedSegment) and isinstance(s2, DoubleStrandedSegment) and A.on_fwd_strand == False: continue
if isinstance(s1, DoubleStrandedSegment) and isinstance(s2, DoubleStrandedSegment) and A.on_fwd_strand == False: continue
@@ -2039,7 +2040,9 @@ class SegmentModel(ArbdModel):
# print( n1.idx,n2.idx,n3.idx,n4.idx,k,t0,sep )
pot = self.get_dihedral_potential(k,t0)
self.add_dihedral( n1,n2,n3,n4, pot )
-
+
+ # ## remove duplicate potentials; ## TODO ensure that they aren't added twice in the first place?
+ # self.remove_duplicate_terms()
def walk_through_helices(segment, direction=1, processed_segments=None):
"""
@@ -2182,7 +2185,7 @@ class SegmentModel(ArbdModel):
if locs is None: continue
# for pos, is_fwd in locs:
for l in locs:
- print("Tracing",l)
+ # print("Tracing",l)
# TODOTODO
pos = seg.contour_to_nt_pos(l.address, round_nt=True)
is_fwd = l.on_fwd_strand