r11163 - /1.3/maths_fns/rdc.py -- April 28, 2010 - 16:35

Author: bugman
Date: Wed Apr 28 16:35:26 2010
New Revision: 11163

URL: http://svn.gna.org/viewcvs/relax?rev=11163&view=rev
Log:
Huge bug fix for the 'population' N-state model.

The probabilities of each state was being set to that of the last state.  
This is fatal for this
model.  Earlier code did not have this problem, but it is unknown when this 
bug was introduced.


Modified:
    1.3/maths_fns/rdc.py

Modified: 1.3/maths_fns/rdc.py
URL: 
http://svn.gna.org/viewcvs/relax/1.3/maths_fns/rdc.py?rev=11163&r1=11162&r2=11163&view=diff
==============================================================================
--- 1.3/maths_fns/rdc.py (original)
+++ 1.3/maths_fns/rdc.py Wed Apr 28 16:35:26 2010
@@ -72,21 +72,17 @@
     # No weights given.
     if weights == None:
         pc = 1.0 / N
-
-    # Loop over the structures.
+        weights = [pc] * N
+
+    # Missing last weight.
+    if len(weights) < N:
+        pN = 1.0 - sum(weights, axis=0)
+        weights = weights.tolist()
+        weights.append(pN)
+
+    # Back-calculate the RDC.
     for c in xrange(N):
-        # The given weights.
-        if weights != None:
-            # Missing last weight.
-            if len(weights) < N: 
-                pc = 1.0 - sum(weights, axis=0)
-
-            # Given weight.
-            else:
-                pc = weights[c]
-
-        # Back-calculate the RDC.
-        val = val + pc * (vect[c, 0]**2 - vect[c, 2]**2)*A[0] + (vect[c, 
1]**2 - vect[c, 2]**2)*A[1] + 2.0*vect[c, 0]*vect[c, 1]*A[2] + 2.0*vect[c, 
0]*vect[c, 2]*A[3] + 2.0*vect[c, 1]*vect[c, 2]*A[4]
+        val = val + weights[c] * (vect[c, 0]**2 - vect[c, 2]**2)*A[0] + 
(vect[c, 1]**2 - vect[c, 2]**2)*A[1] + 2.0*vect[c, 0]*vect[c, 1]*A[2] + 
2.0*vect[c, 0]*vect[c, 2]*A[3] + 2.0*vect[c, 1]*vect[c, 2]*A[4]
 
     # Return the average RDC.
     return dj * val
@@ -151,21 +147,17 @@
     # No weights given.
     if weights == None:
         pc = 1.0 / N
-
-    # Loop over the structures c.
+        weights = [pc] * N
+
+    # Missing last weight.
+    if len(weights) < N:
+        pN = 1.0 - sum(weights, axis=0)
+        weights = weights.tolist()
+        weights.append(pN)
+
+    # Back-calculate the RDC.
     for c in xrange(N):
-        # The given weights.
-        if weights != None:
-            # Missing last weight.
-            if len(weights) < N: 
-                pc = 1.0 - sum(weights, axis=0)
-
-            # Given weight.
-            else:
-                pc = weights[c]
-
-        # Back-calculate the RDC.
-        val = val + pc * dot(vect[c], dot(A, vect[c]))
+        val = val + weights[c] * dot(vect[c], dot(A, vect[c]))
 
     # Return the average RDC.
     return dj * val
@@ -217,21 +209,17 @@
     # No weights given.
     if weights == None:
         pc = 1.0 / N
-
-    # Loop over the structures.
+        weights = [pc] * N
+
+    # Missing last weight.
+    if len(weights) < N:
+        pN = 1.0 - sum(weights, axis=0)
+        weights = weights.tolist()
+        weights.append(pN)
+
+    # Back-calculate the RDC gradient element.
     for c in xrange(N):
-        # The given weights.
-        if weights != None:
-            # Missing last weight.
-            if len(weights) < N: 
-                pc = 1.0 - sum(weights, axis=0)
-
-            # Given weight.
-            else:
-                pc = weights[c]
-
-        # Back-calculate the RDC gradient element.
-        grad = grad + pc * dot(vect[c], dot(dAi_dAmn, vect[c]))
+        grad = grad + weights[c] * dot(vect[c], dot(dAi_dAmn, vect[c]))
 
     # Return the average RDC gradient element.
     return dj * grad