Hi Ed. The problem is, that 2-4 lines of code in a lib file can trigger a math domain error.
To figure out what each function is "sad" about, one would need a little inspection per lib file, and the code would get full of checkings and handlings. That slows it down... It should be, that in 95 percent of the time, the minimise function is just happy and calculating. Only at boundaries values, the math domain errors can get triggered. Rather that 5% of the time, it calculates to much, that 95% is full of checking. Best Troels 2014-05-19 10:18 GMT+02:00 Edward d'Auvergne <[email protected]>: > Hi Troels, > > There is actually a better way that is just as fast. That is to catch > the value prior to the divide by zero, or what ever operation creates > the Inf values. Then you can fill back_calc with 1e100 and return. > That way the test will pass when running with --numpy-raise, and the > user will not see all the numpy warning messages. Having the check > earlier might even be a tiny bit faster. > > Regards, > > Edward > > > > > On 19 May 2014 00:51, <[email protected]> wrote: >> Author: tlinnet >> Date: Mon May 19 00:51:20 2014 >> New Revision: 23220 >> >> URL: http://svn.gna.org/viewcvs/relax?rev=23220&view=rev >> Log: >> Huge speed-up of model B14. >> >> task #7793: (https://gna.org/task/?7793) Speed-up of dispersion models. >> >> Time test for systemtests: >> >> test_baldwin_synthetic >> 2.626s -> 1.990s >> >> test_baldwin_synthetic_full >> 18.326s -> 13.742s >> >> 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/b14.py >> >> Modified: branches/disp_speed/lib/dispersion/b14.py >> URL: >> http://svn.gna.org/viewcvs/relax/branches/disp_speed/lib/dispersion/b14.py?rev=23220&r1=23219&r2=23220&view=diff >> ============================================================================== >> --- branches/disp_speed/lib/dispersion/b14.py (original) >> +++ branches/disp_speed/lib/dispersion/b14.py Mon May 19 00:51:20 2014 >> @@ -110,8 +110,7 @@ >> """ >> >> # Python module imports. >> -import numpy >> -from numpy import arccosh, array, cos, cosh, in1d, log, nonzero, sin, sinh, >> sqrt, power >> +from numpy import arccosh, array, cos, cosh, isfinite, log, power, sin, >> sinh, sqrt, sum >> >> # Repetitive calculations (to speed up calculations). >> g_fact = 1/sqrt(2) >> @@ -216,70 +215,31 @@ >> # Real. The v_1c in paper. >> v1c = F0 * cosh(E0) - F2 * cos(E2) >> >> - # Catch devision with zero in y, when v3 = 0. v3 is 0, when v1c = 1. >> - # If no 1.0, perform normally. >> - if not in1d(1.0, v1c): >> - # Exact result for v2v3. >> - v3 = sqrt(v1c**2 - 1.) >> - >> - y = power( (v1c - v3) / (v1c + v3), ncyc) >> - >> - Tog = 0.5 * (1. + y) + (1. - y) * v5 / (2. * v3 * N ) >> - >> - # Find where Tog has negative values. >> - neg_index = nonzero(Tog.real < 0.0)[0] >> - >> - # Do normal calculation >> - if len(neg_index) == 0: >> - ## -1/Trel * log(LpreDyn). >> - # Rpre = (r20a + r20b + kex) / 2.0 >> - >> - ## Carver and Richards (1972) >> - # R2eff_CR72 = Rpre - inv_tcpmg * ncyc * arccosh(v1c.real) >> - >> - ## Baldwin final. >> - # Estimate R2eff. relax_time = Trel = 1/inv_tcpmg. >> - # R2eff = R2eff_CR72 - inv_tcpmg * log(Tog.real) >> - >> - # Fastest calculation. >> - R2eff = (r20a + r20b + kex) / 2.0 - inv_tcpmg * ( ncyc * >> arccosh(v1c.real) + log(Tog.real) ) >> - >> - # Loop over the time points, back calculating the R2eff values. >> - for i in range(num_points): >> - >> - # Put values back. >> - back_calc[i] = R2eff[i] >> - >> - else: >> - # Loop over each point. >> - for i in range(num_points): >> - >> - # Return large value >> - if i in neg_index: >> - back_calc[i] = 1e100 >> - >> - else: >> - v3 = sqrt(v1c[i]**2 - 1.) >> - y = power( (v1c[i] - v3) / (v1c[i] + v3), ncyc[i]) >> - Tog = 0.5 * (1. + y) + (1. - y) * v5[i] / (2. * v3 * N ) >> - R2eff = (r20a + r20b + kex) / 2.0 - inv_tcpmg * ( >> ncyc[i] * arccosh(v1c[i].real) + log(Tog.real) ) >> - back_calc[i] = R2eff >> - >> - # This section is for catching math domain errors. >> - else: >> - # Find index where >> - one_indexes = nonzero(v1c == 1.0)[0] >> - >> - # Loop over each point. >> - for i in range(num_points): >> - >> - # Return large value >> - if i in one_indexes: >> - back_calc[i] = 1e100 >> - >> - else: >> - v3 = sqrt(v1c[i]**2 - 1.) >> - y = power( (v1c[i] - v3) / (v1c[i] + v3), ncyc[i]) >> - Tog = 0.5 * (1. + y) + (1. - y) * v5[i] / (2. * v3 * N ) >> - R2eff = (r20a + r20b + kex) / 2.0 - inv_tcpmg * ( ncyc[i] >> * arccosh(v1c[i].real) + log(Tog.real) ) >> - back_calc[i] = R2eff >> + # Exact result for v2v3. >> + v3 = sqrt(v1c**2 - 1.) >> + >> + y = power( (v1c - v3) / (v1c + v3), ncyc) >> + >> + Tog = 0.5 * (1. + y) + (1. - y) * v5 / (2. * v3 * N ) >> + >> + ## -1/Trel * log(LpreDyn). >> + # Rpre = (r20a + r20b + kex) / 2.0 >> + >> + ## Carver and Richards (1972) >> + # R2eff_CR72 = Rpre - inv_tcpmg * ncyc * arccosh(v1c.real) >> + >> + ## Baldwin final. >> + # Estimate R2eff. relax_time = Trel = 1/inv_tcpmg. >> + # R2eff = R2eff_CR72 - inv_tcpmg * log(Tog.real) >> + >> + # Fastest calculation. >> + R2eff = (r20a + r20b + kex) / 2.0 - inv_tcpmg * ( ncyc * >> arccosh(v1c.real) + log(Tog.real) ) >> + >> + # 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): >> + back_calc[i] = R2eff[i] >> >> >> _______________________________________________ >> 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 _______________________________________________ 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
