Source Code for Module lib.dispersion.lm63

  1  ############################################################################### 
  2  #                                                                             # 
  3  # Copyright (C) 2009 Sebastien Morin                                          # 
  4  # Copyright (C) 2013-2014 Edward d'Auvergne                                   # 
  5  # Copyright (C) 2014 Troels E. Linnet                                         # 
  6  #                                                                             # 
  7  # This file is part of the program relax (http://www.nmr-relax.com).          # 
  8  #                                                                             # 
  9  # This program is free software: you can redistribute it and/or modify        # 
 10  # it under the terms of the GNU General Public License as published by        # 
 11  # the Free Software Foundation, either version 3 of the License, or           # 
 12  # (at your option) any later version.                                         # 
 13  #                                                                             # 
 14  # This program is distributed in the hope that it will be useful,             # 
 15  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
 16  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
 17  # GNU General Public License for more details.                                # 
 18  #                                                                             # 
 19  # You should have received a copy of the GNU General Public License           # 
 20  # along with this program.  If not, see <http://www.gnu.org/licenses/>.       # 
 21  #                                                                             # 
 22  ############################################################################### 
 23   
 24  # Module docstring. 
 25  r"""The Luz and Meiboom (1963) 2-site fast exchange U{LM63<http://wiki.nmr-relax.com/LM63>} model. 
 26   
 27  Description 
 28  =========== 
 29   
 30  This module is for the function, gradient and Hessian of the U{LM63<http://wiki.nmr-relax.com/LM63>} model. 
 31   
 32   
 33  References 
 34  ========== 
 35   
 36  The model is named after the reference: 
 37   
 38      - Luz, S. and Meiboom S., (1963)  Nuclear Magnetic Resonance study of protolysis of trimethylammonium ion in aqueous solution - order of reaction with respect to solvent, I{J. Chem. Phys.}, B{39}, 366-370 (U{DOI: 10.1063/1.1734254<http://dx.doi.org/10.1063/1.1734254>}). 
 39   
 40   
 41  Equations 
 42  ========= 
 43   
 44  The equation used is:: 
 45   
 46                    phi_ex   /     4 * nu_cpmg         /     kex     \ \  
 47      R2eff = R20 + ------ * | 1 - -----------  * tanh | ----------- | | , 
 48                     kex     \         kex             \ 4 * nu_cpmg / / 
 49   
 50  where:: 
 51   
 52      phi_ex = pA * pB * delta_omega^2 , 
 53   
 54  kex is the chemical exchange rate constant, pA and pB are the populations of states A and B, and delta_omega is the chemical shift difference between the two states. 
 55   
 56   
 57  Links 
 58  ===== 
 59   
 60  More information on the LM63 model can be found in the: 
 61   
 62      - U{relax wiki<http://wiki.nmr-relax.com/LM63>}, 
 63      - U{relax manual<http://www.nmr-relax.com/manual/The_LM63_2_site_fast_exchange_CPMG_model.html>}, 
 64      - U{relaxation dispersion page of the relax website<http://www.nmr-relax.com/analyses/relaxation_dispersion.html#LM63>}. 
 65  """ 
 66   
 67  # Python module imports. 
 68  from numpy import isfinite, min, sum, tanh 
 69  from numpy.ma import fix_invalid, masked_where 
 70   
 71   
 72 -def r2eff_LM63(r20=None, phi_ex=None, kex=None, cpmg_frqs=None, back_calc=None): 
 73      """Calculate the R2eff values for the LM63 model. 
 74   
 75      See the module docstring for details. 
 76   
 77   
 78      @keyword r20:           The R20 parameter value (R2 with no exchange). 
 79      @type r20:              numpy float array of rank [NE][NS][NM][NO][ND] 
 80      @keyword phi_ex:        The phi_ex parameter value (pA * pB * delta_omega^2). 
 81      @type phi_ex:           numpy float array of rank [NE][NS][NM][NO][ND] 
 82      @keyword kex:           The kex parameter value (the exchange rate in rad/s). 
 83      @type kex:              float 
 84      @keyword cpmg_frqs:     The CPMG nu1 frequencies. 
 85      @type cpmg_frqs:        numpy float array of rank [NE][NS][NM][NO][ND] 
 86      @keyword back_calc:     The array for holding the back calculated R2eff values.  Each element corresponds to one of the CPMG nu1 frequencies. 
 87      @type back_calc:        numpy float array of rank [NE][NS][NM][NO][ND] 
 88      """ 
 89   
 90      # Flag to tell if values should be replaced if phi_ex is zero. 
 91      t_phi_ex_zero = False 
 92   
 93      # Catch divide with zeros (to avoid pointless mathematical operations). 
 94      if kex == 0.0: 
 95          back_calc[:] = r20 
 96          return 
 97   
 98      # Catch zeros (to avoid pointless mathematical operations). 
 99      # This will result in no exchange, returning flat lines. 
100      if min(phi_ex) == 0.0: 
101          t_phi_ex_zero = True 
102          mask_phi_ex_zero = masked_where(phi_ex == 0.0, phi_ex) 
103   
104      # Repetitive calculations (to speed up calculations). 
105      rex = phi_ex / kex 
106      kex_4 = 4.0 / kex 
107   
108      # Calculate R2eff. 
109      back_calc[:] = r20 + rex * (1.0 - kex_4 * cpmg_frqs * tanh(kex / (4.0 * cpmg_frqs))) 
110   
111      # Replace data in array. 
112      # If phi_ex is zero. 
113      if t_phi_ex_zero: 
114          back_calc[mask_phi_ex_zero.mask] = r20[mask_phi_ex_zero.mask] 
115   
116      # Catch errors, taking a sum over array is the fastest way to check for 
117      # +/- inf (infinity) and nan (not a number). 
118      if not isfinite(sum(back_calc)): 
119          # Replaces nan, inf, etc. with fill value. 
120          fix_invalid(back_calc, copy=False, fill_value=1e100) 
121