Re: CST branch -- October 27, 2010 - 20:46

Hi,
sorry for a late response. Your suggestions seem reasonable. I would add just a two comments.
1) We should keep as before not only ri_prime but probably also chi2, dchi2 and d2chi2 as well.
.
2)The consideration of Rex as an additional interaction is a one of the
suitable approaches and it fits the general idea. We have to consider
that it requires to add new interaction into the list of the
interactions ( __init__ function in class Mf, file maths_fns/mf.py). 
It might be done just within the __init__ function (class Mf, file
maths_fns/mf.py). It will be based on the presence of Rex term within the loaded
variable param_types. The other possibility is to do it before calling this __init__ function (number of interactions is one of its arguments). 

Then
I address to discussion one additional problem related to the patch in
the preparation. Within the  functions related to the derivations:
dfunc_mf, dfunc_local_tm, dfunc_all, dfunc_diff (class Mf,
maths_fns/mf.py) there is a small complication. It will be illustrated using  "dfunc_mf" function as an example. The original code is
following:

        data = "">
        # Calculate the spectral density gradient components.
        if data.calc_djw_comps:
            data.calc_djw_comps(data, params)

        # Diffusion tensor correlation times.

        self.diff_data.calc_dti(data, self.diff_data)

        # Loop over the gradient.
        for j in xrange(data.total_num_params):
            # Calculate the spectral density gradients.
            if data.calc_djw[j]:

                data.djw = data.calc_djw[j](data, params, j)
            else:
                data.djw = data.djw * 0.0

            # Calculate the relaxation gradient components.
            data.create_dri_comps(data, params)

            # Calculate the R1, R2, and sigma_noe gradients.
            data.dri_prime[j] = data.create_dri_prime[j](data)

            # Loop over the relaxation values and modify the NOE gradients.
            data.dri[j] = data.dri_prime[j]

            for m in xrange(data.num_ri):
                if data.create_dri[m]:
                    data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j)

            # Calculate the chi-squared gradient.

            data.dchi2[j] = dchi2_element(data.relax_data, data.ri, data.dri[j], data.errors)

If  the loop over interactions is incorporated, it is necessary to do it before
"data.calc_djw_comps" is calculated. However then we get into the
conflict with calculation of "modification of NOE gradients" and also
calculation of "dchi2". These quantities can be calculated only when
contributions from all interactions has been already evaluated (it is similar to already discussed ri_prime). In this case it is complicated by the presence of loop over fitted
parameters. The possible solution would be to make the loop over the
interactions the inner one. However this seems to require repeating the loop over
interactions twice within the code:

        data = "">
        # Loop over the interactions
        for k in xrange(self.num_interactions[0]):

            # Calculate the spectral density gradient components.
            if data[k].calc_djw_comps:
                data[k].calc_djw_comps(data[k], params)

            # Diffusion tensor correlation times.
            self.diff_data.calc_dti(data[k], self.diff_data)

        # Loop over the gradient.
        for j in xrange(data.total_num_params):
            # Loop over the interactions
            for k in xrange(self.num_interactions[0])
                # Calculate the spectral density gradients.
                if data[k].calc_djw[j]:
                    data[k].djw = data[k].calc_djw[j](data[k], params, j)
                else:
                    data[k].djw = data[k].djw * 0.0

                    # Calculate the relaxation gradient components.
                    data[k].create_dri_comps(data[k], params)

                # Calculate the R1, R2, and sigma_noe gradients.
                data[k].dri_prime[j] = data[k].create_dri_prime[j](data[k])                              # later data[k].dri_prime[j] will be replaced by data.dri_prime[j]
                                                                                                                               # but you suggest not to include both changes into the following patch together


            # Loop over the relaxation values and modify the NOE gradients.
            data.dri[j] = data.dri_prime[j]
            for m in xrange(data.num_ri):
                if data.create_dri[m]:
                    data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j)

                # Calculate the chi-squared gradient.
                data.dchi2[j] = dchi2_element(data.relax_data, data.ri, data.dri[j], data.errors)

Additionally I would notify that variable "total_num_params" needs to be kept at the position self.data[0] (similar to your suggestion in your last mail, concerning ri_prime)
What do you think about this?
Regards
Pavel

On 22 October 2010 19:46, Edward d'Auvergne <edward@xxxxxxxxxxxxx> wrote:

Hi,

This is a good idea.  First I would suggest putting this into 2
patches, as this is easier to check and apply.  Don't worry about
breaking the code at this point - it is in your own special branch so
as long as it works in the end, you can smash it into a million pieces
and put it back together again if you like.  Commits to the repository
are better when they are many small ones, as they can be more easily
managed.  For example if something is found to be designed not
ideally, or there is a fatal bug, that patch/commit can be reverted
and then new code can be worked on.  And small patches make it much
easier for the other relax developers to read and catch potential
hidden bugs or design issues.

This is an intriguing problem, as the data flow hits a fork here.
ri_prime is the correct target for the merging of the data from all
the interactions, as this needs to occur before the calculation of the
NOE using sigma_noe and the R1.  The relaxation rates R1, R2, and
sigma_noe add.  However the NOE does not.  I just thought I'd explain
the logic for others to follow ;)

I would suggest we store ri_prime somewhere else.  What I can do is to
apply a patch for the first change to accommodate for the multiple
relaxation interactions.  Then I could make a change myself, creating
a special Python object for 'data[i]'.  We can store the ri_prime data
this special object.  Essentially, it will be like the current code
base.  We could have specific interaction data in say:

self.data[10][0].jw  (or self.data[0].jw if data = "">
This is new.  But as before we could have:

self.data[10].ri_prime (or data.ri_prime if data = "">
So the forking can all be managed within the self.data[i] data
structures.  What do you think?  Also, we will have to work out how to
modify or replace func_ri_prime and func_ri_prime_rex in the
maths_fns.ri_prime module.  This is where the merging of data streams
currently occurs.  Maybe Rex needs to be considered as its own
interaction, as this is added to produce ri_prime?

Regards,

Edward


P. S.  I just talked to Kathleen Hall and she seemed very interested
in what you are doing here!

On 22 October 2010 18:26, Pavel Kaderavek <pavel.kaderavek@xxxxxxxxx> wrote:
> Hi,
> I am continuing in the discussion started in my post
> https://mail.gna.org/public/relax-devel/2010-09/msg00020.html
> It covers changes of functions func_mf.py, func_local_tm , func_diff,
> func_all and their equivalents for a first and second derivatives (class Mf,
> file maths_fns/mf.py).
>
> I would like to include into next patch also treatment of the fact, that it
> is necessary to sum together contributions of all interactions. It seems to
> me that the most suitable way is to do that by the modification just revised
> functions (func_mf.py, func_local_tm ...)
> I would suggest to initialize ri_prime before loop over interactions and
> then step by step add the contributions into the ri_prime:
>
>         data = ""                                # the cases
> when `i` is replaced by `0` were discussed in your last mail
>         ri_prime=0
>
>
>         for j in xrange(self.num_interactions[i]):
>
>               ...
>               ...
>
>               # Calculate the R1, R2, and sigma_noe values.
>               ri_prime = ri_prime + data[j].create_ri_prime(data[j])
>
>          data[0].ri_prime = ri_prime
>
>
> When the loop over interactions is finished the accumulated relaxation rate
> is copied into the data storage of the first interaction. Then it is
> possible to call functions, where total ri_prime is needed:
>
>         # Calculate the NOE values.
>         data[0].ri = data[0].ri_prime * 1.0
>         for m in xrange(data[0].num_ri):
>             if data[0].create_ri[m]:
>                 data[0].create_ri[m](data[0], m, data[0].remap_table[m],
> data[0].get_r1, params)
>
>         # Calculate the chi-squared value.
>         data[0].chi2 = chi2(data[0].relax_data, data[0].ri, data[0].errors)
>
> Regards,
> Pavel
>
> On 19 October 2010 13:49, Edward d'Auvergne <edward@xxxxxxxxxxxxx> wrote:
>>
>> Hi,
>>
>> Sorry for the delay, I just came back from a 2 week holiday.  This is
>> correct, the func_mf, func_local_tm, etc. methods are working on a
>> single spin.  This is stored in self.data[0].  The other functions
>> work on multiple spin data located in self.data[0], self.data[1], etc.
>>
>> Regards,
>>
>> Edward
>>
>>
>>
>> On 14 October 2010 19:38, Pavel Kaderavek <pavel.kaderavek@xxxxxxxxx>
>> wrote:
>> > Hi,
>> >
>> > I would like to announce a small clarification of my previous mail about
>> > changes in functions func_mf, func_local_tm, func_diff, func_all, and
>> > their
>> > derivatives (class Mf, file maths_fns/mf.py). Loop suggested in my last
>> > post:
>> >
>> > https://mail.gna.org/public/relax-devel/2010-09/msg00020.html
>> >
>> > is valid just for the functions func_diff, func_all, and their
>> > equivalents
>> > for a first and second derivatives.
>> > While for the functions  func_mf, func_local_tm, and corresponding
>> > derivatives the index of self.num_interactions should be set to `0`
>> > instead
>> > of index `i`
>> >
>> >         for j in xrange(self.num_interactions[0]):
>> >
>> > (The rest of the loop remains the same)
>> > That comes from the fact that these functions (func_mf, func_local_tm,
>> > ...
>> > ) are called for each spin separately. As it is indicated by the
>> > preceding
>> > statement:
>> >
>> >         data = ""> >> >
>> > Regards,
>> > Pavel
>> >
>> > On 29 September 2010 10:53, Edward d'Auvergne <edward@xxxxxxxxxxxxx>
>> > wrote:
>> >>
>> >> Hi,
>> >>
>> >> This is the perfect approach.  It will abstract the calculations so
>> >> that we will not need to touch many of maths_fns modules.  With this
>> >> code in place, I would aim at then making the test suite pass again by
>> >> having the correct data structures pass into maths_fns.mf.  The last
>> >> step would be to input CSA tensors and the multi-dipole interactions
>> >> via user functions.  If you make a patch for all of the func_*(),
>> >> dfunc_*(), and d2func_*() methods, I can check and apply it quickly.
>> >>
>> >> Regards,
>> >>
>> >> Edward
>> >>
>> >>
>> >> On 28 September 2010 12:07, Pavel Kaderavek <pavel.kaderavek@xxxxxxxxx>
>> >> wrote:
>> >> > Hi,
>> >> > we were thinking about next necessary changes in the CST branch.
>> >> >
>> >> > Now we need to break through the problem, which implies the fact we
>> >> > split
>> >> > the relaxation rate calculation into contributions of individual
>> >> > interactions. Each of them has its own data class to store its
>> >> > parameters
>> >> > (so far called data[i][j], where the [i] was a spin index and [j] was
>> >> > the
>> >> > interaction index).
>> >> >
>> >> > It seems to us, that it the best way to deal with it is to edit
>> >> > functions:
>> >> > func_mf, func_local_tm , func_diff, func_all and their equivalents
>> >> > for a
>> >> > first and second derivatives (defined in mf.py file).
>> >> >
>> >> > Within these functions the calculations of direction cosines,
>> >> > diffusion
>> >> > tensor weight calculations, components of the spectral densities and
>> >> > so
>> >> > on
>> >> > are performed. All these must be calculated for each interaction
>> >> > separately,
>> >> > because each interaction has its own data storage (which replaced
>> >> > previously
>> >> > used one data class container for the whole IS spin system).
>> >> >
>> >> > Instead of:
>> >> >
>> >> >
>> >> >         # Direction cosine calculations.
>> >> >         if self.diff_data.calc_di:
>> >> >             self.diff_data.calc_di(data, self.diff_data)
>> >> >
>> >> >         # Diffusion tensor weight calculations.
>> >> >         self.diff_data.calc_ci(data, self.diff_data)
>> >> >
>> >> >         # Diffusion tensor correlation times.
>> >> >         self.diff_data.calc_ti(data, self.diff_data)
>> >> >
>> >> >         # Calculate the components of the spectral densities.
>> >> >         if data.calc_jw_comps:
>> >> >             data.calc_jw_comps(data, params)
>> >> >
>> >> >         # Calculate the R1, R2, and sigma_noe values.
>> >> >         data.ri_prime = data.create_ri_prime(data)
>> >> >
>> >> >
>> >> > we would suggest to introduce a loop over interacations
>> >> >         for j in xrange(self.num_interactions[i]):
>> >> >             # Direction cosine calculations.
>> >> >             if self.diff_data.calc_di:
>> >> >                 self.diff_data.calc_di(data[j], self.diff_data)
>> >> >
>> >> >             # Diffusion tensor weight calculations.
>> >> >             self.diff_data.calc_ci(data[j], self.diff_data)
>> >> >
>> >> >             # Diffusion tensor correlation times.
>> >> >             self.diff_data.calc_ti(data[j], self.diff_data)
>> >> >
>> >> >             # Calculate the components of the spectral densities.
>> >> >             if data.calc_jw_comps:
>> >> >                 data.calc_jw_comps(data[j], params)
>> >> >
>> >> >             # Calculate the R1, R2, and sigma_noe components.
>> >> >             data.ri_prime = data.create_ri_prime(data[j])
>> >> >
>> >> >
>> >> > It must be accompanied in the next step by a change the ri_prime
>> >> > function so
>> >> > that it just calculate only a product of specific interaction
>> >> > constant
>> >> > and
>> >> > corresponding linear combination of spectral densities. While the
>> >> > final
>> >> > sumation over all interactions should be done in a separate step.
>> >> >
>> >> > Moreover it will be also necessary to distinguish within the function
>> >> > setup_equation the type of equation used for contribution of
>> >> > individual
>> >> > interactions according to their type.
>> >> >
>> >> > Best
>> >> > Pavel
>> >> >
>> >> >
>> >> > On 10 September 2010 15:43, Edward d'Auvergne <edward@xxxxxxxxxxxxx>
>> >> > wrote:
>> >> >>
>> >> >> Hi Pavel,
>> >> >>
>> >> >> I missed it in the patches, but there were tab characters '\t'
>> >> >> causing
>> >> >> problems.  These are now fixed.  relax requires that a tab is
>> >> >> replaced
>> >> >> by 4 spaces.  I have also added you to the copyright notices
>> >> >> (http://svn.gna.org/viewcvs/relax?view=rev&revision=11543) of the
>> >> >> files you have modified.
>> >> >>
>> >> >> Regards,
>> >> >>
>> >> >> Edward
>> >> >>
>> >> >> On 10 September 2010 15:36, Edward d'Auvergne <edward@xxxxxxxxxxxxx>
>> >> >> wrote:
>> >> >> > Hi,
>> >> >> >
>> >> >> > I've carefully checked the patches and committed them with the
>> >> >> > messages you wrote.  Sorry again for the delays.  It should be
>> >> >> > faster
>> >> >> > now that I am no longer in the tropical wilderness of far north
>> >> >> > Australia.
>> >> >> >
>> >> >> > Regards,
>> >> >> >
>> >> >> > Edward
>> >> >> >
>> >> >> >
>> >> >> > On 6 September 2010 13:23, Edward d'Auvergne
>> >> >> > <edward@xxxxxxxxxxxxx>
>> >> >> > wrote:
>> >> >> >> Hi Pavel,
>> >> >> >>
>> >> >> >> Sorry for the delayed response.  I was at the ICMRBS conference
>> >> >> >> in
>> >> >> >> Australia and then travelled through the tropical north end of
>> >> >> >> Australia afterwards.  I came back yesterday out of the remote
>> >> >> >> wilderness and can soon start looking at this patches.
>> >> >> >>
>> >> >> >> Regards,
>> >> >> >>
>> >> >> >> Edward
>> >> >> >>
>> >> >> >>
>> >> >> >> On 31 August 2010 18:00, Pavel Kaderavek
>> >> >> >> <pavel.kaderavek@xxxxxxxxx>
>> >> >> >> wrote:
>> >> >> >>> Hi,
>> >> >> >>>
>> >> >> >>> some time ago, we submitted two patches regarding CST branch. We
>> >> >> >>> are
>> >> >> >>> not
>> >> >> >>> sure if we should wait for some additional comment from your
>> >> >> >>> side,
>> >> >> >>> or
>> >> >> >>> we can
>> >> >> >>> continue with introducing further changes in the code.
>> >> >> >>> Next step would be a splitting of the relaxation equation so
>> >> >> >>> that
>> >> >> >>> contribution to the relaxation due to the individual types of
>> >> >> >>> interaction
>> >> >> >>> (dipole-dipole, CSA) can be calculated separately.
>> >> >> >>>
>> >> >> >>>
>> >> >> >>> Regars,
>> >> >> >>>
>> >> >> >>> Pavel, Petr
>> >> >> >>>
>> >> >> >>
>> >> >> >
>> >> >
>> >> >
>> >
>> >
>
>