[yt-users] Setting Units

[John Regan]

HI Britton,

I've pretty much got everything working now as it should but there one
thing that is throwing me. In the above script k7 was fairly straight
forward but suppose I want to calculate the rates for something like:

######################################################
import yt
from yt.utilities.physical_constants import kboltz
import numpy as np

def _k17(field, data):
    T = data["Temperature"]
    print "T.shape = ", T.shape
    k17rate = np.zeros(len(T))
    i = 0
    for t in T:
        if t < 1e4:
            k17rate[i] = 1.0e-8*pow(t, -0.4)
        else:
            k17rate[i] = 4.0e-4*pow(t, -1.4)*exp(-15100.0/t)
        i = i + 1
    return k17rate #cm^3 s^-1


yt.add_field("k17rate", units="cm**3/s", function=_k17)
Filename = "RD0000/RD0000"

ds = yt.load(Filename)
sp = ds.sphere('max', (100, 'kpc'))
print sp["k17rate"]
######################################################

This time the function is messy. There is probably a nice python way of
doing the calculation
without looping - maybe something similar to k17rate[T < 1e4] =
1.0e-8*np.power(T, -0.4) ?

Either way this fails because the array that gets passed to k17 has a shape
of (16, 16, 16) and so everything falls
apart. The _k17 function actually gets called 5 times in total - firstly
with the array of shape (16, 16, 16) and then
with flat arrays but it fails on the first non-flat array anyway.

Just wondering if you know what's going on here (why the function gets
called 5 times and why the first array in a 3D array).
Any solutions would of course also be appreciated!

Cheers,
John







On Fri, Sep 26, 2014 at 2:41 PM, John Regan <johnanthonyregan at gmail.com>
wrote:

> Yip that works! Thanks Britton!
>
> On Fri, Sep 26, 2014 at 2:35 PM, Britton Smith <brittonsmith at gmail.com>
> wrote:
>
>> Hi John,
>>
>> I think you should be able to wipe out the units being used in the
>> function and add your own by doing:
>>
>> return data.ds.arr(krate.d, "cm**3/s")
>>
>> krate.d gives you just the data array without units.
>>
>> Britton
>>
>> On Fri, Sep 26, 2014 at 7:05 AM, John Regan <johnanthonyregan at gmail.com>
>> wrote:
>>
>>> Hi All,
>>>
>>> I'm having a little trouble with Units in YT-3.0.
>>>
>>> I want to calculate the rates for some collisional dissociations as used
>>> in Enzo. The rates are fits so the dimensions are not going to be as
>>> required and YT doesn't like that!
>>>
>>> An example script showing the problem is shown below:
>>>
>>> #####################################################
>>> import yt
>>> from yt.utilities.physical_constants import kboltz
>>> import numpy as np
>>>
>>> def _k7(field, data): #H + e -> HM + gamma
>>>     Tergs = data["Temperature"]*kboltz
>>>     TeV = Tergs.convert_to_units('eV')
>>>     krate = 6.77e-15*np.power(TeV, 0.8779)
>>>     return krate #cm^3 s^-1
>>>
>>> yt.add_field(("gas", "k7rate"), units="cm**3/s", function=_k7)
>>>
>>> Filename = "....../RD0000/RD0000"
>>>
>>> ds = yt.load(Filename)
>>> sp = ds.sphere('max', (100, 'kpc'))
>>> print sp["k7rate"]
>>> #######################################################
>>>
>>> The script fails with the error:
>>> yt.utilities.exceptions.YTUnitConversionError: Unit dimensionalities do
>>> not match. Tried to convert between eV**(8779/10000) (dim
>>> (length)**(8779/5000)*(mass)**(8779/10000)/(time)**(8779/5000)) and cm**3/s
>>> (dim (length)**3/(time)).
>>>
>>> which makes sense because the _k7() function is just a fit based on
>>> temperature. Is there a way to "set" the units of k7 to be cm^3 s^-1?
>>>
>>> Cheers,
>>> John
>>>
>>>
>>>
>>>
>>>
>>> _______________________________________________
>>> yt-users mailing list
>>> yt-users at lists.spacepope.org
>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>
>>>
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20140926/41dcb364/attachment.htm>