mailr23219 - /branches/disp_speed/lib/dispersion/tsmfk01.py


Others Months | Index by Date | Thread Index
>>   [Date Prev] [Date Next] [Thread Prev] [Thread Next]

Header


Content

Posted by tlinnet on May 19, 2014 - 00:51:
Author: tlinnet
Date: Mon May 19 00:51:18 2014
New Revision: 23219

URL: http://svn.gna.org/viewcvs/relax?rev=23219&view=rev
Log:
Speed-up of model TSMFK01.

task #7793: (https://gna.org/task/?7793) Speed-up of dispersion models.

This is won by not checking single values in the R2eff array for math domain
errors, but calculating all steps, and in one single round check for finite 
values.
If just one non-finite value is found, the whole array is returned with a 
large
penalty of 1e100.

This makes all calculations be the fastest numpy array way.

Modified:
    branches/disp_speed/lib/dispersion/tsmfk01.py

Modified: branches/disp_speed/lib/dispersion/tsmfk01.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/disp_speed/lib/dispersion/tsmfk01.py?rev=23219&r1=23218&r2=23219&view=diff
==============================================================================
--- branches/disp_speed/lib/dispersion/tsmfk01.py       (original)
+++ branches/disp_speed/lib/dispersion/tsmfk01.py       Mon May 19 00:51:18 
2014
@@ -67,7 +67,7 @@
 """
 
 # Python module imports.
-from math import sin
+from numpy import sin, isfinite, sum
 
 
 def r2eff_TSMFK01(r20a=None, dw=None, k_AB=None, tcp=None, back_calc=None, 
num_points=None):
@@ -90,24 +90,20 @@
     @type num_points:       int
     """
 
-    # Loop over the time points, back calculating the R2eff values.
+    # Denominator.
+    denom = dw * tcp
+
+    # The numerator.
+    numer = sin(denom)
+
+    # Calculate R2eff.
+    R2eff = r20a + k_AB - k_AB * numer / denom
+
+    # Catch errors, taking a sum over array is the fastest way to check for
+    # +/- inf (infinity) and nan (not a number).
+    if not isfinite(sum(R2eff)):
+        R2eff = array([1e100]*num_points)
+
+    # Parse back the value to update the back_calc class object.
     for i in range(num_points):
-        # Denominator.
-        denom = dw * tcp[i]
-
-        # The numerator.
-        numer = sin(denom)
-
-        # Catch zeros (to avoid pointless mathematical operations).
-        if numer == 0.0:
-            back_calc[i] = r20a + k_AB
-            continue
-
-        # Avoid divide by zero.
-        if denom == 0.0:
-            back_calc[i] = 1e100
-            continue
-
-        # The full formula.
-        else:
-            back_calc[i] = r20a + k_AB - k_AB * numer / denom
+        back_calc[i] = R2eff[i]




Related Messages


Powered by MHonArc, Updated Mon May 19 01:00:03 2014