Interesting! If the magnetisation propagation can also be shifted out of the loops, then the loops will disappear and the numeric models will become far quicker. With a matrix exponential, maybe this is possible. I'll think about this for a while.
Regards, Edward On 19 June 2014 20:52, <[email protected]> wrote: > Author: tlinnet > Date: Thu Jun 19 20:52:57 2014 > New Revision: 24171 > > URL: http://svn.gna.org/viewcvs/relax?rev=24171&view=rev > Log: > Moved the costly calculation of the matrix exponential out of for loops. > > It was the numpy.eig and numpy.inv which was draining power. > > This speeds up model NS R1rho 2site, by a factor 4X: > BEFORE: > Single: > ncalls tottime percall cumtime percall filename:lineno(function) > 1 0.000 0.000 32.552 32.552 <string>:1(<module>) > 1 0.002 0.002 32.552 32.552 pf_nsr1rho2site:530(single) > Cluster: > ncalls tottime percall cumtime percall filename:lineno(function) > 1 0.000 0.000 33.307 33.307 <string>:1(<module>) > 1 0.008 0.008 33.307 33.307 pf_nsr1rho2site:554(cluster) > > AFTER: > Single: > ncalls tottime percall cumtime percall filename:lineno(function) > 1 0.000 0.000 8.286 8.286 <string>:1(<module>) > 1 0.002 0.002 8.286 8.286 pf_nsr1rho2site:530(single) > Cluster: > ncalls tottime percall cumtime percall filename:lineno(function) > 1 0.000 0.000 8.223 8.223 <string>:1(<module>) > 1 0.007 0.007 8.223 8.223 pf_nsr1rho2site:554(cluster) > > 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_r1rho_2site.py > > Modified: branches/disp_spin_speed/lib/dispersion/ns_r1rho_2site.py > URL: > http://svn.gna.org/viewcvs/relax/branches/disp_spin_speed/lib/dispersion/ns_r1rho_2site.py?rev=24171&r1=24170&r2=24171&view=diff > ============================================================================== > --- branches/disp_spin_speed/lib/dispersion/ns_r1rho_2site.py (original) > +++ branches/disp_spin_speed/lib/dispersion/ns_r1rho_2site.py Thu Jun 19 > 20:52:57 2014 > @@ -56,7 +56,7 @@ > # relax module imports. > from lib.dispersion.ns_matrices import rr1rho_3d, rr1rho_3d_rankN > from lib.float import isNaN > -from lib.linear_algebra.matrix_exponential import matrix_exponential > +from lib.linear_algebra.matrix_exponential import matrix_exponential, > matrix_exponential_rankN > > > def ns_r1rho_2site(M0=None, matrix=None, r1rho_prime=None, omega=None, > offset=None, r1=0.0, pA=None, dw=None, kex=None, spin_lock_fields=None, > relax_time=None, inv_relax_time=None, back_calc=None, num_points=None): > @@ -106,6 +106,9 @@ > # The matrix that contains all the contributions to the evolution, i.e. > relaxation, exchange and chemical shift evolution. > R_mat = rr1rho_3d_rankN(R1=r1, r1rho_prime=r1rho_prime, pA=pA, pB=pB, > dw=dw, omega=omega, offset=offset, w1=spin_lock_fields, k_AB=k_AB, k_BA=k_BA, > relax_time=relax_time) > > + # This matrix is a propagator that will evolve the magnetization with > the matrix R. > + Rexpo_mat = matrix_exponential_rankN(R_mat) > + > # Loop over spins. > for si in range(NS): > # Loop over the spectrometer frequencies. > @@ -135,19 +138,16 @@ > > # Loop over the time points, back calculating the R2eff > values. > for j in range(num_points_i): > - # The matrix that contains all the contributions to the > evolution, i.e. relaxation, exchange and chemical shift evolution. > - R_mat_i = R_mat[0, si, mi, oi, j] > - > # The following lines rotate the magnetization previous > to spin-lock into the weff frame. > theta = atan2(spin_lock_fields_i[j], dA) > M0[0] = sin(theta) # The A state initial X > magnetisation. > M0[2] = cos(theta) # The A state initial Z > magnetisation. > > # This matrix is a propagator that will evolve the > magnetization with the matrix R. > - Rexpo = matrix_exponential(R_mat_i) > + Rexpo_i = Rexpo_mat[0, si, mi, oi, j] > > # Magnetization evolution. > - MA = dot(M0, dot(Rexpo, M0)) > + MA = dot(M0, dot(Rexpo_i, M0)) > > # The next lines calculate the R1rho using a two-point > approximation, i.e. assuming that the decay is mono-exponential. > if MA <= 0.0 or isNaN(MA): > > > _______________________________________________ > 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
