r23451 - /branches/disp_speed/lib/dispersion/ns_cpmg_2site_3d.py -- May 27, 2014 - 00:08

Author: tlinnet
Date: Tue May 27 00:08:08 2014
New Revision: 23451

URL: http://svn.gna.org/viewcvs/relax?rev=23451&view=rev
Log:
Critical fix for the math domain catching of model NS CPMG 2site 3D.

This was discovered with the added 8 unit tests demonstrating edge case 'no 
Rex' failures.

This follows from the ideas in the post 
http://article.gmane.org/gmane.science.nmr.relax.devel/5858.
This is related to: task #7793: (https://gna.org/task/?7793) Speed-up of 
dispersion models.

This is to implement catching of math domain errors, before they occur.

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

Modified: branches/disp_speed/lib/dispersion/ns_cpmg_2site_3d.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/disp_speed/lib/dispersion/ns_cpmg_2site_3d.py?rev=23451&r1=23450&r2=23451&view=diff
==============================================================================
--- branches/disp_speed/lib/dispersion/ns_cpmg_2site_3d.py      (original)
+++ branches/disp_speed/lib/dispersion/ns_cpmg_2site_3d.py      Tue May 27 
00:08:08 2014
@@ -103,6 +103,12 @@
     @type power:            numpy int16, rank-1 array
     """
 
+    # Catch parameter values that will result in no exchange, returning flat 
R2eff = R20 lines (when kex = 0.0, k_AB = 0.0).
+    if dw == 0.0 or pA == 1.0 or k_AB == 0.0:
+        for i in range(num_points):
+            back_calc[i] = r20a
+        return
+
     # The matrix R that contains all the contributions to the evolution, 
i.e. relaxation, exchange and chemical shift evolution.
     R = rcpmg_3d(R1A=r10a, R1B=r10b, R2A=r20a, R2B=r20b, pA=pA, pB=pB, 
dw=dw, k_AB=k_AB, k_BA=k_BA)
 
@@ -123,6 +129,8 @@
         # The next lines calculate the R2eff using a two-point 
approximation, i.e. assuming that the decay is mono-exponential.
         Mx = fabs(Mint[1] / pA)
         if Mx <= 0.0 or isNaN(Mx):
-            back_calc[i] = 1e99
+            for i in range(num_points):
+                back_calc[i] = r20a
+            return
         else:
             back_calc[i]= -inv_tcpmg * log(Mx)