Chemical kinetics parameters

From Spinach Documentation Wiki
Jump to: navigation, search

Spinach has a very general chemical kinetics module that can handle arbitrary reaction networks, the only restriction being that the corresponding differential equations must be linear and must have the following general form:

[math]\frac{d}{dt}\left( \begin{matrix} [\text{A}] \\ [\text{B}] \\ [\text{C}] \\ \vdots \\ \end{matrix} \right)=\mathbf{K}\left( \begin{matrix} [\text{A}] \\ [\text{B}] \\ [\text{C}] \\ \vdots \\ \end{matrix} \right)[/math]

where [math]\mathbf{K}[/math] is the reaction rate matrix. For example:

[math]\begin{matrix} \text{A }\xrightarrow{{{k}_{1+}}}\text{ B }\underset{{{k}_{2-}}}{\overset{{{k}_{2+}}}{\longleftrightarrow}}\text{ C} \\ \frac{d}{dt}\left( \begin{matrix} \text{ }\!\![\!\!\text{ A }\!\!]\!\!\text{ } \\ \text{ }\!\![\!\!\text{ B }\!\!]\!\!\text{ } \\ \text{ }\!\![\!\!\text{ C }\!\!]\!\!\text{ } \\ \end{matrix} \right)=\left( \begin{matrix} -{{k}_{1+}} & 0 & 0 \\ {{k}_{1+}} & -{{k}_{2+}} & {{k}_{2-}} \\ 0 & {{k}_{2+}} & -{{k}_{2-}} \\ \end{matrix} \right)\left( \begin{matrix} \text{ }\!\![\!\!\text{ A }\!\!]\!\!\text{ } \\ \text{ }\!\![\!\!\text{ B }\!\!]\!\!\text{ } \\ \text{ }\!\![\!\!\text{ C }\!\!]\!\!\text{ } \\ \end{matrix} \right) \\ \end{matrix}[/math]

Spinach expects you to supply this matrix and the initial concentrations of the molecules, as well as to say which spins in the spin system specification belong to which molecule. Spinach can also assist in determining the steady-state concentrations (see equilibrate.m) if you do not know them.

Networks of chemical reactions

If a number of molecules are involved in a network of chemical reactions, all of the molecules should be specified in the same input stream and then Spinach should be told which spins belong to which molecule. Example:

    % Isotopes
    sys.isotopes={'1H','13C','15N','1H',13C','15N'};
    
    % Chemical shifts
    inter.zeeman.scalar={1.0, 20.0, 15.0, 1.5, 25.0, 16.0};

    % Spins 1,2,3 are molecule A; spins 4,5,6 are molecule B
    inter.chem.parts={[1 2 3],[4 5 6]};
    
    % Kinetic rate matrix (Hz)
    inter.chem.rates=[-0.1 0.2; 0.1 -0.2];

    % Initial concentrations (arbitrary units)
    inter.chem.concs=[1.0, 2.0] 

In the general case, the the parameters, supplied to create.m at the spin system specification stage, must be:

Reaction kinetics parameters
Variable name Variable type Content
inter.chem.parts a cell array of vectors containing integers individual vectors in the cell array must contain the numbers of spins that belong to each of the molecules in the chemical reaction, for example {[1 2],[3 4]} indicates that spins 1 and 2 belong to the first molecule and spins 3 and 4 belong to the second molecule.
inter.chem.rates a matrix of real numbers chemical reaction rate matrix between the molecules identified in inter.chem.parts
inter.chem.concs a vector of non-negative real numbers initial concentrations of the molecules identified in inter.chem.parts

Note that all molecules must have the same number of spins, must have those spins specified in the same order, and must have the same basis set. Only the Liouville space formalism using the spherical tensor basis set is supported for this option.

Spin state fluxes

There are rare cases when the magnetization is transferred through some indirect mechanism within a topologically fixed molecule, for example, a proton swap between two N-H groups of the same molecule. Because protons are fundamentally identical, this essentially amounts to swapping the spin state (together with all associated correlations to all external spins) between the two protons. In this case, the user should supply the matrix of magnetization flux rates between spins, for example:

    % Isotopes
    sys.isotopes={'15N','1H','1H'};

    % Chemical shift tensors
    inter.zeeman.eigs{1}=[-50 -50 100];
    inter.zeeman.euler{1}=[0.1 0.2 0.3];
    inter.zeeman.eigs{2}=[-5 -5 10]+1;
    inter.zeeman.euler{2}=[0.2 0.3 0.1];
    inter.zeeman.eigs{3}=[-5 -5 10]+5;
    inter.zeeman.euler{3}=[0.3 0.1 0.2];

    % J-couplings
    inter.coupling.scalar=cell(3,3);
    inter.coupling.scalar{1,2}=50;

    % Coordinates
    inter.coordinates={[0.00 0.00 0.00];
                       [1.05 0.00 0.00];
                       [0.00 3.50 0.00]};

    % Swapping the protons and their correlations
    inter.chem.flux_rate=zeros(3,3);
    inter.chem.flux_rate(2,3)=1000;
    inter.chem.flux_rate(3,2)=1000;
    inter.chem.flux_type='intramolecular';

In the general case, the the parameters, supplied to create.m at the spin system specification stage, must be:

Magnetisation flux parameters
Variable name Variable type Content
inter.chem.flux_rate nspins × nspins matrix of real numbers magnetization transport rates, in Hz
inter.chem.flux_type string: 'intermolecular', 'intramolecular' controls the fate of correlations with external spins, which are damped in the 'intermolecular' case and transported in the 'intramolecular' case.

The flux type parameter requires further explanation. When we swap a pair of spins in a molecule, their coherences with any observer spins should in general be preserved. This is what the 'intramolecular' option does. There are cases, however, when the exchange happens through some external medium, such as the solvent, or the end point of the jump is a different molecule. This would lead to the emergence of inter-molecular coherences between identical molecules – a situation that the density matrix formalism is ill-equipped to accommodate. Because such coherences are non-observable, and the spin is very unlikely to jump back to the exact same molecule that it had originally come from, the coherences in question are counted as lost. The corresponding contribution to the kinetics superoperator effectively causes relaxation. This is what happens when the 'intermolecular' option is specified.

Spin-selective radical pair recombination

This exotic process occurs in radical pairs involved in CIDNP, where the two-electron singlet state recombines (and forms diamagnetic products) faster than the two-electron triplet state. Because the rate of the singlet-triplet interconversion depends on the hyperfine couplings present, the result is the sorting of nuclei by orientation (some nuclear spin states accumulate in diamagnetic products, others are left behind in the unreacted radicals).

The kinetic modelling of the process is reviewed in the recent paper by Peter Hore and co-authors (http://dx.doi.org/10.1016/j.cplett.2011.03.082). The exponential model, the Haberkorn model and the Jones-Hore model are currently supported in Spinach. An example of the specification appears below:

    % Particles present
    sys.isotopes ={'E','E','1H'};

    % g-factors and chemical shifts
    inter.zeeman.scalar={2.0023 2.0024 1.0};

    % Hyperfine couplings
    inter.coupling.scalar=cell(3,3);
    inter.coupling.scalar{2,3}=1e7;

    % Radical pair kinetics model
    inter.chem.rp_theory='haberkorn';

    % Electron identification
    inter.chem.rp_electrons=[1 2];

    % Singlet and triplet recombination rates 
    inter.chem.rp_rates=[1e7 0];

In the general case, the the parameters, supplied to create.m at the spin system specification stage, must be:

Magnetisation flux parameters
Variable name Variable type Content
inter.chem.rp_theory string: 'exponential', 'haberkorn', 'jones-hore' radical pair recombination theory
inter.chem.rp_electrons 1 x 2 vector of integers rlectron spins involved in the recombination reaction.
inter.chem.rp_rates 1 x 2 vector of reals singlet and triplet recombination rates in a radical pair, in Hz

The theory is only applicable to systems with exactly two unpaired electrons in Liouville space formalisms.


Version 1.8, authors: Ilya Kuprov, Hannah Hogben, Peter Hore