From 3cfec3cef2923ebd9fbba6b79ce2620a0444a42b Mon Sep 17 00:00:00 2001
From: Chris Maffeo <cmaffeo2@illinois.edu>
Date: Thu, 20 Sep 2018 21:22:33 -0500
Subject: [PATCH] Updated method for spring constant
---
mrdna/segmentmodel.py | 14 +++++++-------
1 file changed, 7 insertions(+), 7 deletions(-)
diff --git a/mrdna/segmentmodel.py b/mrdna/segmentmodel.py
index 8f1b305..f62513c 100644
--- a/mrdna/segmentmodel.py
+++ b/mrdna/segmentmodel.py
@@ -15,6 +15,8 @@ from scipy import interpolate
from .model.CanonicalNucleotideAtoms import canonicalNtFwd, canonicalNtRev, seqComplement
from .model.CanonicalNucleotideAtoms import enmTemplateHC, enmTemplateSQ, enmCorrectionsHC
+from .model.spring_from_lp import k_angle as angle_spring_from_lp
+
# import pdb
"""
TODO:
@@ -2069,10 +2071,7 @@ class SegmentModel(ArbdModel):
kT = 0.58622522 # kcal/mol
if b1.type_.name[0] == "D" and b2.type_.name[0] == "D" and b3.type_.name[0] == "D":
- ## <cos(q)> = exp(-s/Lp) = integrate( x^4 exp(-A x^2) / 2, {x, 0, pi} ) / integrate( x^2 exp(-A x^2), {x, 0, pi} )
- ## <cos(q)> ~ 1 - 3/4A
- ## where A = k_spring / (2 kT)
- k = 1.5 * kT * (1.0 / (1-np.exp(-float(sep)/147))) * 0.00030461742; # kcal_mol/degree^2
+ k = angle_spring_from_lp(sep,147)
if local_twist:
## TODO optimize this paramter
k *= 0.5 # halve because orientation beads have similar springs
@@ -2082,7 +2081,7 @@ class SegmentModel(ArbdModel):
else:
## TODO: get correct number from ssDNA model
- k = 1.5 * kT * (1.0 / (1-np.exp(-float(sep)/3))) * 0.00030461742; # kcal_mol/degree^2
+ k = angle_spring_from_lp(sep,2)
angle = self.get_angle_potential(k,180)
parent.add_angle( b1, b2, b3, angle )
@@ -2126,7 +2125,7 @@ class SegmentModel(ArbdModel):
parent = self._getParent( b1, b2 )
""" Add heuristic 90 degree potential to keep orientation bead orthogonal """
- k = (1.0/2) * 1.5 * kT * (1.0 / (1-np.exp(-float(sep)/147))) * 0.00030461742; # kcal_mol/degree^2
+ k = 0.5*angle_spring_from_lp(sep,147)
pot = self.get_angle_potential(k,90)
parent.add_angle(o1,b1,b2, pot)
parent.add_angle(b1,b2,o2, pot)
@@ -2146,7 +2145,8 @@ class SegmentModel(ArbdModel):
parent.add_dihedral(o1,b1,b2,o2, pot)
def k_angle(sep):
- return 1.5 * kT * (1.0 / (1-np.exp(-float(sep)/147))) * 0.00030461742; # kcal_mol/degree^2
+ return 0.5*angle_spring_from_lp(sep,147)
+
def k_xover_angle(sep):
return 0.5 * k_angle(sep)
--
GitLab