bug in ortho_ray and derived fields?
Hi I'm experiencing some weird behavior with ortho_ray() when used with derived fields. I suspect this is a bug, but I'm new to yt, so maybe I'm doing something incorrectly. This is for a 1D Enzo simulation, and the only fields in the simulation dataset are Density, x-velocity, and TotalEnergy. I defined two derived fields, InternalEnergy and Pressure: ### define InternalEnergy field def _InternalEnergy(field, data): return data['TotalEnergy'] - 0.5*data['x-velocity']*data['x-velocity'] add_field('InternalEnergy', function=_InternalEnergy, units=r'\rm{erg}/\rm{g}') ### define Pressure field def _Pressure(field, data): return (data.pf['Gamma'] - 1.0) * data['Density'] * data['InternalEnergy'] add_field('Pressure', function=_Pressure, units=r'\rm{dyne}/\rm{cm}^{2}') Next I want to extract a ray of Density, x-velocity, Pressure, and InternalEnergy along the only dimension. I figured ortho_ray is the right tool. ray = pf.h.ortho_ray(0, [0.5, 0.5], ['x-velocity', 'Density', 'InternalEnergy', 'Pressure']) print ray.keys() This produces: ['Pressure', 'x', 'x-velocity', 'InternalEnergy', 'Density'] Perfect, everything as desired. However, if I change the order of the input fields, some disappear from the ray. Examples: (1) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['Density', 'x-velocity', 'Pressure', 'InternalEnergy']) print ray.keys() ['x', 'x-velocity', 'Pressure', 'Density'] -- 'InternalEnergy' is missing. (2) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['Pressure', 'Density', 'x-velocity']) print ray.keys() ['Pressure'] -- 'Density' and 'x-velocity' are missing. (3) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['InternalEnergy', 'Density', 'x-velocity']) print ray.keys() ['InternalEnergy', 'Density'] -- 'x-velocity' is missing. The pattern I'm seeing here is that whenever 'InternalEnergy' is specified as an input only one field after it makes it into the ray. When 'Pressure' is specified, no further fields make it. I don't imagine this behavior is intended, but maybe it's due to some user error on my part? Thanks in advance for any clarifications. Mike -- ********************************************************************* * * * Dr. Michael Kuhlen Theoretical Astrophysics Center * * email: mqk@astro.berkeley.edu UC Berkeley * * cell phone: (831) 588-1468 601 Campbell Hall * * skype username: mikekuhlen Berkeley, CA 94720 * * * *********************************************************************
Hi Michael, This is indeed a bug, and unfortunately I don't think I'll be able to completely track it down until Wednesday or so. Until then, there are two other alternatives: 1. Don't specify the initial list of fields, simply access them -- typically they only need to be specified if the data object is expensive to read or expensive to generate fields for. 2. Use a 3D data object, which should also work on 1- and 2-D data. I'll be in touch when I figure it out. Best, Matt On Fri, Apr 16, 2010 at 4:28 PM, Michael Kuhlen <mqk@astro.berkeley.edu> wrote:
Hi
I'm experiencing some weird behavior with ortho_ray() when used with derived fields. I suspect this is a bug, but I'm new to yt, so maybe I'm doing something incorrectly.
This is for a 1D Enzo simulation, and the only fields in the simulation dataset are Density, x-velocity, and TotalEnergy. I defined two derived fields, InternalEnergy and Pressure:
### define InternalEnergy field def _InternalEnergy(field, data): return data['TotalEnergy'] - 0.5*data['x-velocity']*data['x-velocity']
add_field('InternalEnergy', function=_InternalEnergy, units=r'\rm{erg}/\rm{g}')
### define Pressure field def _Pressure(field, data): return (data.pf['Gamma'] - 1.0) * data['Density'] * data['InternalEnergy']
add_field('Pressure', function=_Pressure, units=r'\rm{dyne}/\rm{cm}^{2}')
Next I want to extract a ray of Density, x-velocity, Pressure, and InternalEnergy along the only dimension. I figured ortho_ray is the right tool.
ray = pf.h.ortho_ray(0, [0.5, 0.5], ['x-velocity', 'Density', 'InternalEnergy', 'Pressure']) print ray.keys()
This produces: ['Pressure', 'x', 'x-velocity', 'InternalEnergy', 'Density']
Perfect, everything as desired.
However, if I change the order of the input fields, some disappear from the ray. Examples:
(1) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['Density', 'x-velocity', 'Pressure', 'InternalEnergy']) print ray.keys()
['x', 'x-velocity', 'Pressure', 'Density'] -- 'InternalEnergy' is missing.
(2) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['Pressure', 'Density', 'x-velocity']) print ray.keys()
['Pressure'] -- 'Density' and 'x-velocity' are missing.
(3) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['InternalEnergy', 'Density', 'x-velocity']) print ray.keys()
['InternalEnergy', 'Density'] -- 'x-velocity' is missing.
The pattern I'm seeing here is that whenever 'InternalEnergy' is specified as an input only one field after it makes it into the ray. When 'Pressure' is specified, no further fields make it.
I don't imagine this behavior is intended, but maybe it's due to some user error on my part?
Thanks in advance for any clarifications.
Mike
-- ********************************************************************* * * * Dr. Michael Kuhlen Theoretical Astrophysics Center * * email: mqk@astro.berkeley.edu UC Berkeley * * cell phone: (831) 588-1468 601 Campbell Hall * * skype username: mikekuhlen Berkeley, CA 94720 * * * *********************************************************************
_______________________________________________ yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Hi Matt Thanks for your reply.
1. Don't specify the initial list of fields, simply access them -- typically they only need to be specified if the data object is expensive to read or expensive to generate fields for.
Ah, I didn't realize that's how it worked. But if I don't specify the fields, then they don't show up under ray.keys() at all... ray = pf.h.ortho_ray(0, [0.5, 0.5]) print ray.keys() ['x'] I see that all the fields, including the derived ones, are there and are accessible, but for some reason they don't show up as defined keys until I access them. Is that by design? If so, how do I iterate over the fields? Before I was doing something like this: for key in ray.keys(): do stuff with ray[key] But that won't work now.
2. Use a 3D data object, which should also work on 1- and 2-D data.
Actually I'm not explicitly setting what kind of data object to use. I simply did a 'pf = load(filename)' and it decided to go with EnzoHierarchy1D. How would I enforce a 3D data object, and why would that be better? Note that all of this is with the trunk (development) version. With the stable version, I couldn't get 1D simulation data to work at all. pf.h.find_max("Density") gave this error: File "/home/mqk/yt/src/yt-1.6/yt/lagos/BaseDataTypes.py", line 1664, in _get_data_from_grid return grid[field][pointI].ravel() IndexError: index (1) out of range (0<=index<0) in dimension 2 Mike On 04/16/2010 04:22 PM, Matthew Turk wrote:
Hi Michael,
This is indeed a bug, and unfortunately I don't think I'll be able to completely track it down until Wednesday or so. Until then, there are two other alternatives:
1. Don't specify the initial list of fields, simply access them -- typically they only need to be specified if the data object is expensive to read or expensive to generate fields for. 2. Use a 3D data object, which should also work on 1- and 2-D data.
I'll be in touch when I figure it out. Best,
Matt
On Fri, Apr 16, 2010 at 4:28 PM, Michael Kuhlen <mqk@astro.berkeley.edu> wrote:
Hi
I'm experiencing some weird behavior with ortho_ray() when used with derived fields. I suspect this is a bug, but I'm new to yt, so maybe I'm doing something incorrectly.
This is for a 1D Enzo simulation, and the only fields in the simulation dataset are Density, x-velocity, and TotalEnergy. I defined two derived fields, InternalEnergy and Pressure:
### define InternalEnergy field def _InternalEnergy(field, data): return data['TotalEnergy'] - 0.5*data['x-velocity']*data['x-velocity']
add_field('InternalEnergy', function=_InternalEnergy, units=r'\rm{erg}/\rm{g}')
### define Pressure field def _Pressure(field, data): return (data.pf['Gamma'] - 1.0) * data['Density'] * data['InternalEnergy']
add_field('Pressure', function=_Pressure, units=r'\rm{dyne}/\rm{cm}^{2}')
Next I want to extract a ray of Density, x-velocity, Pressure, and InternalEnergy along the only dimension. I figured ortho_ray is the right tool.
ray = pf.h.ortho_ray(0, [0.5, 0.5], ['x-velocity', 'Density', 'InternalEnergy', 'Pressure']) print ray.keys()
This produces: ['Pressure', 'x', 'x-velocity', 'InternalEnergy', 'Density']
Perfect, everything as desired.
However, if I change the order of the input fields, some disappear from the ray. Examples:
(1) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['Density', 'x-velocity', 'Pressure', 'InternalEnergy']) print ray.keys()
['x', 'x-velocity', 'Pressure', 'Density'] -- 'InternalEnergy' is missing.
(2) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['Pressure', 'Density', 'x-velocity']) print ray.keys()
['Pressure'] -- 'Density' and 'x-velocity' are missing.
(3) ray = pf.h.ortho_ray(0, [0.5, 0.5], ['InternalEnergy', 'Density', 'x-velocity']) print ray.keys()
['InternalEnergy', 'Density'] -- 'x-velocity' is missing.
The pattern I'm seeing here is that whenever 'InternalEnergy' is specified as an input only one field after it makes it into the ray. When 'Pressure' is specified, no further fields make it.
I don't imagine this behavior is intended, but maybe it's due to some user error on my part?
Thanks in advance for any clarifications.
Mike
-- ********************************************************************* * * * Dr. Michael Kuhlen Theoretical Astrophysics Center * * email: mqk@astro.berkeley.edu UC Berkeley * * cell phone: (831) 588-1468 601 Campbell Hall * * skype username: mikekuhlen Berkeley, CA 94720 * * * *********************************************************************
_______________________________________________ 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
-- ********************************************************************* * * * Dr. Michael Kuhlen Theoretical Astrophysics Center * * email: mqk@astro.berkeley.edu UC Berkeley * * cell phone: (831) 588-1468 601 Campbell Hall * * skype username: mikekuhlen Berkeley, CA 94720 * * * *********************************************************************
Hi Michael,
I see that all the fields, including the derived ones, are there and are accessible, but for some reason they don't show up as defined keys until I access them. Is that by design? If so, how do I iterate over the fields? Before I was doing something like this:
Yup, it's by design -- it loads on demand.
But that won't work now.
You can access the full list of fields as: pf.h.field_list # fields in the files pf.h.derived_field_list # fields it thinks it can make
Actually I'm not explicitly setting what kind of data object to use. I simply did a 'pf = load(filename)' and it decided to go with EnzoHierarchy1D. How would I enforce a 3D data object, and why would that be better?
Oh, I meant a 3D data object like a region, sphere, etc etc.
Note that all of this is with the trunk (development) version. With the stable version, I couldn't get 1D simulation data to work at all. pf.h.find_max("Density") gave this error:
File "/home/mqk/yt/src/yt-1.6/yt/lagos/BaseDataTypes.py", line 1664, in _get_data_from_grid return grid[field][pointI].ravel() IndexError: index (1) out of range (0<=index<0) in dimension 2
Ah, yes. I think I forgot to backport a fix. I will take a look at this shortly. Thanks for the heads up! -Matt
Hi Michael, I believe I have identified the problem, and it was twofold. One issue was that the path-parameter (in an ortho-ray, this would be x, y, or z, but in a non-ortho ray it would be t) was not being pre-pended to the list of fields being collected. The other issue was in how the derived fields work. Because often derived fields will make use of a number of other fields -- in your case this was density and x-velocity, for instance -- the derived field machinery checks to do two things. The first is that if at all possible, generate the field inside the data object, rather than generating it inside the grids and then re-selecting the cells that go to the data object. The other is that any fields that are created in the intermediate steps inside the data object are deleted. This causes a bit of a problem if two different nested iterators iterate over the same list; in this case, it was the fields_to_get list, which was being changed by one iterator at a higher stack level and still being iterated over at a lower stack level. I've changed this in r1694 in trunk. Tonight I'll see what I can do about backporting things to 1.6, but I believe that we're almost at our set of features for 1.7, so it might not happen. Thanks very much for the heads up on this! -Matt On Fri, Apr 16, 2010 at 4:47 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi Michael,
I see that all the fields, including the derived ones, are there and are accessible, but for some reason they don't show up as defined keys until I access them. Is that by design? If so, how do I iterate over the fields? Before I was doing something like this:
Yup, it's by design -- it loads on demand.
But that won't work now.
You can access the full list of fields as:
pf.h.field_list # fields in the files pf.h.derived_field_list # fields it thinks it can make
Actually I'm not explicitly setting what kind of data object to use. I simply did a 'pf = load(filename)' and it decided to go with EnzoHierarchy1D. How would I enforce a 3D data object, and why would that be better?
Oh, I meant a 3D data object like a region, sphere, etc etc.
Note that all of this is with the trunk (development) version. With the stable version, I couldn't get 1D simulation data to work at all. pf.h.find_max("Density") gave this error:
File "/home/mqk/yt/src/yt-1.6/yt/lagos/BaseDataTypes.py", line 1664, in _get_data_from_grid return grid[field][pointI].ravel() IndexError: index (1) out of range (0<=index<0) in dimension 2
Ah, yes. I think I forgot to backport a fix. I will take a look at this shortly. Thanks for the heads up!
-Matt
participants (2)
-
Matthew Turk -
Michael Kuhlen