r24286 - in /branches/disp_spin_speed: lib/dispersion/ns_r1rho_2site.py target_functions/relax_disp.py -- June 24, 2014 - 16:35

Author: tlinnet
Date: Tue Jun 24 16:35:01 2014
New Revision: 24286

URL: http://svn.gna.org/viewcvs/relax?rev=24286&view=rev
Log:
Speeded up ns r1rho 2site, by preforming the evolution matrices, and the M0 
matrix in the init part of the target function.

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
    branches/disp_spin_speed/target_functions/relax_disp.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=24286&r1=24285&r2=24286&view=diff
==============================================================================
--- branches/disp_spin_speed/lib/dispersion/ns_r1rho_2site.py   (original)
+++ branches/disp_spin_speed/lib/dispersion/ns_r1rho_2site.py   Tue Jun 24 
16:35:01 2014
@@ -50,8 +50,7 @@
 """
 
 # Python module imports.
-from math import atan2
-from numpy import array, cos, dot, float64, log, multiply, sin, sum
+from numpy import array, cos, dot, einsum, float64, log, multiply, sin, sum
 
 # relax module imports.
 from lib.float import isNaN
@@ -236,6 +235,9 @@
     # This matrix is a propagator that will evolve the magnetization with 
the matrix R.
     Rexpo_mat = matrix_exponential_rank_NE_NS_NM_NO_ND_x_x(R_mat)
 
+    # Magnetization evolution.
+    Rexpo_M0_mat = einsum('...ij,...jk', Rexpo_mat, M0)
+
     # Loop over spins.
     for si in range(NS):
         # Loop over the spectrometer frequencies.
@@ -245,23 +247,19 @@
                 # Extract number of points.
                 num_points_i = num_points[0, si, mi, oi]
 
-                # Repetitive calculations (to speed up calculations).
-                # Offset of spin-lock from A.
-                dA = omega[0, si, mi, oi, 0] - offset[0, si, mi, oi, 0]
-
                 # Loop over the time points, back calculating the R2eff 
values.
                 for j in range(num_points_i):
-                    # The following lines rotate the magnetization previous 
to spin-lock into the weff frame.
-                    theta = atan2(spin_lock_fields[0, si, mi, oi, j], dA)
-                    M0[0] = sin(theta)    # The A state initial X 
magnetisation.
-                    M0[2] = cos(theta)    # The A state initial Z 
magnetisation.
+                    # Extract the preformed matrix that rotate the 
magnetization previous to spin-lock into the weff frame.
+                    M0_i= M0[0, si, mi, oi, j, :, 0]
 
                     # This matrix is a propagator that will evolve the 
magnetization with the matrix R.
                     Rexpo_i = Rexpo_mat[0, si, mi, oi, j]
 
+                    # Extract from the pre-formed Magnetization evolution 
matrix.
+                    Rexpo_M0_mat_i = Rexpo_M0_mat[0, si, mi, oi, j, :, 0]
+
                     # Magnetization evolution.
-                    MA = dot(Rexpo_i, M0)
-                    MA = dot(M0, MA)
+                    MA = dot(M0_i, Rexpo_M0_mat_i)
 
                     # 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):

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=24286&r1=24285&r2=24286&view=diff
==============================================================================
--- branches/disp_spin_speed/target_functions/relax_disp.py     (original)
+++ branches/disp_spin_speed/target_functions/relax_disp.py     Tue Jun 24 
16:35:01 2014
@@ -27,7 +27,7 @@
 # Python module imports.
 from copy import deepcopy
 from math import pi
-from numpy import add, array, asarray, complex64, dot, float64, int16, max, 
multiply, ones, sqrt, sum, tile, zeros
+from numpy import add, array, arctan2, asarray, cos, complex64, dot, 
float64, int16, max, multiply, ones, sin, sqrt, sum, tile, zeros
 from numpy.ma import masked_equal
 
 # relax module imports.
@@ -404,7 +404,17 @@
             self.M0 = zeros(7, float64)
             self.M0[0] = 0.5
         if model in [MODEL_NS_R1RHO_2SITE]:
-            self.M0 = zeros(6, float64)
+            # Offset of spin-lock from A.
+            da_mat = self.chemical_shifts - self.offset
+            # The following lines rotate the magnetization previous to 
spin-lock into the weff frame.
+            theta_mat = arctan2(self.spin_lock_omega1, da_mat)
+            M0_0 = zeros([6, 1], float64)
+            M0_0[0, 0] = 1
+            M0_sin = multiply.outer( sin(theta_mat), M0_0 )
+            M0_2 = zeros([6, 1], float64)
+            M0_2[2, 0] = 1
+            M0_cos = multiply.outer( cos(theta_mat), M0_2 )
+            self.M0 = M0_sin + M0_cos
         if model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]:
             self.M0 = zeros(9, float64)