Hi Troels,
I have spotted a bug in this change. The reason it is a bug is
because of the 'MMQ CR72' equations
(http://www.nmr-relax.com/manual/MMQ_CR72_model.html). A value of dwH
== 0.0 does not result in no exchange if dw != 0.0. The converse is
also true, there is exchange when dw == 0.0 if dwH != 0.0. So you
could modify the check to be:
# Catch parameter values that will result in no exchange,
returning flat R2eff = R20 lines (when kex = 0.0, k_AB = 0.0).
if (dw == 0.0 and dwH == 0.0) or pA == 1.0 or k_AB == 0.0:
return array([r20]*num_points)
You would need to modify the unit tests to match.
Regards,
Edward
On 26 May 2014 23:09, <[email protected]> wrote:
> Author: tlinnet
> Date: Mon May 26 23:09:43 2014
> New Revision: 23439
>
> URL: http://svn.gna.org/viewcvs/relax?rev=23439&view=rev
> Log:
> Critical fix for the math domain catching of model MMQ CR72.
>
> This was discovered with the added 9 unit tests demonstrating edge case 'no
> Rex' failures.
>
> This follows from the ideas in the post
> http://article.gmane.org/gmane.science.nmr.relax.devel/5858.
> This is related to: task #7793: (https://gna.org/task/?7793) Speed-up of
> dispersion models.
>
> This is to implement catching of math domain errors, before they occur.
>
> Modified:
> branches/disp_speed/lib/dispersion/mmq_cr72.py
>
> Modified: branches/disp_speed/lib/dispersion/mmq_cr72.py
> URL:
> http://svn.gna.org/viewcvs/relax/branches/disp_speed/lib/dispersion/mmq_cr72.py?rev=23439&r1=23438&r2=23439&view=diff
> ==============================================================================
> --- branches/disp_speed/lib/dispersion/mmq_cr72.py (original)
> +++ branches/disp_speed/lib/dispersion/mmq_cr72.py Mon May 26 23:09:43
> 2014
> @@ -84,6 +84,10 @@
> @type power: numpy int16, rank-1 array
> """
>
> + # Catch parameter values that will result in no exchange, returning flat
> R2eff = R20 lines (when kex = 0.0, k_AB = 0.0).
> + if dw == 0.0 or pA == 1.0 or k_AB == 0.0 or dwH == 0.0:
> + return array([r20]*num_points)
> +
> # Repetitive calculations (to speed up calculations).
> dw2 = dw**2
> r20_kex = r20 + kex/2.0
> @@ -126,9 +130,7 @@
> # Catch math domain error of cosh(val > 710).
> # This is when etapos > 710.
> if max(etapos) > 700:
> - R2eff = array([1e100]*num_points)
> -
> - return R2eff
> + return array([r20]*num_points)
>
> # The full eta - values.
> etaneg = etaneg_part / cpmg_frqs
> @@ -154,4 +156,4 @@
> if not isfinite(sum(R2eff)):
> R2eff = array([1e100]*num_points)
>
> - return R2eff
> + return R2eff
>
>
> _______________________________________________
> relax (http://www.nmr-relax.com)
>
> This is the relax-commits mailing list
> [email protected]
>
> To unsubscribe from this list, get a password
> reminder, or change your subscription options,
> visit the list information page at
> https://mail.gna.org/listinfo/relax-commits
_______________________________________________
relax (http://www.nmr-relax.com)
This is the relax-devel mailing list
[email protected]
To unsubscribe from this list, get a password
reminder, or change your subscription options,
visit the list information page at
https://mail.gna.org/listinfo/relax-devel
