Modified program for chemical exchange

[katarkon]

I have made some modification in Stefan Stoll's exchange program viewtopic.php?f=8&t=309. The major changes are follows:
- separate linewidth for each site
- automatic calculation of each site population (manual definition is possible too)
- multicomponent simulation (spectra without exchange are possible too)
- natural isotope mixtures are supported
- gaussian line broadening (now it's single for whole spectrum, use Sys.lwg parameter)
- tumbling effect compensation (via Sys.lb parameter and Sys.lbm matrix)
- frequency drift compensation (via Sys.mwFreq parameter)
- alternative definition for g, A and k to limit the number of parameters for fitting (or alternatively to define the restrictions)
- temperature support (now via simple 1/T scalling factor, more precise calculation may be issued if necessary)
- alternative definition of k via activation entрalpy, entropy and temperature
- more precise line position calculations via Breit-Rabi or second-order perturbation methods (equivalent nuclei not supported yet)
If you use this code, please cite original Stefan Stoll's work http://dx.doi.org/10.1021/jp3104358 and my work http://dx.doi.org/10.1016/j.molstruc.2017.07.001 where modified code was firstly used.
Alternative definition for g, A and k is follow. Sys.g, Sys.A and Sys.k are lists of values. Additional parameters Sys.gm, Sys.Am and Sys.km are required. They should be a integer numbers meaning the index of corresponding value. I.e. Sys.Am(i,j)=2 means second value from Sys.A. So, one value may be used several times. For Sys.km a special value 0 (zero) is also defined. It means the absence of corresponding exchange.
Activation entрalpy Sys.dH and entropy Sys.dS may be used instead Sys.k, the Exp.Temperature is additionally required.
Tumbling effect is defined by two values Sys.lb and Sys.lbm. Sys.lb means max line broadening. Sys.lbm contains values between 0 and 1 meaning fractional broadening for each line of each site. I.e. Sys.lbm(i,j)=0.5 means that j-th line of i-th site have a half of max broadening. At least one of values shold be 1 and 0. So, Sys.lbm should be defined once before fitting and only Sys.lb should be varied.
For nuclei with large hyperfine splitting two alternative methods are proposed. Use Opt.CalcMethod='breitrabi' for Breit-Rabi treatment or Opt.CalcMethod='perturb2' for second-order perturbation treatment. Equivalent nuclei are not supported. Cross-nuclear terms are neglected for both methods, so only system with one large hyperfine splitting will be treated correctly. The matrix diagonalization method is still under development. You may use it by defining Opt.CalcMethod='matrix', but it may stuck or give unexpected reults. Use it for own risk.
Also, simple wrapper exchange2D for temperature-sliced spectra is issued. It uses parameters Sys.dH and Sys.dS for activation enthalpy and entropy and set of temperatures in Exp.Temperature parameter to calculate the rate constants and returns a sequence of the spectra as concatenated "megaspectrum". The wrapper is attended for using together with esfit function. The parameters Sys.lw, Sys.lb and Sys.mwFreq generally should be arrays of values corresponding to each temperature.

[katarkon]

Some brief instructions for new features.
Alternative way for definition of the rate constants allow to restrict the number of varying parameters and define the exchange model with restrictions. For example, let's build the matrix for three-site exchange A <-> B <-> A' with two rate constants A->B and B->A (or B->A' which equals B->A). In common way it will be Sys.k=[0 k1 0;k2 0 k2;0 k1 0];. The function esfit will vary four values of Sys.k, and their equivalence will be lost. In alternative way Sys.k=[k1 k2];, Sys.km=[0 1 0;2 0 2;0 1 0]. Only two parameters need to be varied by this way. The parameters Sys.g and Sys.A may be defined in similar manner:

Code: Select all

Sys.A=[A1 A2 A3];
Sys.Am(1,:)=[1 2];
Sys.Am(2,:)=[3 3];
Sys.Am(3,:)=[2 1];

Also, activation enthalpy and entropy may be used for calculation of rate constants together with Exp.Temperature parameter. Corresponding parameters Sys.dH (in J/mol) and Sys.dS (in J/mol*K). They should have the same size as desired Sys.k. The best way is use they together with Sys.km as noted above.
Arbitrary line broadening may be defined by two parameters Sys.lb and Sys.lbm. Sys.lb means maximal broadening in the spectrum (i.e. difference between max and min linewidths in the spectrum). Sys.lbm is a matrix containing the coefficients for each line of each site in the range 0...1. Result linewidth calculates by equation lw(iSite,iLine)=Sys.lw(iSite)+Sys.lb(iSite,iLine). This allows vary only two parametersSys.lw and Sys.lb. Sys.lbm values may be adjusted manually or calucated if anisotropic parameters are known:

Code: Select all

lw=fastmotion(Sys,B0,Sys.tcorr);
lw=lw-min(lw);
lw=lw/max(lw);

This procedure should be made once for each site. B0 value is the mean of the range, Sys.tcorr is arbitrary but reasonable value i.e. 1e-9.
For isotopes mixture separate Sys.lbm matrix for each isotopologue is required. It may be defined as cell array so Sys.lbm{i} is Sys.lbm for i-th isotopologue. If all of isotopes have the same spin (like Cu, Cl, Br), single value of Sys.lb may be used.
Sys.mwFreq is alternative way to define of MW frequency. In some cases it may be adjusted during EPR experiment (e.g. temperature sweep) and not match with recorded value for some slices. Such definition allow varying of MW frequency during the fitting. Sys.g and Sys.mwFreqshould not be varied together, of course.
Alternative methods for calculation of lines positions may be used via definition of Opt.CalcMethod parameter. It may have the values 'perturb1' (default), 'perturb2', 'breitrabi' or 'matrix' coresponding to first or second order perturbation, Breit-Rabi or matrix diagonalization correspondingly. Only first method works with equivalent nuclei. Equivcouple representation for Breit-Rabi and second-order perturbation methods does not released yet (and most probably never will be released). Each nucleus should be defined separately. Matrix diagonalization method now is under development and may be used only for your own risk. But it have no alternatives if several nuclei with large hyperfine coupling are present.
Additional wrapper exchange2D is attended for temperature-sliced spectra treatment with esfit(). It calculates concatenated 'long' spectrum for temperatures defined in Exp.Temperature array. Also, Sys.dH and Sys.dS should be used instead Sys.k. Because line width changes with temperature, Sys.lw generally should be an array with the same size as Exp.Temperature. Sys.lb and Sys.mwFreq should be defined in the same manner if necessary. A little example how to use exchange2D:

Code: Select all

[x, y, par] = eprload('filename.par');
Exp.Temperature=x{2}; %temperature list
x=x{1}; %field values
y=reshape(y,1,[]); %long concatenated spectrum
Sys=...
... % spin system, parameters, options, etc
esfit('exchange2D',y,Sys,Vary,Exp,Opt,FitOpt);