Re: CST branchRe: CST branchHi, I think we should avoid optimising the bond length and the CSA tensor, if the full tensor is used or if there are multiple bond lengths for one spin, at least for now. This could be caught in the setup method and a RelaxError thrown saying that the optimisation of these things is not yet implemented. That will allow for the current models and the new code to work together. Optimisation of the CSA tensor or multiple dipole-dipole distances is too much effort for now, and there probably isn't enough information in the current data to extract this. What do you think? Cheers, Edward On 12 April 2011 15:28, Pavel Kaderavek <pavel.kaderavek@xxxxxxxxx> wrote:
Hi, it seems possible, but I think that it means that there will be for each type of the interaction an unique (d,d2)ri_comps function. It will be defined in the setup_equation (mf.py, class Mf). Currently, the (d,d2)ri_comps function is selected just based on the fact whether the values of physical quantities (internuclei distance, csa) are optimized or not. Now it will be selected again according to the same criteria and moreover according the type of interaction. This ensure that inside the function (d,d2)ri_comps the function create_(d,d2)ri_comps will be called with appropriate parameters. Generally, it is probably not very convenient to optimize all internuclei distances, but only the most crucial one (to the bonded hydrogen) would be sufficient. However, the selection should be upon user. The optimalization of CST will be also quite difficult, because 5 parameters would be needed to optimize (3 euler angles + 2 csa values - corresponding to two pseudo chemical shielding tensors which stand for asymmetric chemical shielding tensor - hence I always mention csa1 and csa2 as well as the cross-term, because they are not independent). So far I would consider only the optimalization of the eigenvalues of the tensors, the optimalization of the orientation might be added later. Best regards, Pavel On 22 March 2011 12:23, Edward d'Auvergne <edward@xxxxxxxxxxxxx> wrote:At this point we would like to address a related question. Currently the calculation of physical constant is done in a several steps. First, the physical constant is calculated and the value is stored in the data.dip_const_func or data.csa_const_func (grad, hess). Then, when the relaxation rates are calculated, the physical constant is recalculated by the function create_dip_func or create_csa_func (grad, hess) (method setup_equations in class Mf, maths_fns/mf.py). comp_dip_const_func(data, data.bond_length) comp_csa_const_func(data, data.csa) for i in xrange(data.num_ri): data.dip_comps_func[i] = data.dip_const_func if data.create_dip_func[i]: data.dip_comps_func[i] = data.create_dip_func[i](data.dip_const_func) if data.create_csa_func[i]: data.csa_comps_func[i] = data.create_csa_func[i](data.csa_const_func[data.remap_table[i]]) There is one exception, the dipolar physical constant is not recalculated in the case of calculation R1 relaxation rate, because the function create_dip_func does not exist in this case. We do not see a reason for such a recalculation.The reason is because of the m10 to m39 models built into relax. I have made it possible to optimise the bond length and CSA information. However these models are not stable with the current relaxation data. I do plan on working with these in the future though, so it would be useful to keep them. Note that for models m0 to m9, the data.create_dip_func[i] and data.create_csa_func[i] function pointers are set to None. Therefore for normal model-free analysis the constant is not recalculated.It seems better to us to just change the coefficient in the functions comp_r1_dip_jw, comp_r2_dip_jw, comp_r1_csa_jw, comp_r2_csa_jw (maths_fns/ri_comps.py). I guess, that this design was dedicated to avoid multiple calculation of the same interaction constant for each measured relaxation rate. We would suggest to reach the same effect by this construction: for i in xrange(data.num_ri): if data.const_func[0]: data.const_func[i] = data.const_func[0] else data.create_const_func(data)For models m10 to m39, this construct will not work. The constants are already pre-calculated for models m0 to m9 so this is not needed.Note, the comp_dip_const_func and comp_csa_const_func should be change so that, it is possible to call them just with the argument data (maths_fns/ri_comps.py). Instead of: def comp_dip_const_func(data, bond_length): """Calculate the dipolar constant. ... if bond_length == 0.0: data.dip_const_func = 1e99 else: data.dip_const_func = 0.25 * data.dip_const_fixed * bond_length**-6 It should look like: def comp_dip_const_func(data): """Calculate the dipolar constant. ... if data.internuclei_distance == 0.0: data.const_func = 1e99 else: data.const_func = 0.25 * data.dip_const_fixed * data.internuclei_distance**-6The bond_length arg was designed so that the bond length could either come from a fixed value supplied by the user or from the parameter vector when bond lengths or CSA values are optimised. This behaviour might have to be preserved.data.dip_const_func were renamed to more general data.const_func and instead of bond_length the function directly takes the internuclei distance for the current dipole-dipole interaction. The change of data.dip_const_func to data.const_func later simplify the code design in the maths_fns/ri_prime.py . It will be reduced just to a multiplication of constant and the linear combination of spectral density functions.For models m10 to m39, I'm not sure if this design would work. Could we redesign this in another way in which these complex models are still functional?Then, we would suggest to call the function comp_dip_const_func, comp_csa_const_func ... with full set of possible parameters, i.e. comp_dip_const_func(data,internuclei_distance,csa1,csa2,rex) comp_csa_const_func(data,internuclei_distance,csa1,csa2,rex) ... If we would call the function with just physically relevant arguments then we would have to use yet another condition to decide the type of the interaction in the loop where individual interaction contributions are calculated. The physical quantities irrelevent for the given interaction (for example csa1,csa2,rex for the dipole-dipole interaction) are None from the initialization and will not be used by the function anyway.It would be easier to have a different function for each physical constant, as this is very specific code. For example comp_dip_const_func is only called from functions which require the dipolar constant. For the normal model-free models, this is only called during setup. So there is no need to pass in the csa1 and csa2 constants, or Rex for this specific function. For the CSA constants, do you calculate 2 separate constants? For example would you use: comp_dip_const_func(data[i],internuclei_distance) comp_csa_const_func(data[j],csa1) comp_csa_const_func(data[k],csa2) where data[x] in each case is a different interaction?Moreover, there is an unanswered question about the NOE and the additional dipolar interaction. I am not sure if the suggested design is physically correct, rather not. During the NOE experiment, the protons are saturated in order to reach the steady state. Then a complete set cross relaxation rates between all interacting spin pairs should be taken into the account, not only between the spin of interest and all other interacting nuclei. On the other hand this is probably beyond the aim of the program relax. What do you think about that?This is getting quite complex as you would need to take the cross-correlated relaxation rates between the different relaxation interactions into account, as well as the motion of all spins if they are not directly bonded. Is this needed for the current work? Of course anything is accepted into relax, especially if you would like to probe this area (with a paper in mind), but it has to play nicely with the rest of relax and not be a burden on the relax developers to maintain in the future (as well as not make the current number crunching code slower than it already is). The code would therefore need to be designed in public. So if you would like to tackle such a task, I would first recommend finishing off the cst branch, and then make a new branch for this work.Then, I will design the code as I suggested. So, the sigma_noe will be calculated separately for all dipole-dipole interactions with the central spin (assuming they are isolated spin pairs). The total sigma_noe will be calculated as a sum of all individual sigma_noe contributions.This sounds the most reasonable for a first approximation. In reality the cross rates will be important as well. Cheers, Edward