Wow, this is impressive! And if you come up with a multi-rank replacement for the square_matrix_power() function, you'll be able to finally eliminate all the looping and make this numeric model hugely faster. Although we can't test it directly, I would not be surprised if this would then be by far the fastest implementation in the field. Though we don't have all of the tricks implemented for the numeric models - off-resonance effects, the correction factor Andy talked about, etc. - this will still be huge.
Regards, Edward On 24 June 2014 14:58, <[email protected]> wrote: > Author: tlinnet > Date: Tue Jun 24 14:58:07 2014 > New Revision: 24280 > > URL: http://svn.gna.org/viewcvs/relax?rev=24280&view=rev > Log: > Speeded up model NS CPMG 2site star, by moving the forming of the propagator > matrix out of the for loops, and preform it. > > Task #7807 (https://gna.org/task/index.php?7807): Speed-up of dispersion > models for Clustered analysis. > > Modified: > branches/disp_spin_speed/lib/dispersion/ns_cpmg_2site_star.py > > Modified: branches/disp_spin_speed/lib/dispersion/ns_cpmg_2site_star.py > URL: > http://svn.gna.org/viewcvs/relax/branches/disp_spin_speed/lib/dispersion/ns_cpmg_2site_star.py?rev=24280&r1=24279&r2=24280&view=diff > ============================================================================== > --- branches/disp_spin_speed/lib/dispersion/ns_cpmg_2site_star.py > (original) > +++ branches/disp_spin_speed/lib/dispersion/ns_cpmg_2site_star.py Tue > Jun 24 14:58:07 2014 > @@ -57,7 +57,7 @@ > """ > > # Python module imports. > -from numpy import add, array, conj, dot, fabs, float64, isfinite, log, min, > multiply, sum > +from numpy import add, array, conj, dot, einsum, fabs, float64, isfinite, > log, min, multiply, sum > from numpy.ma import fix_invalid, masked_where > > # relax module imports. > @@ -221,8 +221,15 @@ > # The matrix R that contains all the contributions to the evolution, > i.e. relaxation, exchange and chemical shift evolution. > R_mat, cR2_mat, Rr_mat, Rex_mat, RCS_mat = rcpmg_star_rankN(R2A=r20a, > R2B=r20b, dw=dw, k_AB=k_AB, k_BA=k_BA, tcp=tcp) > > + # The the essential evolution matrix. > + # This matrix is a propagator that will evolve the magnetization with > the matrix R for a delay tcp. > eR_mat = matrix_exponential_rank_NE_NS_NM_NO_ND_x_x(R_mat) > ecR2_mat = matrix_exponential_rank_NE_NS_NM_NO_ND_x_x(cR2_mat) > + > + # Preform the matrix. > + # This is the propagator for an element of [delay tcp; 180 deg pulse; 2 > times delay tcp; 180 deg pulse; delay tau], i.e. for 2 times tau-180-tau. > + prop_2_mat = evolution_matrix_mat = einsum('...ij,...jk', eR_mat, > ecR2_mat) > + prop_2_mat = evolution_matrix_mat = einsum('...ij,...jk', prop_2_mat, > eR_mat) > > # Loop over the spins > for si in range(NS): > @@ -236,16 +243,11 @@ > # Extract the values from the higher dimensional arrays. > power_si_mi_di = int(power[0, si, mi, 0, di]) > > - # This matrix is a propagator that will evolve the > magnetization with the matrix R for a delay tcp. > - eR_tcp = eR_mat[0, si, mi, 0, di] > - ecR2_tcp = ecR2_mat[0, si, mi, 0, di] > - > # This is the propagator for an element of [delay tcp; 180 > deg pulse; 2 times delay tcp; 180 deg pulse; delay tau], i.e. for 2 times > tau-180-tau. > - prop_2 = dot(eR_tcp, ecR2_tcp) > - prop_2 = dot(prop_2, eR_tcp) > + prop_2_i = prop_2_mat[0, si, mi, 0, di] > > # Now create the total propagator that will evolve the > magnetization under the CPMG train, i.e. it applies the above > tau-180-tau-tau-180-tau so many times as required for the CPMG frequency > under consideration. > - prop_total = square_matrix_power(prop_2, power_si_mi_di) > + prop_total = square_matrix_power(prop_2_i, power_si_mi_di) > > # Now we apply the above propagator to the initial > magnetization vector - resulting in the magnetization that remains after the > full CPMG pulse train. It is called M of t (t is the time after the CPMG > train). > Moft = dot(prop_total, M0) > > > _______________________________________________ > 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
