<div dir="ltr">Hi John,<div><br></div><div>Following up on Kacper's comment, which is what I was going to suggest :):<br><div><br></div><div>The 'data' passed to a field is actually a chunk of data, and the chunking method depends on how it was implemented in the frontend.  That said, the derived field should be somewhat agnostic to 'data's size/shape.  In your case, you get around this by saying you want k17rate to be np.zeros_like(T), which means to have the same shape as T.</div></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Fri, Sep 26, 2014 at 8:38 AM, John Regan <span dir="ltr"><<a href="mailto:johnanthonyregan@gmail.com" target="_blank">johnanthonyregan@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div><div><div><div><div><div><div><div><div>HI Britton, <br><br></div>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:<br><br>######################################################<br>import yt<br>from yt.utilities.physical_constants import kboltz<br>import numpy as np<br><br>def _k17(field, data):<br>    T = data["Temperature"]<br>    print "T.shape = ", T.shape<br>    k17rate = np.zeros(len(T))<br>    i = 0<br>    for t in T:<br>        if t < 1e4:<br>            k17rate[i] = 1.0e-8*pow(t, -0.4)<br>        else:<br>            k17rate[i] = 4.0e-4*pow(t, -1.4)*exp(-15100.0/t)<br>        i = i + 1<br>    return k17rate #cm^3 s^-1<br><br><br>yt.add_field("k17rate", units="cm**3/s", function=_k17)<br>Filename = "RD0000/RD0000"<br><br>ds = yt.load(Filename)<br>sp = ds.sphere('max', (100, 'kpc'))<br>print sp["k17rate"]<br>######################################################<br><br></div>This time the function is messy. There is probably a nice python way of doing the calculation <br></div>without looping - maybe something similar to k17rate[T < 1e4] = 1.0e-8*np.power(T, -0.4) ?<br><br></div>Either way this fails because the array that gets passed to k17 has a shape of (16, 16, 16) and so everything falls <br></div>apart. The _k17 function actually gets called 5 times in total - firstly with the array of shape (16, 16, 16) and then<br></div>with flat arrays but it fails on the first non-flat array anyway. <br><br></div>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).<br></div>Any solutions would of course also be appreciated!<br><br></div>Cheers,<br>John<br><br><br><div><div><div><br><div><div><div><div><div><div><br><br><br></div></div></div></div></div></div></div></div></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Fri, Sep 26, 2014 at 2:41 PM, John Regan <span dir="ltr"><<a href="mailto:johnanthonyregan@gmail.com" target="_blank">johnanthonyregan@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr">Yip that works! Thanks Britton!<br></div><div><div><div class="gmail_extra"><br><div class="gmail_quote">On Fri, Sep 26, 2014 at 2:35 PM, Britton Smith <span dir="ltr"><<a href="mailto:brittonsmith@gmail.com" target="_blank">brittonsmith@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr">Hi John,<div><br></div><div>I think you should be able to wipe out the units being used in the function and add your own by doing:</div><div><br></div><div>return data.ds.arr(krate.d, "cm**3/s")</div><div><br></div><div>krate.d gives you just the data array without units.</div><div><br></div><div>Britton</div></div><div class="gmail_extra"><br><div class="gmail_quote"><div><div>On Fri, Sep 26, 2014 at 7:05 AM, John Regan <span dir="ltr"><<a href="mailto:johnanthonyregan@gmail.com" target="_blank">johnanthonyregan@gmail.com</a>></span> wrote:<br></div></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><div><div dir="ltr"><div><div><div><div><div><div><div>Hi All, <br><br></div>I'm having a little trouble with Units in YT-3.0.<br><br></div>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!<br><br></div>An example script showing the problem is shown below:<br><br>#####################################################<br>import yt<br>from yt.utilities.physical_constants import kboltz<br>import numpy as np<br><br>def _k7(field, data): #H + e -> HM + gamma<br>    Tergs = data["Temperature"]*kboltz<br>    TeV = Tergs.convert_to_units('eV')<br>    krate = 6.77e-15*np.power(TeV, 0.8779)<br>    return krate #cm^3 s^-1<br><br>yt.add_field(("gas", "k7rate"), units="cm**3/s", function=_k7)<br><br>Filename = "....../RD0000/RD0000"<br><br>ds = yt.load(Filename)<br>sp = ds.sphere('max', (100, 'kpc'))<br>print sp["k7rate"]<br>#######################################################<br><br></div>The script fails with the error:<br>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)).<br><br></div>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?<br><br></div>Cheers,<br></div>John<br><br><br><div><div><div><div><br><div><div><br></div></div></div></div></div></div></div>
<br></div></div>_______________________________________________<br>
yt-users mailing list<br>
<a href="mailto:yt-users@lists.spacepope.org" target="_blank">yt-users@lists.spacepope.org</a><br>
<a href="http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org" target="_blank">http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org</a><br>
<br></blockquote></div><br></div>
<br>_______________________________________________<br>
yt-users mailing list<br>
<a href="mailto:yt-users@lists.spacepope.org" target="_blank">yt-users@lists.spacepope.org</a><br>
<a href="http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org" target="_blank">http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org</a><br>
<br></blockquote></div><br></div>
</div></div></blockquote></div><br></div>
<br>_______________________________________________<br>
yt-users mailing list<br>
<a href="mailto:yt-users@lists.spacepope.org">yt-users@lists.spacepope.org</a><br>
<a href="http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org" target="_blank">http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org</a><br>
<br></blockquote></div><br></div>