[yt-svn] commit/yt: xarthisius: Merged in ngoldbaum/yt (pull request #2213)
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Wed Jun 29 11:06:47 PDT 2016
1 new commit in yt:
https://bitbucket.org/yt_analysis/yt/commits/9b2d8e598c95/
Changeset: 9b2d8e598c95
Branch: yt
User: xarthisius
Date: 2016-06-29 18:06:34+00:00
Summary: Merged in ngoldbaum/yt (pull request #2213)
MKS code unit system fixes. Closes #1229
Affected #: 11 files
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d tests/tests.yaml
--- a/tests/tests.yaml
+++ b/tests/tests.yaml
@@ -2,7 +2,7 @@
local_artio_000:
- yt/frontends/artio/tests/test_outputs.py
- local_athena_000:
+ local_athena_001:
- yt/frontends/athena
local_chombo_000:
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/data_objects/static_output.py
--- a/yt/data_objects/static_output.py
+++ b/yt/data_objects/static_output.py
@@ -43,6 +43,7 @@
ParameterFileStore, \
NoParameterShelf, \
output_type_registry
+from yt.units.dimensions import current_mks
from yt.units.unit_object import Unit, unit_system_registry
from yt.units.unit_registry import UnitRegistry
from yt.fields.derived_field import \
@@ -237,15 +238,19 @@
self.no_cgs_equiv_length = False
self._create_unit_registry()
- self._parse_parameter_file()
- self.set_units()
- self._setup_coordinate_handler()
create_code_unit_system(self)
if unit_system == "code":
unit_system = str(self)
+ else:
+ unit_system = str(unit_system).lower()
+
self.unit_system = unit_system_registry[unit_system]
+ self._parse_parameter_file()
+ self.set_units()
+ self._setup_coordinate_handler()
+
# Because we need an instantiated class to check the ds's existence in
# the cache, we move that check to here from __new__. This avoids
# double-instantiation.
@@ -865,18 +870,30 @@
self._set_code_unit_attributes()
# here we override units, if overrides have been provided.
self._override_code_units()
+
self.unit_registry.modify("code_length", self.length_unit)
self.unit_registry.modify("code_mass", self.mass_unit)
self.unit_registry.modify("code_time", self.time_unit)
if hasattr(self, 'magnetic_unit'):
- # If we do not have this set, but some fields come in in
- # "code_magnetic", this will allow them to remain in that unit.
+ # if the magnetic unit is in T, we need to recreate the code unit
+ # system as an MKS-like system
+ if current_mks in self.magnetic_unit.units.dimensions.free_symbols:
+
+ if self.unit_system == unit_system_registry[str(self)]:
+ unit_system_registry.pop(str(self))
+ create_code_unit_system(self, current_mks_unit='A')
+ self.unit_system = unit_system_registry[str(self)]
+ elif str(self.unit_system) == 'mks':
+ pass
+ else:
+ self.magnetic_unit = \
+ self.magnetic_unit.to_equivalent('gauss', 'CGS')
self.unit_registry.modify("code_magnetic", self.magnetic_unit)
vel_unit = getattr(
self, "velocity_unit", self.length_unit / self.time_unit)
pressure_unit = getattr(
self, "pressure_unit",
- self.mass_unit / (self.length_unit * self.time_unit)**2)
+ self.mass_unit / (self.length_unit * (self.time_unit)**2))
temperature_unit = getattr(self, "temperature_unit", 1.0)
density_unit = getattr(self, "density_unit", self.mass_unit / self.length_unit**3)
self.unit_registry.modify("code_velocity", vel_unit)
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/frontends/athena/data_structures.py
--- a/yt/frontends/athena/data_structures.py
+++ b/yt/frontends/athena/data_structures.py
@@ -327,11 +327,6 @@
self.dataset.domain_dimensions = \
np.round(self.dataset.domain_width/gdds[0]).astype('int')
- # Need to reset the units in the dataset based on the correct
- # domain left/right/dimensions.
- # DEV: Is this really necessary?
- #self.dataset._set_code_unit_attributes()
-
if self.dataset.dimensionality <= 2 :
self.dataset.domain_dimensions[2] = np.int(1)
if self.dataset.dimensionality == 1 :
@@ -483,19 +478,10 @@
# We set these to cgs for now, but they may be overridden later.
mylog.warning("Assuming 1.0 = 1.0 %s", cgs)
setattr(self, "%s_unit" % unit, self.quan(1.0, cgs))
-
- def set_code_units(self):
- super(AthenaDataset, self).set_code_units()
- mag_unit = getattr(self, "magnetic_unit", None)
- vel_unit = getattr(self, "velocity_unit", None)
- if mag_unit is None:
- self.magnetic_unit = np.sqrt(4*np.pi * self.mass_unit /
- (self.time_unit**2 * self.length_unit))
+ self.magnetic_unit = np.sqrt(4*np.pi * self.mass_unit /
+ (self.time_unit**2 * self.length_unit))
self.magnetic_unit.convert_to_units("gauss")
- self.unit_registry.modify("code_magnetic", self.magnetic_unit)
- if vel_unit is None:
- self.velocity_unit = self.length_unit/self.time_unit
- self.unit_registry.modify("code_velocity", self.velocity_unit)
+ self.velocity_unit = self.length_unit / self.time_unit
def _parse_parameter_file(self):
self._handle = open(self.parameter_filename, "rb")
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/frontends/athena/fields.py
--- a/yt/frontends/athena/fields.py
+++ b/yt/frontends/athena/fields.py
@@ -19,7 +19,7 @@
kboltz, mh
b_units = "code_magnetic"
-pres_units = "code_mass/(code_length*code_time**2)"
+pres_units = "code_pressure"
erg_units = "code_mass * (code_length/code_time)**2"
rho_units = "code_mass / code_length**3"
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/frontends/flash/data_structures.py
--- a/yt/frontends/flash/data_structures.py
+++ b/yt/frontends/flash/data_structures.py
@@ -253,12 +253,9 @@
self.time_unit = self.quan(1.0, "s")
self.velocity_unit = self.quan(1.0, "cm/s")
self.temperature_unit = self.quan(temperature_factor, "K")
- self.unit_registry.modify("code_magnetic", self.magnetic_unit)
-
+
def set_code_units(self):
super(FLASHDataset, self).set_code_units()
- self.unit_registry.modify("code_temperature",
- self.temperature_unit.value)
def _find_parameter(self, ptype, pname, scalar = False):
nn = "/%s %s" % (ptype,
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/frontends/gdf/data_structures.py
--- a/yt/frontends/gdf/data_structures.py
+++ b/yt/frontends/gdf/data_structures.py
@@ -26,6 +26,10 @@
GridIndex
from yt.data_objects.static_output import \
Dataset
+from yt.units.dimensions import \
+ dimensionless as sympy_one
+from yt.units.unit_object import \
+ Unit
from yt.utilities.exceptions import \
YTGDFUnknownGeometry
from yt.utilities.lib.misc_utilities import \
@@ -211,12 +215,22 @@
current_unit = h5f["/dataset_units/%s" % unit_name]
value = current_unit.value
unit = current_unit.attrs["unit"]
+ # need to convert to a Unit object and check dimensions
+ # because unit can be things like
+ # 'dimensionless/dimensionless**3' so naive string
+ # comparisons are insufficient
+ unit = Unit(unit, registry=self.unit_registry)
+ if unit_name.endswith('_unit') and unit.dimensions is sympy_one:
+ un = unit_name[:-5]
+ un = un.replace('magnetic', 'magnetic_field', 1)
+ unit = self.unit_system[un]
+ setattr(self, unit_name, self.quan(value, unit))
setattr(self, unit_name, self.quan(value, unit))
if unit_name in h5f["/field_types"]:
if unit_name in self.field_units:
mylog.warning("'field_units' was overridden by 'dataset_units/%s'"
% (unit_name))
- self.field_units[unit_name] = unit.decode('utf8')
+ self.field_units[unit_name] = str(unit)
else:
self.length_unit = self.quan(1.0, "cm")
self.mass_unit = self.quan(1.0, "g")
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/units/tests/test_unit_systems.py
--- a/yt/units/tests/test_unit_systems.py
+++ b/yt/units/tests/test_unit_systems.py
@@ -15,7 +15,7 @@
from yt.units.unit_systems import UnitSystem
from yt.units import dimensions
from yt.convenience import load
-from yt.testing import assert_almost_equal, requires_file
+from yt.testing import assert_almost_equal, assert_allclose, requires_file
from yt.config import ytcfg
def test_unit_systems():
@@ -46,7 +46,7 @@
"velocity_divergence": "1/Myr",
"density_gradient_x": "Msun/kpc**4"}
test_units["code"] = {"density": "code_mass/code_length**3",
- "kinetic_energy": "code_mass/(code_length*code_time**2)",
+ "kinetic_energy": "code_pressure",
"velocity_magnitude": "code_velocity",
"velocity_divergence": "code_velocity/code_length",
"density_gradient_x": "code_mass/code_length**4"}
@@ -96,4 +96,44 @@
v1 = dd_cgs[field].in_base(us)
v2 = dd[field]
assert_almost_equal(v1.v, v2.v)
- assert str(v2.units) == test_units[us][field]
\ No newline at end of file
+ assert str(v2.units) == test_units[us][field]
+
+wdm = 'WDMerger_hdf5_chk_1000/WDMerger_hdf5_chk_1000.hdf5'
+ at requires_file(wdm)
+def test_tesla_magnetic_unit():
+ ytcfg["yt", "skip_dataset_cache"] = "True"
+ for us in ['cgs', 'mks', 'code']:
+ ds = load(wdm, unit_system=us,
+ units_override={'magnetic_unit': (1.0, 'T')})
+ ad = ds.all_data()
+ dens = ad['density']
+ magx = ad['magx']
+ magnetic_field_x = ad['magnetic_field_x']
+
+ if us == 'cgs':
+ assert str(dens.units) == 'g/cm**3'
+ assert str(magx.units) == 'code_magnetic'
+ assert magx.uq == ds.quan(1e4, 'G')
+ assert str(magnetic_field_x.units) == 'gauss'
+ assert_allclose(magx.value, magnetic_field_x.value/1e4)
+ assert_allclose(
+ magnetic_field_x.to_equivalent('T', 'SI').value,
+ magnetic_field_x.value/1e4)
+
+ if us == 'mks':
+ assert str(dens.units) == 'kg/m**3'
+ assert str(magx.units) == 'code_magnetic'
+ assert magx.uq == ds.quan(1, 'T')
+ assert str(magnetic_field_x.units) == 'T'
+ assert_allclose(magx.value, magnetic_field_x.value)
+ assert_allclose(magnetic_field_x.to_equivalent('G', 'CGS').value,
+ magnetic_field_x.value*1e4)
+
+ if us == 'code':
+ assert str(dens.units) == 'code_mass/code_length**3'
+ assert str(magx.units) == 'code_magnetic'
+ assert magx.uq == ds.quan(1, 'T')
+ assert str(magnetic_field_x.units) == 'code_magnetic'
+ assert_allclose(magx.value, magnetic_field_x.value)
+ assert_allclose(magnetic_field_x.to_equivalent('G', 'CGS').value,
+ magnetic_field_x.value*1e4)
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/units/unit_object.py
--- a/yt/units/unit_object.py
+++ b/yt/units/unit_object.py
@@ -26,6 +26,8 @@
base_dimensions, temperature, \
dimensionless, current_mks, \
angle
+from yt.units.equivalencies import \
+ equivalence_registry
from yt.units.unit_lookup_table import \
unit_prefixes, prefixable_units, latex_prefixes, \
default_base_units
@@ -432,6 +434,25 @@
return False
return True
+ def list_equivalencies(self):
+ """
+ Lists the possible equivalencies associated with this unit object
+ """
+ for k, v in equivalence_registry.items():
+ if self.has_equivalent(k):
+ print(v())
+
+ def has_equivalent(self, equiv):
+ """
+ Check to see if this unit object as an equivalent unit in *equiv*.
+ """
+ try:
+ this_equiv = equivalence_registry[equiv]()
+ except KeyError:
+ raise KeyError("No such equivalence \"%s\"." % equiv)
+ old_dims = self.dimensions
+ return old_dims in this_equiv.dims
+
def get_base_equivalent(self, unit_system="cgs"):
"""
Create and return dimensionally-equivalent units in a specified base.
@@ -652,6 +673,8 @@
my_dims = dimensions.expand()
if my_dims is dimensionless:
return ""
+ if my_dims in base_units:
+ return base_units[my_dims]
for factor in my_dims.as_ordered_factors():
dim = list(factor.free_symbols)[0]
unit_string = str(base_units[dim])
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/units/unit_registry.py
--- a/yt/units/unit_registry.py
+++ b/yt/units/unit_registry.py
@@ -93,11 +93,13 @@
"in this registry." % symbol)
del self.lut[symbol]
+ if symbol in self.unit_objs:
+ del self.unit_objs[symbol]
def modify(self, symbol, base_value):
"""
- Change the base value of a dimension. Useful for adjusting code units
- after parsing parameters."
+ Change the base value of a unit symbol. Useful for adjusting code units
+ after parsing parameters.
"""
if symbol not in self.lut:
@@ -106,9 +108,17 @@
"in this registry." % symbol)
if hasattr(base_value, "in_base"):
- base_value = base_value.in_base().value
+ new_dimensions = base_value.units.dimensions
+ base_value = base_value.in_base('cgs-ampere')
+ base_value = base_value.value
+ else:
+ new_dimensions = self.lut[symbol][1]
- self.lut[symbol] = (float(base_value), ) + self.lut[symbol][1:]
+ self.lut[symbol] = ((float(base_value), new_dimensions) +
+ self.lut[symbol][2:])
+
+ if symbol in self.unit_objs:
+ del self.unit_objs[symbol]
def keys(self):
"""
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/units/unit_systems.py
--- a/yt/units/unit_systems.py
+++ b/yt/units/unit_systems.py
@@ -79,10 +79,9 @@
def __getitem__(self, key):
if isinstance(key, string_types):
- if key not in self._dims:
- self._dims.append(key)
key = getattr(dimensions, key)
- if key not in self.units_map:
+ um = self.units_map
+ if key not in um or um[key].dimensions is not key:
units = _get_system_unit_string(key, self.units_map)
self.units_map[key] = Unit(units, registry=self.registry)
return self.units_map[key]
@@ -109,10 +108,16 @@
repr += " %s: %s\n" % (key, self.units_map[dim])
return repr
-def create_code_unit_system(ds):
- code_unit_system = UnitSystem(str(ds), "code_length", "code_mass", "code_time",
- "code_temperature", registry=ds.unit_registry)
+def create_code_unit_system(ds, current_mks_unit=None):
+ code_unit_system = UnitSystem(
+ str(ds), "code_length", "code_mass", "code_time", "code_temperature",
+ current_mks_unit=current_mks_unit, registry=ds.unit_registry)
code_unit_system["velocity"] = "code_velocity"
+ if current_mks_unit:
+ code_unit_system["magnetic_field_mks"] = "code_magnetic"
+ else:
+ code_unit_system["magnetic_field_cgs"] = "code_magnetic"
+ code_unit_system["pressure"] = "code_pressure"
cgs_unit_system = UnitSystem("cgs", "cm", "g", "s")
cgs_unit_system["energy"] = "erg"
@@ -147,3 +152,14 @@
planck_unit_system = UnitSystem("planck", "l_pl", "m_pl", "t_pl", temperature_unit="T_pl")
planck_unit_system["energy"] = "E_pl"
planck_unit_system["charge_cgs"] = "q_pl"
+
+
+cgs_ampere_unit_system = UnitSystem('cgs-ampere', 'cm', 'g', 's',
+ current_mks_unit='A')
+cgs_ampere_unit_system["energy"] = "erg"
+cgs_ampere_unit_system["specific_energy"] = "erg/g"
+cgs_ampere_unit_system["pressure"] = "dyne/cm**2"
+cgs_ampere_unit_system["force"] = "dyne"
+cgs_ampere_unit_system["magnetic_field_cgs"] = "gauss"
+cgs_ampere_unit_system["charge_cgs"] = "esu"
+cgs_ampere_unit_system["current_cgs"] = "statA"
diff -r 3215234b6e95186ebf3a1858937e0e476f2b260f -r 9b2d8e598c9587b5508dc0c9e5530c32fccc329d yt/units/yt_array.py
--- a/yt/units/yt_array.py
+++ b/yt/units/yt_array.py
@@ -596,18 +596,21 @@
The unit that you wish to convert to.
equiv : string
The equivalence you wish to use. To see which equivalencies are
- supported for this unitful quantity, try the :meth:`list_equivalencies`
- method.
+ supported for this unitful quantity, try the
+ :meth:`list_equivalencies` method.
Examples
--------
>>> a = yt.YTArray(1.0e7,"K")
>>> a.to_equivalent("keV", "thermal")
"""
- unit_quan = YTQuantity(1.0, unit, registry=self.units.registry)
+ conv_unit = Unit(unit, registry=self.units.registry)
this_equiv = equivalence_registry[equiv]()
- if self.has_equivalent(equiv) and (unit_quan.has_equivalent(equiv) or this_equiv._one_way):
- new_arr = this_equiv.convert(self, unit_quan.units.dimensions, **kwargs)
+ oneway_or_equivalent = (
+ conv_unit.has_equivalent(equiv) or this_equiv._one_way)
+ if self.has_equivalent(equiv) and oneway_or_equivalent:
+ new_arr = this_equiv.convert(
+ self, conv_unit.dimensions, **kwargs)
if isinstance(new_arr, tuple):
try:
return YTArray(new_arr[0], new_arr[1]).in_units(unit)
@@ -623,21 +626,14 @@
Lists the possible equivalencies associated with this YTArray or
YTQuantity.
"""
- for k,v in equivalence_registry.items():
- if self.has_equivalent(k):
- print(v())
+ self.units.list_equivalencies()
def has_equivalent(self, equiv):
"""
Check to see if this YTArray or YTQuantity has an equivalent unit in
*equiv*.
"""
- try:
- this_equiv = equivalence_registry[equiv]()
- except KeyError:
- raise KeyError("No such equivalence \"%s\"." % equiv)
- old_dims = self.units.dimensions
- return old_dims in this_equiv.dims
+ return self.units.has_equivalent(equiv)
def ndarray_view(self):
"""
Repository URL: https://bitbucket.org/yt_analysis/yt/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the yt-svn
mailing list