1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26 """The numerical fit of 2-site Bloch-McConnell equations for CPMG-type experiments, the U{NS CPMG 2-site 3D<http://wiki.nmr-relax.com/NS_CPMG_2-site_3D>} and U{NS CPMG 2-site 3D full<http://wiki.nmr-relax.com/NS_CPMG_2-site_3D_full>} models.
27
28 Description
29 ===========
30
31 The function uses an explicit matrix that contains relaxation, exchange and chemical shift terms. It does the 180deg pulses in the CPMG train. The approach of Bloch-McConnell can be found in chapter 3.1 of Palmer, A. G. Chem Rev 2004, 104, 3623-3640. This function was written, initially in MATLAB, in 2010.
32
33
34 Code origin
35 ===========
36
37 This is the model of the numerical solution for the 2-site Bloch-McConnell equations. It originates as optimization function number 1 from the fitting_main_kex.py script from Mathilde Lescanne, Paul Schanda, and Dominique Marion (see U{http://thread.gmane.org/gmane.science.nmr.relax.devel/4138}, U{https://web.archive.org/web/https://gna.org/task/?7712#comment2} and U{https://web.archive.org/web/https://gna.org/support/download.php?file_id=18262}).
38
39
40 Links
41 =====
42
43 More information on the NS CPMG 2-site 3D model can be found in the:
44
45 - U{relax wiki<http://wiki.nmr-relax.com/NS_CPMG_2-site_3D>},
46 - U{relax manual<http://www.nmr-relax.com/manual/reduced_NS_2_site_3D_CPMG_model.html>},
47 - U{relaxation dispersion page of the relax website<http://www.nmr-relax.com/analyses/relaxation_dispersion.html#NS_CPMG_2-site_3D>}.
48
49 More information on the NS CPMG 2-site 3D full model can be found in the:
50
51 - U{relax wiki<http://wiki.nmr-relax.com/NS_CPMG_2-site_3D_full>},
52 - U{relax manual<http://www.nmr-relax.com/manual/full_NS_2_site_3D_CPMG_model.html>},
53 - U{relaxation dispersion page of the relax website<http://www.nmr-relax.com/analyses/relaxation_dispersion.html#NS_CPMG_2-site_3D_full>}.
54 """
55
56
57 from math import log
58 from numpy import dot
59
60
61 from lib.dispersion.ns_matrices import rcpmg_3d
62 from lib.float import isNaN
63 from lib.linear_algebra.matrix_exponential import matrix_exponential
64
65
66 -def r2eff_ns_cpmg_2site_3D(r180x=None, M0=None, r10a=0.0, r10b=0.0, r20a=None, r20b=None, pA=None, pB=None, dw=None, k_AB=None, k_BA=None, inv_tcpmg=None, tcp=None, back_calc=None, num_points=None, power=None):
67 """The 2-site numerical solution to the Bloch-McConnell equation.
68
69 This function calculates and stores the R2eff values.
70
71
72 @keyword r180x: The X-axis pi-pulse propagator.
73 @type r180x: numpy float64, rank-2, 7D array
74 @keyword M0: This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations.
75 @type M0: numpy float64, rank-1, 7D array
76 @keyword r10a: The R1 value for state A.
77 @type r10a: float
78 @keyword r10b: The R1 value for state B.
79 @type r10b: float
80 @keyword r20a: The R2 value for state A in the absence of exchange.
81 @type r20a: float
82 @keyword r20b: The R2 value for state B in the absence of exchange.
83 @type r20b: float
84 @keyword pA: The population of state A.
85 @type pA: float
86 @keyword pB: The population of state B.
87 @type pB: float
88 @keyword dw: The chemical exchange difference between states A and B in rad/s.
89 @type dw: float
90 @keyword k_AB: The rate of exchange from site A to B (rad/s).
91 @type k_AB: float
92 @keyword k_BA: The rate of exchange from site B to A (rad/s).
93 @type k_BA: float
94 @keyword inv_tcpmg: The inverse of the total duration of the CPMG element (in inverse seconds).
95 @type inv_tcpmg: float
96 @keyword tcp: The tau_CPMG times (1 / 4.nu1).
97 @type tcp: numpy rank-1 float array
98 @keyword back_calc: The array for holding the back calculated R2eff values. Each element corresponds to one of the CPMG nu1 frequencies.
99 @type back_calc: numpy rank-1 float array
100 @keyword num_points: The number of points on the dispersion curve, equal to the length of the tcp and back_calc arguments.
101 @type num_points: int
102 @keyword power: The matrix exponential power array.
103 @type power: numpy int16, rank-1 array
104 """
105
106
107 if dw == 0.0 or pA == 1.0 or k_AB == 0.0:
108 for i in range(num_points):
109 back_calc[i] = r20a
110 return
111
112
113 R = rcpmg_3d(R1A=r10a, R1B=r10b, R2A=r20a, R2B=r20b, pA=pA, pB=pB, dw=dw, k_AB=k_AB, k_BA=k_BA)
114
115
116 for i in range(num_points):
117
118 Mint = M0
119
120
121 Rexpo = matrix_exponential(R*tcp[i])
122
123
124 for j in range(2*power[i]):
125 Mint = dot(Rexpo, Mint)
126 Mint = dot(r180x, Mint)
127 Mint = dot(Rexpo, Mint)
128
129
130 Mx = Mint[1] / pA
131 if Mx <= 0.0 or isNaN(Mx):
132 back_calc[i] = r20a
133 else:
134 back_calc[i]= -inv_tcpmg * log(Mx)
135