r23278 - in /branches/disp_speed: lib/dispersion/mmq_cr72.py target_functions/relax_disp.py -- May 20, 2014 - 22:29

Author: tlinnet
Date: Tue May 20 22:29:54 2014
New Revision: 23278

URL: http://svn.gna.org/viewcvs/relax?rev=23278&view=rev
Log:
Math-domain catching for model MMQ CR72.

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

This is to implement catching of math domain errors, before they occur.
These can be found via the --numpy-raise function to the systemtests.

To make the code look clean, the class object "back_calc" is no longer
being updated per time point, but is updated in the relax_disp target 
function in
one go.

Modified:
    branches/disp_speed/lib/dispersion/mmq_cr72.py
    branches/disp_speed/target_functions/relax_disp.py

Modified: branches/disp_speed/lib/dispersion/mmq_cr72.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/disp_speed/lib/dispersion/mmq_cr72.py?rev=23278&r1=23277&r2=23278&view=diff
==============================================================================
--- branches/disp_speed/lib/dispersion/mmq_cr72.py      (original)
+++ branches/disp_speed/lib/dispersion/mmq_cr72.py      Tue May 20 22:29:54 
2014
@@ -47,10 +47,10 @@
 """
 
 # Python module imports.
-from numpy import arccosh, array, cos, cosh, isfinite, log, sin, sqrt, sum
+from numpy import arccosh, array, cos, cosh, isfinite, log, max, sin, sqrt, 
sum
 
 
-def r2eff_mmq_cr72(r20=None, pA=None, pB=None, dw=None, dwH=None, kex=None, 
k_AB=None, k_BA=None, cpmg_frqs=None, inv_tcpmg=None, tcp=None, 
back_calc=None, num_points=None, power=None):
+def r2eff_mmq_cr72(r20=None, pA=None, pB=None, dw=None, dwH=None, kex=None, 
k_AB=None, k_BA=None, cpmg_frqs=None, inv_tcpmg=None, tcp=None, 
num_points=None, power=None):
     """The CR72 model extended to MMQ CPMG data.
 
     This function calculates and stores the R2eff values.
@@ -78,9 +78,7 @@
     @type inv_tcpmg:        float
     @keyword tcp:           The tau_CPMG times (1 / 4.nu1).
     @type tcp:              numpy rank-1 float array
-    @keyword back_calc:     The array for holding the back calculated R2eff 
values.  Each element corresponds to one of the CPMG nu1 frequencies.
-    @type back_calc:        numpy rank-1 float array
-    @keyword num_points:    The number of points on the dispersion curve, 
equal to the length of the tcp and back_calc arguments.
+    @keyword num_points:    The number of points on the dispersion curve, 
equal to the length of the tcp.
     @type num_points:       int
     @keyword power:         The matrix exponential power array.
     @type power:            numpy int16, rank-1 array
@@ -122,8 +120,17 @@
     etapos_part = eta_scale * sqrt(Psi + sqrt_psi2_zeta2)
     etaneg_part = eta_scale * sqrt(-Psi + sqrt_psi2_zeta2)
 
-    # The full eta+/- values.
+    # The full eta+ values.
     etapos = etapos_part / cpmg_frqs
+
+    # Catch math domain error of cosh(val > 710).
+    # This is when etapos > 710.
+    if max(etapos) > 700:
+        R2eff = array([1e100]*num_points)
+
+        return R2eff
+
+    # The full eta - values.
     etaneg = etaneg_part / cpmg_frqs
 
     # The mD value.
@@ -147,6 +154,4 @@
     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):
-        back_calc[i] = R2eff[i]
+    return R2eff

Modified: branches/disp_speed/target_functions/relax_disp.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/disp_speed/target_functions/relax_disp.py?rev=23278&r1=23277&r2=23278&view=diff
==============================================================================
--- branches/disp_speed/target_functions/relax_disp.py  (original)
+++ branches/disp_speed/target_functions/relax_disp.py  Tue May 20 22:29:54 
2014
@@ -1335,7 +1335,7 @@
                         aliased_dwH = dw_frq
 
                     # Back calculate the R2eff values.
-                    r2eff_mmq_cr72(r20=R20[r20_index], pA=pA, pB=pB, 
dw=aliased_dw, dwH=aliased_dwH, kex=kex, k_AB=k_AB, k_BA=k_BA, 
cpmg_frqs=self.cpmg_frqs[ei][mi][0], inv_tcpmg=self.inv_relax_times[ei][mi], 
tcp=self.tau_cpmg[ei][mi], back_calc=self.back_calc[ei][si][mi][0], 
num_points=self.num_disp_points[ei][si][mi][0], power=self.power[ei][mi])
+                    self.back_calc[ei][si][mi][0] = 
r2eff_mmq_cr72(r20=R20[r20_index], pA=pA, pB=pB, dw=aliased_dw, 
dwH=aliased_dwH, kex=kex, k_AB=k_AB, k_BA=k_BA, 
cpmg_frqs=self.cpmg_frqs[ei][mi][0], inv_tcpmg=self.inv_relax_times[ei][mi], 
tcp=self.tau_cpmg[ei][mi], num_points=self.num_disp_points[ei][si][mi][0], 
power=self.power[ei][mi])
 
                     # For all missing data points, set the back-calculated 
value to the measured values so that it has no effect on the chi-squared 
value.
                     for di in range(self.num_disp_points[ei][si][mi][0]):