diff --git a/mrdna/segmentmodel.py b/mrdna/segmentmodel.py index c4f125c098dcf275ff7cdf3aa0253b458a7796aa..e679a4527bc263eab6fd5ae4e8bcfebc96ac3a50 100644 --- a/mrdna/segmentmodel.py +++ b/mrdna/segmentmodel.py @@ -2461,6 +2461,15 @@ class SegmentModel(ArbdModel): parent.add_dihedral(o1,b1,b2,o2, pot) + def count_crossovers(beads): + count = 0 + for b in beads: + for l in b.locations: + if l.connection is not None: + if l.connection.type_ == "crossover": + count += 1 + return count + def add_local_crossover_strand_orientation_potential(b1,b2, b1_on_fwd_strand): """ Adds a dihedral angle potential so bead b2 at opposite @@ -2559,7 +2568,13 @@ class SegmentModel(ArbdModel): ## TODO?: Check length-dependence of this potential if a is not None: - k = k_xover_angle( dists[b][a]+dists[c][d] ) + k_intrinsic = 0.00086 + k = [1/k for k in (4*k_xover_angle( dists[b][a] ), + k_intrinsic, + 4*k_xover_angle( dists[c][d] ))] + k = sum(k)**-1 + k = k * count_crossovers((b1,b2))/4 + pot = self.get_dihedral_potential(k,t0) self.add_dihedral( a,b,c,d, pot ) ... @@ -2612,7 +2627,8 @@ class SegmentModel(ArbdModel): return """ Add bond potential """ - pot = self.get_bond_potential(4,18.5) + k = 1.0 * count_crossovers((b1,b2)) + pot = self.get_bond_potential(k,18.5) self.add_bond(b1,b2, pot) """ Add parallel helices potential, possibly """