[yt-users] yt-3.0 derived field

Nathan Goldbaum nathan12343 at gmail.com
Sun Jan 4 19:49:44 PST 2015 

Hi Junhwan,



On Sun, Jan 4, 2015 at 7:20 PM, Junhwan Choi (최준환) <choi.junhwan at gmail.com>
wrote:

> Hi yt user,
>
> I start to use yt-3.x version for the gadget analyses.
> Although I have previously used the yt-2.6x version for the enzo
> analyses, I find that there are many differences between two.
> [Now, I find that new yt is significantly improved than previous one
> but I need some adjustment..]
> Here, I have a question regarding to derived field.
> Since, my gadegt outputs only have InterenalEnergy, I need to make a
> derived field for temperature as follow:


There are a couple things wrong here.  First, is this data in the Gadget
binary format or is it in the HDF5 format?  If it's the former, the "Gas"
particle type should have all of the fields you're referencing.  That means
that where you have data["gas","ElectronAbundance"], you should instead
have data["Gas", "ElectronAbundance"] (note the change from "g" to "G") .
If this is Gadget HDF5 data, instead of "Gas", you should be referencing
"PartType0".

The "gas" field type represents fluid quantities defined on a space-filling
mesh.  The "Gas" (or "PartType0") field type represents quantities only
defined at the location of SPH particles. You should never mix particle and
mesh fields, since they're not defined at the same locations in space.  To
get around this, yt provides a number of internal facilities for creating a
mesh field from a particle field using a smoothing or deposition
operation.  For example, the ("gas", "density") field is generated by yt
using an SPH smoothing operation for SPH datasets.  We will need to do
something similar here in the end to create a ("gas", "temperature") field
from the internal ("Gas", "InternalEnergy") field.


>

==============
> @derived_field(name=("gas","temperature"), units="K")
> def _temperature(field, data):
>     BOLTZMANN   = 1.3806e-16
>     PROTONMASS  = 1.6726e-24
>     GAMMA       = 5.0 / 3.0
>     H_MASSFRAC = 0.76
>     mu = 4.0 / (3.0 * H_MASSFRAC + 1.0 + 4.0 * H_MASSFRAC *
> data["gas","ElectronAbundance"])
>     return data["gas","InternalEnergy"]*(GAMMA-1)*mu*PROTONMASS/BOLTZMANN
>

Next, In order to make the units work out correctly, you're going to need
to ensure that all of the quantities in the expression have correct units.

There are two things we'll have to fix here in order for that to work.

First, it seems the InternalEnergy field is being read in with
dimensionless units.  I don't think that's correct and it's something we'll
need to fix in yt.

Looking at the gadget frontend, it seems the InternalEnergy field is indeed
hard-coded to have dimensionless units (see line 32 in
yt/frontends/sph/fields.py).  I think this is a bug.  Do you happen to know
what units Gadget writes the InternalEnergy field out to disk in?  That
section of the SPH frontend should be rewritten to use the correct units.

Next, the hard-coded constants you're using should be quantities with units
rather than floats.  In addition, for consistency with the rest of yt you
should use the physical constant definitions provided in the
yt.utilities.physical_constants module.  If you import from here, you'll
get quantities with units for free.

Last, as I alluded to above, once you create a derived Temperature field
defined for all of the SPH particles, you will also need to create the
("gas", "temperature") field via an SPH smoothing operation. There's an
example of how to do this in the OWLS frontend, inside the
setup_particle_fields function in yt/frontends/owls/fields.py.  In fact,
some of the logic in the OWLS frontend that sets up extra smoothed fields
could probably be better suited for living in the generic gadget frontend,
which the OWLS frontend subclasses.

In the end, I think wanting a proper temperature field set up for all
gadget datasets is probably something that we should have set up in the
gadget frontend itself, obviating the need for future users in your
footsteps to need to define a derived field.  If you'd like, you could file
a bug about this issue and one of us might take a crack at setting this up
once and for all. Doing this yourself would probably also be a good way to
get more familiarity with the internals of the gadget frontend.

Sorry that I don't have easy solutions for you here.

-Nathan



> ================
>
> and try to make phase diagram and got the following error message:
> .....
> yt : [INFO     ] 2015-01-04 19:56:35,791 Max Value is 4.03205e-21 at
> 30150.5327224731445312 31066.5464401245117188 29380.5456161499023438
> Traceback (most recent call last):
>   File "phase.py", line 32, in <module>
>     weight_field=None)
>   File "/Users/jhchoi/common/src/yt/yt/visualization/profile_plotter.py",
> line 699, in __init__
>     fractional=fractional)
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/profiles.py", line
> 1331, in create_profile
>     for f, l in zip(bin_fields, logs)]
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/derived_quantities.py",
> line 482, in __call__
>     rv = super(Extrema, self).__call__(fields, non_zero)
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/derived_quantities.py",
> line 56, in __call__
>     sto.result = self.process_chunk(ds, *args, **kwargs)
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/derived_quantities.py",
> line 490, in process_chunk
>     fd = data[field]
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/data_containers.py",
> line 240, in __getitem__
>     self.get_data(f)
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/data_containers.py",
> line 661, in get_data
>     self._generate_fields(fields_to_generate)
>   File "/Users/jhchoi/common/src/yt/yt/data_objects/data_containers.py",
> line 687, in _generate_fields
>     fd.convert_to_units(fi.units)
>   File "/Users/jhchoi/common/src/yt/yt/units/yt_array.py", line 416,
> in convert_to_units
>     new_units = self._unit_repr_check_same(units)
>   File "/Users/jhchoi/common/src/yt/yt/units/yt_array.py", line 402,
> in _unit_repr_check_same
>     self.units, self.units.dimensions, units, units.dimensions)
> yt.utilities.exceptions.YTUnitConversionError: Unit dimensionalities
> do not match. Tried to convert between dimensionless (dim 1) and K
> (dim (temperature)).
>
> It looks that new yt become a bit picky on the units.
> In this derived field, the unit contribution from BOLTZMANN and
> PROTONMASS is not recognized and it causes error, I think.
> How can I resolve this error?
> Moreover, how can I make derived field without sensitive unit matching?
> Sometimes, I used derived field when I need to convert log units for
> visualization, and when dealing with the dimensionless field.
>
> Thank you for help,
> Junhwan
> _______________________________________________
> 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/20150104/38d4968a/attachment.html>

More information about the yt-users mailing list