Re: r24191 - /branches/disp_spin_speed/lib/dispersion/ns_mmq_2site.py -- June 23, 2014

Great work!  You can take this much, much further though and move the
M1_M2, M2_M1, M1_M2_star, M2_M1_star, M1_M2_M2_M1_star, and
M2_M1_M1_M2_star dot products out of the loop.  For this model, that
would be a huge win!  Then for this model, just like the other numeric
models, the only thing keeping the slow looping alive is the
magentisation evolution - which might itself be shifted out of the
loop in the future.

Regards,

Edward


On 20 June 2014 09:38,  <tlinnet@xxxxxxxxxxxxx> wrote:

Author: tlinnet
Date: Fri Jun 20 09:38:30 2014
New Revision: 24191

URL: http://svn.gna.org/viewcvs/relax?rev=24191&view=rev
Log:
Moved the costly calculation of matrix_exponential of M1 and M2 out of for 
loop, in model ns_mmq_2site_mq.

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_mmq_2site.py

Modified: branches/disp_spin_speed/lib/dispersion/ns_mmq_2site.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/disp_spin_speed/lib/dispersion/ns_mmq_2site.py?rev=24191&r1=24190&r2=24191&view=diff
==============================================================================
--- branches/disp_spin_speed/lib/dispersion/ns_mmq_2site.py     (original)
+++ branches/disp_spin_speed/lib/dispersion/ns_mmq_2site.py     Fri Jun 20 
09:38:30 2014
@@ -55,7 +55,7 @@

 # relax module imports.
 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
 from lib.linear_algebra.matrix_power import square_matrix_power


@@ -212,8 +212,22 @@
     NS, NM, NO = num_points.shape

     # Populate the m1 and m2 matrices (only once per function call for 
speed).
+    # D+ matrix component.
     m1_mat = populate_matrix_rankN(R20A=R20A, R20B=R20B, dw=-dw - dwH, 
k_AB=k_AB, k_BA=k_BA, tcp=tcp)
+    # Z- matrix component.
     m2_mat = populate_matrix_rankN(R20A=R20A, R20B=R20B, dw=dw - dwH, 
k_AB=k_AB, k_BA=k_BA, tcp=tcp)
+
+    # The M1 and M2 matrices.
+    # Equivalent to D+.
+    M1_mat = matrix_exponential_rankN(m1_mat)
+    # Equivalent to Z-.
+    M2_mat = matrix_exponential_rankN(m2_mat)
+
+    # The complex conjugates M1* and M2*
+    # Equivalent to D+*.
+    M1_mat_star = conj(M1_mat)
+    # Equivalent to Z-*.
+    M2_mat_star = conj(M2_mat)

     # Loop over spins.
     for si in range(NS):
@@ -221,51 +235,29 @@
         for mi in range(NM):
             # Loop over offsets:
             for oi in range(NO):
-
-                r20a_i = R20A[si, mi, oi, 0]
-                r20b_i = R20B[si, mi, oi, 0]
-                dw_i = dw[si, mi, oi, 0]
-                dwH_i = dwH[si, mi, oi, 0]
                 num_points_i = num_points[si, mi, oi]
-
-                # Populate the m1 and m2 matrices (only once per function 
call for speed).
-                populate_matrix(matrix=m1, R20A=r20a_i, R20B=r20b_i, 
dw=-dw_i - dwH_i, k_AB=k_AB, k_BA=k_BA)     # D+ matrix component.
-                populate_matrix(matrix=m2, R20A=r20a_i, R20B=r20b_i, 
dw=dw_i - dwH_i, k_AB=k_AB, k_BA=k_BA)    # Z- matrix component.

                 # Loop over the time points, back calculating the R2eff 
values.
                 for i in range(num_points_i):
-                    m1_mat_i = m1_mat[si, mi, oi, i]
-                    m2_mat_i = m2_mat[si, mi, oi, i]
-
-                    diff_m1 = abs(sum(m1*tcp[si, mi, oi, i] - m1_mat_i))
-                    if diff_m1 > 1.0e-06:
-                        print diff_m1
-                        print m1*tcp[si, mi, oi, i]
-                        print m1_mat_i
-                        print asd
-
-                    diff_m2 = abs(sum(m2*tcp[si, mi, oi, i] - m2_mat_i))
-                    if diff_m2 > 1.0e-06:
-                        print diff_m2
-                        print m2*tcp[si, mi, oi, i]
-                        print m2_mat_i
-                        print asd
-
                     # The M1 and M2 matrices.
-                    M1 = matrix_exponential(m1_mat_i)    # Equivalent to 
D+.
-                    M2 = matrix_exponential(m2_mat_i)    # Equivalent to 
Z-.
+                    # Equivalent to D+.
+                    M1_i = M1_mat[0, si, mi, oi, i]
+                    # Equivalent to Z-.
+                    M2_i = M1_mat[0, si, mi, oi, i]

                     # The complex conjugates M1* and M2*
-                    M1_star = conj(M1)    # Equivalent to D+*.
-                    M2_star = conj(M2)    # Equivalent to Z-*.
+                    # Equivalent to D+*.
+                    M1_star_i = M1_mat_star[0, si, mi, oi, i]
+                    # Equivalent to Z-*.
+                    M2_star_i = M1_mat_star[0, si, mi, oi, i]

                     # Repetitive dot products (minimised for speed).
-                    M1_M2 = dot(M1, M2)
-                    M2_M1 = dot(M2, M1)
+                    M1_M2 = dot(M1_i, M2_i)
+                    M2_M1 = dot(M2_i, M1_i)
                     M1_M2_M2_M1 = dot(M1_M2, M2_M1)
                     M2_M1_M1_M2 = dot(M2_M1, M1_M2)
-                    M1_M2_star = dot(M1_star, M2_star)
-                    M2_M1_star = dot(M2_star, M1_star)
+                    M1_M2_star = dot(M1_star_i, M2_star_i)
+                    M2_M1_star = dot(M2_star_i, M1_star_i)
                     M1_M2_M2_M1_star = dot(M1_M2_star, M2_M1_star)
                     M2_M1_M1_M2_star = dot(M2_M1_star, M1_M2_star)



_______________________________________________
relax (http://www.nmr-relax.com)

This is the relax-commits mailing list
relax-commits@xxxxxxx

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

Re: r24191 - /branches/disp_spin_speed/lib/dispersion/ns_mmq_2site.py

Header

Content

Related Messages