Hi, Again I don't think that the denominator can be zero. Did you see this in the test suite?
Regards, Edward On 20 May 2014 22:29, <[email protected]> wrote: > Author: tlinnet > Date: Tue May 20 22:29:52 2014 > New Revision: 23277 > > URL: http://svn.gna.org/viewcvs/relax?rev=23277&view=rev > Log: > Math-domain catching for model M61. > > 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/m61.py > branches/disp_speed/target_functions/relax_disp.py > > Modified: branches/disp_speed/lib/dispersion/m61.py > URL: > http://svn.gna.org/viewcvs/relax/branches/disp_speed/lib/dispersion/m61.py?rev=23277&r1=23276&r2=23277&view=diff > ============================================================================== > --- branches/disp_speed/lib/dispersion/m61.py (original) > +++ branches/disp_speed/lib/dispersion/m61.py Tue May 20 22:29:52 2014 > @@ -64,10 +64,10 @@ > """ > > # Python module imports. > -from numpy import array, isfinite, sum > +from numpy import abs, array, isfinite, min, sum > > > -def r1rho_M61(r1rho_prime=None, phi_ex=None, kex=None, > spin_lock_fields2=None, back_calc=None, num_points=None): > +def r1rho_M61(r1rho_prime=None, phi_ex=None, kex=None, > spin_lock_fields2=None, num_points=None): > """Calculate the R2eff values for the M61 model. > > See the module docstring for details. > @@ -81,9 +81,7 @@ > @type kex: float > @keyword spin_lock_fields2: The R1rho spin-lock field strengths squared > (in rad^2.s^-2). > @type spin_lock_fields2: numpy rank-1 float array > - @keyword back_calc: The array for holding the back calculated > R1rho values. Each element corresponds to one of the spin-lock fields. > - @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 spin_lock_fields and back_calc arguments. > + @keyword num_points: The number of points on the dispersion > curve, equal to the length of the spin_lock_fields. > @type num_points: int > """ > > @@ -96,6 +94,13 @@ > # Denominator. > denom = kex2 + spin_lock_fields2 > > + # Catch math domain error of dividing with 0. > + # This is when denom=0. > + if min(abs(denom)) == 0: > + R1rho = array([1e100]*num_points) > + > + return R1rho > + > # R1rho calculation. > R1rho = r1rho_prime + numer / denom > > @@ -104,6 +109,4 @@ > if not isfinite(sum(R1rho)): > R1rho = array([1e100]*num_points) > > - # Parse back the value to update the back_calc class object. > - for i in range(num_points): > - back_calc[i] = R1rho[i] > + return R1rho > > 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=23277&r1=23276&r2=23277&view=diff > ============================================================================== > --- branches/disp_speed/target_functions/relax_disp.py (original) > +++ branches/disp_speed/target_functions/relax_disp.py Tue May 20 22:29:52 > 2014 > @@ -1163,7 +1163,7 @@ > phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2 > > # Back calculate the R2eff values. > - r1rho_M61(r1rho_prime=R20[r20_index], phi_ex=phi_ex_scaled, > kex=kex, spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][0], > back_calc=self.back_calc[0][si][mi][0], > num_points=self.num_disp_points[0][si][mi][0]) > + self.back_calc[0][si][mi][0] = > r1rho_M61(r1rho_prime=R20[r20_index], phi_ex=phi_ex_scaled, kex=kex, > spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][0], > num_points=self.num_disp_points[0][si][mi][0]) > > # 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[0][si][mi][0]): > > > _______________________________________________ > relax (http://www.nmr-relax.com) > > This is the relax-commits mailing list > [email protected] > > To unsubscribe from this list, get a password > reminder, or change your subscription options, > visit the list information page at > https://mail.gna.org/listinfo/relax-commits _______________________________________________ relax (http://www.nmr-relax.com) This is the relax-devel mailing list [email protected] To unsubscribe from this list, get a password reminder, or change your subscription options, visit the list information page at https://mail.gna.org/listinfo/relax-devel
