Hi Edward.
There is a very very smart trick for that.
Masked arrays!
One can:
mask_dw_t = False
# Test if dw is zero.
if allclose(dw, zeros(dw.shape)):
mask_dw_t = True
mask_dw = ma.masked_values(dw, 0.0)
# Calculate R2eff.
R2eff = r20_kex - cpmg_frqs * arccosh( fact )
# Replace data in array.
if mask_dw_t:
R2eff[mask_dw] = r20a
back_calc[:] = R2eff
2014-06-10 13:52 GMT+02:00 Edward d'Auvergne <[email protected]>:
> Hi,
>
> I think the 'no Rex' tests and what R2eff is set to needs to be
> thought out carefully. We have two cases:
>
> - The spin independent parameters pA and kex. These will affect
> the R2eff values for all spins of the cluster.
>
> - The spin dependent parameter dw. This will effect the R2eff
> values for individual spins, and catching this and dealing with it
> correctly might be much harder for a spin cluster.
>
> The first case is simplest, catch single pA and kex values. So numpy
> operations are not needed there. For the second, I'm not sure how we
> handle that. But more importantly we also need to work out what to do
> with the R2eff values in that case. How do we catch dw = 0.0 for one
> spin but not the others, and then set the R2eff values to R20 for just
> that spin? But also let the R2eff values be calculated for all other
> spins. Maybe it is better to not catch dw values and to deal with
> that value later in the lib.dispersion functions?
>
> Regards,
>
> Edward
>
>
>
>
> On 10 June 2014 12:05, Troels Emtekær Linnet <[email protected]> wrote:
>> Here are some hints:
>>
>> http://stackoverflow.com/questions/10580676/comparing-two-numpy-arrays-for-equality
>>
>> It would be the timing of:
>>
>> numpy.allclose
>> http://docs.scipy.org/doc/numpy/reference/generated/numpy.allclose.html
>>
>> numpy.isclose
>> http://docs.scipy.org/doc/numpy/reference/generated/numpy.isclose.html#numpy.isclose
>>
>> numpy.ma
>> http://docs.scipy.org/doc/numpy/reference/generated/numpy.ma.masked_values.html#numpy.ma.masked_values
>>
>>
>> 2014-06-10 11:57 GMT+02:00 Troels Emtekær Linnet <[email protected]>:
>>> Hi Edward.
>>>
>>> I thought the same.
>>>
>>> One could send in a ones and zeros numpy array in function.
>>> Set them to standard to None.
>>>
>>> If they are None, then create them.
>>>
>>> One can see on the timings, that
>>> r2eff_CR72
>>> numeric.py:2056(allclose)
>>>
>>> are the most time consuming.
>>> Maybe there is a faster way for allclose.
>>> Masking ???
>>>
>>>
>>> 1 0.000 0.000 2.084 2.084 <string>:1(<module>)
>>> 1 0.002 0.002 2.084 2.084 profiling_cr72.py:441(cluster)
>>> 1000 0.002 0.000 2.005 0.002 profiling_cr72.py:405(calc)
>>> 1000 0.034 0.000 2.003 0.002
>>> relax_disp.py:995(func_CR72_full)
>>> 1000 0.141 0.000 1.960 0.002
>>> relax_disp.py:524(calc_CR72_chi2)
>>> 1300 1.100 0.001 1.676 0.001 cr72.py:101(r2eff_CR72)
>>> 4300 0.245 0.000 0.507 0.000 numeric.py:2056(allclose)
>>> 3000 0.037 0.000 0.164 0.000 shape_base.py:761(tile)
>>> 17828 0.126 0.000 0.126 0.000 {method 'reduce' of
>>> 'numpy.ufunc' objects}
>>> 4000 0.110 0.000 0.110 0.000 {method 'repeat' of
>>> 'numpy.ndarray' objects}
>>> 8609 0.011 0.000 0.086 0.000 fromnumeric.py:1762(any)
>>>
>>>
>>>
>>>
>>>
>>> 2014-06-10 11:44 GMT+02:00 Edward d'Auvergne <[email protected]>:
>>>
>>>> Hi Troels,
>>>>
>>>> For speed, maybe you should send in the single value parameters
>>>> together with the array versions into the lib.dispersion modules? The
>>>> check:
>>>>
>>>> + if np.allclose(dw, np.zeros(dw.shape)) or np.allclose(pA,
>>>> np.ones(dw.shape)) or np.allclose(kex, np.zeros(dw.shape)):
>>>>
>>>> will be very expensive.
>>>>
>>>> Regards,
>>>>
>>>> Edward
>>>>
>>>>
>>>>
>>>> On 8 June 2014 19:48, <[email protected]> wrote:
>>>> > Author: tlinnet
>>>> > Date: Sun Jun 8 19:48:35 2014
>>>> > New Revision: 23737
>>>> >
>>>> > URL: http://svn.gna.org/viewcvs/relax?rev=23737&view=rev
>>>> > Log:
>>>> > Re-implemented safety checks in lib/dispersion/cr72.py.
>>>> >
>>>> > This is now implemented for both rank-1 float array and of higher
>>>> > dimensions.
>>>> >
>>>> > This makes the unit tests pass for multi dimensional computing.
>>>> >
>>>> > Task #7807 (https://gna.org/task/index.php?7807): Speed-up of dispersion
>>>> > models for Clustered analysis.
>>>> >
>>>> > Modified:
>>>> > branches/disp_spin_speed/lib/dispersion/cr72.py
>>>> >
>>>> > Modified: branches/disp_spin_speed/lib/dispersion/cr72.py
>>>> > URL:
>>>> > http://svn.gna.org/viewcvs/relax/branches/disp_spin_speed/lib/dispersion/cr72.py?rev=23737&r1=23736&r2=23737&view=diff
>>>> >
>>>> > ==============================================================================
>>>> > --- branches/disp_spin_speed/lib/dispersion/cr72.py (original)
>>>> > +++ branches/disp_spin_speed/lib/dispersion/cr72.py Sun Jun 8
>>>> > 19:48:35 2014
>>>> > @@ -128,9 +128,15 @@
>>>> > rank_1 = False
>>>> >
>>>> > # Catch parameter values that will result in no exchange, returning
>>>> > flat R2eff = R20 lines (when kex = 0.0, k_AB = 0.0).
>>>> > + # For rank-1 float array.
>>>> > if rank_1:
>>>> > if dw == 0.0 or pA == 1.0 or kex == 0.0:
>>>> > back_calc[:] = array([r20a]*num_points)
>>>> > + return
>>>> > + # For higher dimensions, return same structure.
>>>> > + else:
>>>> > + if np.allclose(dw, np.zeros(dw.shape)) or np.allclose(pA,
>>>> > np.ones(dw.shape)) or np.allclose(kex, np.zeros(dw.shape)):
>>>> > + back_calc[:] = r20a
>>>> > return
>>>> >
>>>> > # The B population.
>>>> > @@ -165,16 +171,23 @@
>>>> >
>>>> > # Catch math domain error of cosh(val > 710).
>>>> > # This is when etapos > 710.
>>>> > - if rank_1:
>>>> > - if max(etapos) > 700:
>>>> > + if max(etapos) > 700:
>>>> > + if rank_1:
>>>> > back_calc[:] = array([r20a]*num_points)
>>>> > + return
>>>> > + # For higher dimensions, return same structure.
>>>> > + else:
>>>> > + back_calc[:] = r20a
>>>> > return
>>>> >
>>>> > # The arccosh argument - catch invalid values.
>>>> > fact = Dpos * cosh(etapos) - Dneg * cos(etaneg)
>>>> > - if rank_1:
>>>> > - if min(fact) < 1.0:
>>>> > + if min(fact) < 1.0:
>>>> > + if rank_1:
>>>> > back_calc[:] = array([r20_kex]*num_points)
>>>> > + return
>>>> > + else:
>>>> > + back_calc[:] = r20_kex
>>>> > return
>>>> >
>>>> > # Calculate R2eff.
>>>> > @@ -182,8 +195,10 @@
>>>> >
>>>> > # Catch errors, taking a sum over array is the fastest way to check
>>>> > for
>>>> > # +/- inf (infinity) and nan (not a number).
>>>> > - if rank_1:
>>>> > - if not isfinite(sum(R2eff)):
>>>> > + if not isfinite(sum(R2eff)):
>>>> > + if rank_1:
>>>> > R2eff = array([1e100]*num_points)
>>>> > + else:
>>>> > + R2eff = np.ones(R2eff.shape) * 1e100
>>>> >
>>>> > back_calc[:] = 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
>>>
>>>
_______________________________________________
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
