Hello! I'm really sorry for the barrage of questions; hopefully this will be the last one!
I now get the error:  operands could not be broadcast together with shapes (2072685,) (2153824,)
on line: fHI = (recomb * ne) / (gamma_HI + (gamma_c * ne))

I printed the shapes of all of those arrays, and this is what I get:
nH.shape
(16, 16, 16)
ne.shape
(16, 16, 16)
gamma_c.shape
(1,)
recomb.shape
(1,)
nH.shape
(2153824,)
ne.shape
(2153824,)
gamma_c.shape
(2072685,)
recomb.shape
(2072685,)

Firstly, I'm a little confused as to why their dimensions change, but what really matters are the bottom ones.  I know how many gas particles are in my gadget simulation: 2072685.  All of the arrays were created from either data["Gas", "InternalEnergy"] or data["gas", "density"], which should, I would imagine, have the same dimensions: the number of gas particles.  I was just wondering if anyone had seen something like this before? The function in question was posted to a pastebin in an earlier email on this thread, and I have since made all of the arrays YTArrays, as suggested above. Thanks!
-Jared

On Tue, Aug 25, 2015 at 5:58 PM, Jared Coughlin wrote:
Thanks, I'll try that!
-Jared

On Tue, Aug 25, 2015 at 5:55 PM, Andrew Myers wrote:
Hi Jared,

You probably don't want to loop over the entire array in pure python; that could be very slow for large datasets. One thing you can do is use slicing to select the elements of the array that meet your criterion:

indices = data['gas','InternalEnergy'] <= 1e4
data['gas', 'InternalEnergy'][indices] = *something*

that would set all the elements where the InternalEnergy field is less than or equal to 1e4 to *something*. Usually you can express most operations you'd want to do this way without explicitly looping over the elements.

-Andrew

On Tue, Aug 25, 2015 at 2:44 PM, Jared Coughlin wrote:
That seemed to have worked, thanks!

I just have another quick question, if you wouldn't mind? In my function to generate my new field, there's an if statement that looks at the value of the temperature: if data['gas','InternalEnergy'] <= 1e4:
However, I get the error: "The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()." I wasn't sure how to loop over the elements of the array to access a specific element, or if that's even possible? Using a.any() or a.all() don't seem to fit the bill for this particular problem, either. Thanks!
-Jared

On Tue, Aug 25, 2015 at 3:19 PM, Nathan Goldbaum wrote:

On Tue, Aug 25, 2015 at 2:17 PM, Jared Coughlin wrote:
Great, thanks! I tried that, and it gave me the error: "YTQuantity instances must be scalars." I'm not sure what to do about that, since all of these quantities depend on either temperature or density, both of which are arrays. I tried looking at the doc page for YTQuantity, but I didn't see anything that looked like it could solve this problem, though it's possible I missed it.

You need to use YTArray instead. YTQuantity is a subclass of YTArray that only holds scalars. YTArray can hold arrays of data with units.

-Jared

On Mon, Aug 24, 2015 at 8:13 PM, Nathan Goldbaum wrote:

On Monday, August 24, 2015, Cameron Hummels <chummels@gmail.com> wrote:
Hi Jared,

The problem is one of units.  In order to make common arithmetic work like addition, yt requires both arguments going into the addition operator to have the same units.

fHI = (recomb * ne) / (gamma_HI + (gamma_c *ne))

you're adding two things (gamma_HI) and (gamma_c*ne) that have different units--in this case, something with (code_mass/code_length**3) units and something with no units defined as (1).

Just a quick clarification: they need not have the same units, just the same dimensions.  You can add two quantities with units of grams and solar masses (for example) but not grams and liters.

In order to fix this, you can re-assign the appropriate units to different arrays or quantities with the YTArray and YTQuantity classes.  In this case, make sure that the two arguments going into the addition have the same units.  If you want to recast "ne" to have number density units, you can do this:  from yt.units.yt_array import YTQuantity; ne = YTQuantity(ne, 'cm**-3').  For more info on units, check out http://yt-project.org/docs/dev/analyzing/units/index.html

I hope this helps!

Cameron

On Mon, Aug 24, 2015 at 4:12 PM, Jared Coughlin wrote:
Hello! I have a gadget snapshot for which the standard internal energy block has been replaced by one containing the temperature.  I'm trying to calculate a derived field for the neutral hydrogen number density,  but I'm getting an error when I try to do ds.add_field():

However, I get an error about being unable to add quantities with differing units:

The docs say not to do any unit conversion because that is apparently taken care of behind the scenes, so I didn't.  I was just wondering if there was a way around this? Thanks!
-Jared

_______________________________________________
yt-users mailing list
yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org

--
Cameron Hummels
NSF Postdoctoral Fellow
Department of Astronomy
California Institute of Technology

_______________________________________________
yt-users mailing list
yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org

_______________________________________________
yt-users mailing list
yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org

_______________________________________________
yt-users mailing list
yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org

_______________________________________________
yt-users mailing list
yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org

_______________________________________________
yt-users mailing list
yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org