Author: tlinnet Date: Thu Jun 12 13:56:12 2014 New Revision: 23880 URL: http://svn.gna.org/viewcvs/relax?rev=23880&view=rev Log: Large increase in speed for model TSMFK01 by changing target functions to use multidimensional numpy arrays in calculation. This is done by restructuring data into multidimensional arrays of dimension [NE][NS][NM][NO][ND], which are number of spins, number of magnetic field strength, number of offsets, maximum number of dispersion point. The speed comes from using numpy ufunc operations. The new version is 2.4X as fast per spin calculation, and 54X as fast for clustered analysis. The different in timings for 3 spectrometer frequencies, calculated for 1 spin or 100 clustered spins with 1000 iterations are: ---- VERSION 3.2.2 ---- 1 spin: ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 0.262 0.262 <string>:1(<module>) 100 spin: ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 25.391 25.391 <string>:1(<module>) ---- New version --- 1 spin: ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 0.111 0.111 <string>:1(<module>) 100 spin: ncalls tottime percall cumtime percall filename:lineno(function) 1 0.000 0.000 0.468 0.468 <string>:1(<module>) Task #7807 (https://gna.org/task/index.php?7807): Speed-up of dispersion models for Clustered analysis. Modified: branches/disp_spin_speed/target_functions/relax_disp.py Modified: branches/disp_spin_speed/target_functions/relax_disp.py URL: http://svn.gna.org/viewcvs/relax/branches/disp_spin_speed/target_functions/relax_disp.py?rev=23880&r1=23879&r2=23880&view=diff ============================================================================== --- branches/disp_spin_speed/target_functions/relax_disp.py (original) +++ branches/disp_spin_speed/target_functions/relax_disp.py Thu Jun 12 13:56:12 2014 @@ -395,7 +395,7 @@ self.func = self.func_ns_mmq_3site_linear # Setup special numpy array structures, for higher dimensional computation. - if model in [MODEL_B14, MODEL_B14_FULL, MODEL_CR72, MODEL_CR72_FULL]: + if model in [MODEL_B14, MODEL_B14_FULL, MODEL_CR72, MODEL_CR72_FULL, MODEL_TSMFK01]: # Get the shape of back_calc structure. # If using just one field, or having the same number of dispersion points, the shape would extend to that number. # Shape has to be: [ei][si][mi][oi]. @@ -478,7 +478,7 @@ if self.missing[ei][si][mi][oi][di]: self.has_missing = True missing_a[ei][si][mi][oi][di] = 1.0 - if model in [MODEL_B14, MODEL_B14_FULL]: + if model in [MODEL_B14, MODEL_B14_FULL, MODEL_TSMFK01]: self.power_a[ei][si][mi][oi][di] = int(round(self.cpmg_frqs[ei][mi][0][di] * self.relax_times[ei][mi])) self.tau_cpmg_a[ei][si][mi][oi][di] = 0.25 / self.cpmg_frqs[ei][mi][0][di] @@ -1989,29 +1989,22 @@ dw = params[self.end_index[0]:self.end_index[1]] k_AB = params[self.end_index[1]] - # Initialise. - chi2_sum = 0.0 - - # Loop over the spins. - for si in range(self.num_spins): - # Loop over the spectrometer frequencies. - for mi in range(self.num_frq): - # The R20 index. - r20a_index = mi + si*self.num_frq - - # Convert dw from ppm to rad/s. - dw_frq = dw[si] * self.frqs[0][si][mi] - - # Back calculate the R2eff values. - r2eff_TSMFK01(r20a=R20A[r20a_index], dw=dw_frq, dw_orig=dw, k_AB=k_AB, tcp=self.tau_cpmg[0][mi], back_calc=self.back_calc[0][si][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]): - if self.missing[0][si][mi][0][di]: - self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] - - # Calculate and return the chi-squared value. - chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) + # Convert dw from ppm to rad/s. Use the out argument, to pass directly to structure. + multiply( multiply.outer( dw.reshape(self.NE, self.NS), self.nm_no_nd_struct ), self.frqs_struct, out=self.dw_struct ) + + # Reshape R20A and R20B to per experiment, spin and frequency. + self.r20a_struct[:] = multiply.outer( R20A.reshape(self.NE, self.NS, self.NM), self.no_nd_struct ) + + # Back calculate the R2eff values. + r2eff_TSMFK01(r20a=self.r20a_struct, dw=self.dw_struct, dw_orig=dw, k_AB=k_AB, tcp=self.tau_cpmg_a, back_calc=self.back_calc_a, num_points=self.num_disp_points_a) + + # Clean the data for all values, which is left over at the end of arrays. + self.back_calc_a = self.back_calc_a*self.disp_struct + + ## 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. + if self.has_missing: + # Replace with values. + self.back_calc_a[self.mask_replace_blank.mask] = self.values_a[self.mask_replace_blank.mask] # Return the total chi-squared value. - return chi2_sum + return chi2_rankN(self.values_a, self.back_calc_a, self.errors_a)