Hi ed. All the: num_points=self.num_disp_points_a
can also be killed. They are not used. The disp_struct is actually this structure, in higher dimensions? Best Troels 2014-06-13 9:02 GMT+02:00 Edward d'Auvergne <[email protected]>: > Hi Troels, > > Thanks to the lightning quick infrastructure you are putting into > place, we can also simplify the target_functions.relax_disp to > lib.dispersion API. You will notice a lot of code like in this TP02 > model: > > + # Once off parameter conversions. > + pB = 1.0 - pA > > This was originally in lib.dispersion (well at least for some of the > models), but I shifted it into the func_*() target functions to speed > the code up, as then the calculation would happen only once rather > than once for each iteration of that massive looping beast you have > killed. > > So now that the lib.dispersion functions are only called once per > target function call, all of these 'Once off parameter conversions' > can be shifted back into the lib.dispersion functions. Then the > number of arguments for these functions will drop significantly, as > the {pB, k_AB, k_BA} parameters will no longer need to be passed in. > This will be far more significant for the 3-site models where there > are many, many parameter conversions. This will have the added > benefit of simplifying the use of the lib.dispersion modules outside > of the dispersion target functions, for example with the unit testing. > > Cheers, > > Edward > > > > > On 13 June 2014 07:21, <[email protected]> wrote: > > Author: tlinnet > > Date: Fri Jun 13 07:21:02 2014 > > New Revision: 23901 > > > > URL: http://svn.gna.org/viewcvs/relax?rev=23901&view=rev > > Log: > > Replaced the loop structure in target function of TAP03 with numpy > arrays. > > > > This makes the model faster. > > > > 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=23901&r1=23900&r2=23901&view=diff > > > ============================================================================== > > --- branches/disp_spin_speed/target_functions/relax_disp.py > (original) > > +++ branches/disp_spin_speed/target_functions/relax_disp.py Fri Jun > 13 07:21:02 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, MODEL_DPL94, MODEL_TSMFK01]: > > + if model in [MODEL_B14, MODEL_B14_FULL, MODEL_CR72, > MODEL_CR72_FULL, MODEL_DPL94, MODEL_TAP03, 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]. > > @@ -435,10 +435,12 @@ > > self.power_a = ones(self.numpy_array_shape, int16) > > > > # For R1rho data. > > - if model in [MODEL_DPL94]: > > + if model in [MODEL_DPL94, MODEL_TAP03]: > > self.tilt_angles_a = deepcopy(zeros_a) > > self.spin_lock_omega1_squared_a = deepcopy(zeros_a) > > + self.spin_lock_omega1_a = deepcopy(zeros_a) > > self.phi_ex_struct = deepcopy(zeros_a) > > + self.offset_a = deepcopy(zeros_a) > > > > self.frqs_a = deepcopy(zeros_a) > > self.disp_struct = deepcopy(zeros_a) > > @@ -447,6 +449,7 @@ > > # Create special numpy structures. > > # Structure of dw. The full and the outer dimensions > structures. > > self.dw_struct = deepcopy(zeros_a) > > + self.no_nd_struct = ones([self.NO, self.ND], float64) > > self.nm_no_nd_struct = ones([self.NM, self.NO, self.ND], > float64) > > > > # Structure of r20a and r20b. The full and outer dimensions > structures. > > @@ -459,10 +462,11 @@ > > # Expand relax times. > > self.inv_relax_times_a = 1.0 / multiply.outer( > tile(self.relax_times[:,None],(1, 1, self.NS)).reshape(self.NE, self.NS, > self.NM), self.no_nd_struct ) > > > > - if model in [MODEL_DPL94]: > > + if model in [MODEL_DPL94, MODEL_TAP03]: > > self.r1_a = multiply.outer( self.r1.reshape(self.NE, > self.NS, self.NM), self.no_nd_struct ) > > - > > - # Extract the frequencies to numpy array. > > + self.chemical_shifts_a = multiply.outer( > self.chemical_shifts, self.no_nd_struct ) > > + > > + # Extract the frequencies to numpy array. > > self.frqs_a = multiply.outer( > asarray(self.frqs).reshape(self.NE, self.NS, self.NM), self.no_nd_struct ) > > > > # Loop over the experiment types. > > @@ -476,7 +480,7 @@ > > # Extract number of dispersion points. > > num_disp_points = > self.num_disp_points[ei][si][mi][oi] > > > > - if model not in [MODEL_DPL94]: > > + if model not in [MODEL_DPL94, MODEL_TAP03]: > > # Extract cpmg_frqs and num_disp_points > from lists. > > > self.cpmg_frqs_a[ei][si][mi][oi][:num_disp_points] = > self.cpmg_frqs[ei][mi][oi] > > > self.num_disp_points_a[ei][si][mi][oi][:num_disp_points] = > self.num_disp_points[ei][si][mi][oi] > > @@ -498,12 +502,14 @@ > > 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] > > # For R1rho data. > > - if model in [MODEL_DPL94]: > > + if model in [MODEL_DPL94, MODEL_TAP03]: > > > self.disp_struct[ei][si][mi][oi][di] = 1.0 > > > > # Extract the frequencies to numpy > array. > > > self.tilt_angles_a[ei][si][mi][oi][di] = > self.tilt_angles[ei][si][mi][oi][di] > > > self.spin_lock_omega1_squared_a[ei][si][mi][oi][di] = > self.spin_lock_omega1_squared[ei][mi][oi][di] > > + > self.spin_lock_omega1_a[ei][si][mi][oi][di] = > self.spin_lock_omega1[ei][mi][oi][di] > > + self.offset_a[ei][si][mi][oi] = > self.offset[ei][si][mi][oi] > > > > if spin_lock_nu1 != None and > len(spin_lock_nu1[ei][mi][oi]): > > > self.num_disp_points_a[ei][si][mi][oi][di] = num_disp_points > > @@ -1908,6 +1914,49 @@ > > # Once off parameter conversions. > > pB = 1.0 - pA > > > > + # 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 R20 to per experiment, spin and frequency. > > + self.r20_struct[:] = multiply.outer( R20.reshape(self.NE, > self.NS, self.NM), self.no_nd_struct ) > > + > > + # Back calculate the R1rho values. > > + r1rho_TAP03(r1rho_prime=self.r20_struct, > omega=self.chemical_shifts_a, offset=self.offset_a, pA=pA, pB=pB, > dw=self.dw_struct, kex=kex, R1=self.r1_a, > spin_lock_fields=self.spin_lock_omega1_a, > spin_lock_fields2=self.spin_lock_omega1_squared_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_rankN(self.values_a, self.back_calc_a, > self.errors_a) > > + > > + > > + def func_TP02(self, params): > > + """Target function for the Trott and Palmer (2002) R1rho > off-resonance 2-site model. > > + > > + @param params: The vector of parameter values. > > + @type params: numpy rank-1 float array > > + @return: The chi-squared value. > > + @rtype: float > > + """ > > + > > + # Scaling. > > + if self.scaling_flag: > > + params = dot(params, self.scaling_matrix) > > + > > + # Unpack the parameter values. > > + R20 = params[:self.end_index[0]] > > + dw = params[self.end_index[0]:self.end_index[1]] > > + pA = params[self.end_index[1]] > > + kex = params[self.end_index[1]+1] > > + > > + # Once off parameter conversions. > > + pB = 1.0 - pA > > + > > # Initialise. > > chi2_sum = 0.0 > > > > @@ -1924,7 +1973,7 @@ > > # Loop over the offsets. > > for oi in range(self.num_offsets[0][si][mi]): > > # Back calculate the R1rho values. > > - r1rho_TAP03(r1rho_prime=R20[r20_index], > omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], > pA=pA, pB=pB, dw=dw_frq, kex=kex, R1=self.r1[si, mi], > spin_lock_fields=self.spin_lock_omega1[0][mi][oi], > spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], > back_calc=self.back_calc[0][si][mi][oi], > num_points=self.num_disp_points[0][si][mi][oi]) > > + r1rho_TP02(r1rho_prime=R20[r20_index], > omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], > pA=pA, pB=pB, dw=dw_frq, kex=kex, R1=self.r1[si, mi], > spin_lock_fields=self.spin_lock_omega1[0][mi][oi], > spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], > back_calc=self.back_calc[0][si][mi][oi], > num_points=self.num_disp_points[0][si][mi][oi]) > > > > # 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][oi]): > > @@ -1938,58 +1987,6 @@ > > return chi2_sum > > > > > > - def func_TP02(self, params): > > - """Target function for the Trott and Palmer (2002) R1rho > off-resonance 2-site model. > > - > > - @param params: The vector of parameter values. > > - @type params: numpy rank-1 float array > > - @return: The chi-squared value. > > - @rtype: float > > - """ > > - > > - # Scaling. > > - if self.scaling_flag: > > - params = dot(params, self.scaling_matrix) > > - > > - # Unpack the parameter values. > > - R20 = params[:self.end_index[0]] > > - dw = params[self.end_index[0]:self.end_index[1]] > > - pA = params[self.end_index[1]] > > - kex = params[self.end_index[1]+1] > > - > > - # Once off parameter conversions. > > - pB = 1.0 - pA > > - > > - # 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. > > - r20_index = mi + si*self.num_frq > > - > > - # Convert dw from ppm to rad/s. > > - dw_frq = dw[si] * self.frqs[0][si][mi] > > - > > - # Loop over the offsets. > > - for oi in range(self.num_offsets[0][si][mi]): > > - # Back calculate the R1rho values. > > - r1rho_TP02(r1rho_prime=R20[r20_index], > omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], > pA=pA, pB=pB, dw=dw_frq, kex=kex, R1=self.r1[si, mi], > spin_lock_fields=self.spin_lock_omega1[0][mi][oi], > spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], > back_calc=self.back_calc[0][si][mi][oi], > num_points=self.num_disp_points[0][si][mi][oi]) > > - > > - # 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][oi]): > > - if self.missing[0][si][mi][oi][di]: > > - self.back_calc[0][si][mi][oi][di] = > self.values[0][si][mi][oi][di] > > - > > - # Calculate and return the chi-squared value. > > - chi2_sum += chi2(self.values[0][si][mi][oi], > self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) > > - > > - # Return the total chi-squared value. > > - return chi2_sum > > - > > - > > def func_TSMFK01(self, params): > > """Target function for the the Tollinger et al. (2001) 2-site > very-slow exchange model, range of microsecond to second time scale. > > > > > > > > _______________________________________________ > > 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
