[yt-svn] commit/yt: MatthewTurk: Merged in jzuhone/yt-3.x (pull request #1286)
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Mon Nov 3 08:53:19 PST 2014
1 new commit in yt:
https://bitbucket.org/yt_analysis/yt/commits/1e0faa5cb0fc/
Changeset: 1e0faa5cb0fc
Branch: yt
User: MatthewTurk
Date: 2014-11-03 16:53:07+00:00
Summary: Merged in jzuhone/yt-3.x (pull request #1286)
Conversions for different electromagnetic unit systems, Unit equivalencies
Affected #: 16 files
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb doc/source/analyzing/units/2)_Fields_and_unit_conversion.ipynb
--- a/doc/source/analyzing/units/2)_Fields_and_unit_conversion.ipynb
+++ b/doc/source/analyzing/units/2)_Fields_and_unit_conversion.ipynb
@@ -1,7 +1,7 @@
{
"metadata": {
"name": "",
- "signature": "sha256:2faff88abc93fe2bc9d91467db786a8b69ec3ece6783a7055942ecc7c47a0817"
+ "signature": "sha256:c7cfb2db456d127bb633b7eee7ad6fe14290aa622ac62694c7840d80137afaba"
},
"nbformat": 3,
"nbformat_minor": 0,
@@ -36,7 +36,7 @@
"input": [
"import yt\n",
"ds = yt.load('IsolatedGalaxy/galaxy0030/galaxy0030')\n",
- " \n",
+ "\n",
"dd = ds.all_data()\n",
"maxval, maxloc = ds.find_max('density')\n",
"\n",
@@ -222,6 +222,69 @@
"level": 3,
"metadata": {},
"source": [
+ "Electrostatic/Electromagnetic Units"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Electromagnetic units can be a bit tricky, because the units for such quantities in different unit systems can have entirely different dimensions, even if they are meant to represent the same physical quantities. For example, in the SI system of units, current in Amperes is a fundamental unit of measure, so the unit of charge \"coulomb\" is equal to one ampere-second. On the other hand, in the Gaussian/CGS system, there is no equivalent base electromagnetic unit, and the electrostatic charge unit \"esu\" is equal to one $\\mathrm{cm^{3/2}g^{-1/2}s^{-1}}$ (which does not have any apparent physical significance). `yt` recognizes this difference:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "q1 = yt.YTArray(1.0,\"C\") # Coulombs\n",
+ "q2 = yt.YTArray(1.0,\"esu\") # electrostatic units / statcoulomb\n",
+ "\n",
+ "print \"units =\", q1.in_mks().units, \", dims =\", q1.units.dimensions\n",
+ "print \"units =\", q2.in_cgs().units, \", dims =\", q2.units.dimensions"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Under the hood, the `yt` units system has a translation layer that converts between these two systems, without any further effort required. For example:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "from yt.utilities.physical_constants import elementary_charge\n",
+ "\n",
+ "print elementary_charge\n",
+ "elementary_charge_C = elementary_charge.in_units(\"C\")\n",
+ "print elementary_charge_C"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "The electromagnetic unit translations `yt` understands are:\n",
+ "\n",
+ "* Charge: 1 coulomb (C) $\\leftrightarrow$ 0.1c electrostatic unit (esu, Fr)\n",
+ "* Current: 1 ampere (A, C/s) $\\leftrightarrow$ 0.1c statampere (statA, esu/s, Fr) \n",
+ "* Magnetic Field: 1 tesla (T) $\\leftrightarrow 10^4$ gauss (G)\n",
+ "\n",
+ "where \"Fr\" is the franklin, an alternative name for the electrostatic unit, and c is the speed of light. "
+ ]
+ },
+ {
+ "cell_type": "heading",
+ "level": 3,
+ "metadata": {},
+ "source": [
"Working with views and converting to ndarray"
]
},
@@ -324,10 +387,9 @@
"collapsed": false,
"input": [
"from astropy import units as u\n",
- "from yt import YTQuantity, YTArray\n",
"\n",
"x = 42.0 * u.meter\n",
- "y = YTQuantity.from_astropy(x) "
+ "y = yt.YTQuantity.from_astropy(x) "
],
"language": "python",
"metadata": {},
@@ -349,7 +411,7 @@
"collapsed": false,
"input": [
"a = np.random.random(size=10) * u.km/u.s\n",
- "b = YTArray.from_astropy(a)"
+ "b = yt.YTArray.from_astropy(a)"
],
"language": "python",
"metadata": {},
@@ -436,7 +498,7 @@
"collapsed": false,
"input": [
"k1 = kboltz.to_astropy()\n",
- "k2 = YTQuantity.from_astropy(kb)\n",
+ "k2 = yt.YTQuantity.from_astropy(kb)\n",
"print k1 == k2"
],
"language": "python",
@@ -447,7 +509,7 @@
"cell_type": "code",
"collapsed": false,
"input": [
- "c = YTArray.from_astropy(a)\n",
+ "c = yt.YTArray.from_astropy(a)\n",
"d = c.to_astropy()\n",
"print a == d"
],
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb doc/source/analyzing/units/6)_Unit_Equivalencies.ipynb
--- /dev/null
+++ b/doc/source/analyzing/units/6)_Unit_Equivalencies.ipynb
@@ -0,0 +1,179 @@
+{
+ "metadata": {
+ "name": "",
+ "signature": "sha256:cee652d703dd3369d81ebc670882d3734f73d0274aab98823a784d8039355480"
+ },
+ "nbformat": 3,
+ "nbformat_minor": 0,
+ "worksheets": [
+ {
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Some physical quantities are directly related to other unitful quantities by a constant, but otherwise do not have the same units. To facilitate conversions between these quantities, `yt` implements a system of unit equivalencies (inspired by the [AstroPy implementation](http://docs.astropy.org/en/latest/units/equivalencies.html)). The possible unit equivalencies are:\n",
+ "\n",
+ "* `\"thermal\"`: conversions between temperature and energy ($E = k_BT$)\n",
+ "* `\"spectral\"`: conversions between wavelength, frequency, and energy for photons ($E = h\\nu = hc/\\lambda, c = \\lambda\\nu$)\n",
+ "* `\"mass_energy\"`: conversions between mass and energy ($E = mc^2$)\n",
+ "* `\"lorentz\"`: conversions between velocity and Lorentz factor ($\\gamma = 1/\\sqrt{1-(v/c)^2}$)\n",
+ "* `\"schwarzschild\"`: conversions between mass and Schwarzschild radius ($R_S = 2GM/c^2$)\n",
+ "* `\"compton\"`: conversions between mass and Compton wavelength ($\\lambda = h/mc$)\n",
+ "\n",
+ "The following unit equivalencies only apply under conditions applicable for an ideal gas with a constant mean molecular weight $\\mu$ and ratio of specific heats $\\gamma$:\n",
+ "\n",
+ "* `\"number_density\"`: conversions between density and number density ($n = \\rho/\\mu{m_p}$)\n",
+ "* `\"sound_speed\"`: conversions between temperature and sound speed for an ideal gas ($c_s^2 = \\gamma{k_BT}/\\mu{m_p}$)\n",
+ "\n",
+ "A `YTArray` or `YTQuantity` can be converted to an equivalent using `to_equivalent`, where the unit and the equivalence name are provided as arguments:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "import yt\n",
+ "import numpy as np\n",
+ "\n",
+ "ds = yt.load('IsolatedGalaxy/galaxy0030/galaxy0030')\n",
+ "\n",
+ "dd = ds.all_data()\n",
+ "\n",
+ "print dd[\"temperature\"].to_equivalent(\"erg\", \"thermal\")\n",
+ "print dd[\"temperature\"].to_equivalent(\"eV\", \"thermal\")\n",
+ "\n",
+ "# Rest energy of the proton\n",
+ "from yt.utilities.physical_constants import mp\n",
+ "E_p = mp.to_equivalent(\"GeV\", \"mass_energy\")\n",
+ "print E_p"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Equivalencies can go in both directions, without any information required other than the unit you want to convert to:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "from yt.utilities.physical_constants import clight\n",
+ "v = 0.1*clight\n",
+ "g = v.to_equivalent(\"dimensionless\", \"lorentz\")\n",
+ "print g\n",
+ "print g.to_equivalent(\"c\", \"lorentz\")"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "heading",
+ "level": 3,
+ "metadata": {},
+ "source": [
+ "Special Equivalencies"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Some equivalencies can take supplemental information. The `\"number_density\"` equivalence can take a custom mean molecular weight (default is $\\mu = 0.6$):"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "print dd[\"density\"].max()\n",
+ "print dd[\"density\"].to_equivalent(\"cm**-3\", \"number_density\").max()\n",
+ "print dd[\"density\"].to_equivalent(\"cm**-3\", \"number_density\", mu=0.75).max()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "The `\"sound_speed\"` equivalence optionally takes the ratio of specific heats $\\gamma$ and the mean molecular weight $\\mu$ (defaults are $\\gamma$ = 5/3, $\\mu = 0.6$):"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "print dd[\"temperature\"].to_equivalent(\"km/s\", \"sound_speed\").mean()\n",
+ "print dd[\"temperature\"].to_equivalent(\"km/s\", \"sound_speed\", gamma=4./3., mu=0.5).mean()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "These options must be used with caution, and only if you know the underlying data adheres to these assumptions!"
+ ]
+ },
+ {
+ "cell_type": "heading",
+ "level": 3,
+ "metadata": {},
+ "source": [
+ "Determining Valid Equivalencies"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "If a certain equivalence does not exist for a particular unit, then an error will be thrown:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "from yt.utilities.exceptions import YTInvalidUnitEquivalence\n",
+ "\n",
+ "try:\n",
+ " x = v.to_equivalent(\"angstrom\", \"spectral\")\n",
+ "except YTInvalidUnitEquivalence as e:\n",
+ " print e"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "To list the equivalencies available for a given `YTArray` or `YTQuantity`, use the `list_equivalencies` method:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "E_p.list_equivalencies()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ }
+ ],
+ "metadata": {}
+ }
+ ]
+}
\ No newline at end of file
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb doc/source/analyzing/units/index.rst
--- a/doc/source/analyzing/units/index.rst
+++ b/doc/source/analyzing/units/index.rst
@@ -33,6 +33,7 @@
comoving_units_and_code_units
comparing_units_from_different_datasets
units_and_plotting
+ unit_equivalencies
.. note::
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb doc/source/analyzing/units/unit_equivalencies.rst
--- /dev/null
+++ b/doc/source/analyzing/units/unit_equivalencies.rst
@@ -0,0 +1,7 @@
+.. _symbolic_units:
+
+Symbolic units: :code:`yt.units`
+================================
+
+.. notebook:: 6)_Unit_Equivalencies.ipynb
+ :skip_exceptions:
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/analysis_modules/photon_simulator/spectral_models.py
--- a/yt/analysis_modules/photon_simulator/spectral_models.py
+++ b/yt/analysis_modules/photon_simulator/spectral_models.py
@@ -39,7 +39,7 @@
self.emax = emax*units.keV
self.nchan = nchan
self.ebins = np.linspace(emin, emax, nchan+1)*units.keV
- self.de = np.diff(self.ebins)*units.keV
+ self.de = np.diff(self.ebins)
self.emid = 0.5*(self.ebins[1:]+self.ebins[:-1])
def prepare(self):
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/fields/astro_fields.py
--- a/yt/fields/astro_fields.py
+++ b/yt/fields/astro_fields.py
@@ -102,6 +102,15 @@
function=_xray_emissivity,
units="") # add correct units here
+ def _mazzotta_weighting(field, data):
+ # Spectroscopic-like weighting field for galaxy clusters
+ # Only useful as a weight_field for temperature, metallicity, velocity
+ return data["density"]*data["density"]*data["kT"]**-0.25/mh/mh
+
+ registry.add_field((ftype,"mazzotta_weighting"),
+ function=_mazzotta_weighting,
+ units="keV**-0.25*cm**-6")
+
def _sz_kinetic(field, data):
scale = 0.88 * sigma_thompson / mh / clight
vel_axis = data.get_field_parameter("axis")
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/dimensions.py
--- a/yt/units/dimensions.py
+++ b/yt/units/dimensions.py
@@ -19,9 +19,11 @@
time = Symbol("(time)", positive=True)
temperature = Symbol("(temperature)", positive=True)
angle = Symbol("(angle)", positive=True)
+current_mks = Symbol("(current_mks)", positive=True)
dimensionless = sympify(1)
-base_dimensions = [mass, length, time, temperature, angle, dimensionless]
+base_dimensions = [mass, length, time, temperature, angle, current_mks,
+ dimensionless]
#
# Derived dimensions
@@ -44,16 +46,31 @@
power = energy / time
flux = power / area
specific_flux = flux / rate
-charge = (energy * length)**Rational(1, 2) # proper 1/2 power
+number_density = 1/(length*length*length)
+density = mass * number_density
-electric_field = charge / length**2
-magnetic_field = electric_field
+# Gaussian electromagnetic units
+charge_cgs = (energy * length)**Rational(1, 2) # proper 1/2 power
+current_cgs = charge_cgs / time
+electric_field_cgs = charge_cgs / length**2
+magnetic_field_cgs = electric_field_cgs
+
+# SI electromagnetic units
+charge_mks = current_mks * time
+electric_field_mks = force / charge_mks
+magnetic_field_mks = electric_field_mks / velocity
+
+# Since cgs is our default, I'm adding these aliases for backwards-compatibility
+charge = charge_cgs
+electric_field = electric_field_cgs
+magnetic_field = magnetic_field_cgs
solid_angle = angle * angle
derived_dimensions = [rate, velocity, acceleration, jerk, snap, crackle, pop,
- momentum, force, energy, power, charge, electric_field,
- magnetic_field, solid_angle, flux, specific_flux, volume,
- area]
+ momentum, force, energy, power, charge_cgs, electric_field_cgs,
+ magnetic_field_cgs, solid_angle, flux, specific_flux, volume,
+ area, current_cgs, charge_mks, electric_field_mks,
+ magnetic_field_mks]
dimensions = base_dimensions + derived_dimensions
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/equivalencies.py
--- /dev/null
+++ b/yt/units/equivalencies.py
@@ -0,0 +1,166 @@
+"""
+Equivalencies between different kinds of units
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import yt.utilities.physical_constants as pc
+from yt.units.dimensions import temperature, mass, energy, length, rate, \
+ velocity, dimensionless, density, number_density, flux
+from yt.extern.six import add_metaclass
+import numpy as np
+
+equivalence_registry = {}
+
+class RegisteredEquivalence(type):
+ def __init__(cls, name, b, d):
+ type.__init__(cls, name, b, d)
+ if hasattr(cls, "_type_name") and not cls._skip_add:
+ equivalence_registry[cls._type_name] = cls
+
+ at add_metaclass(RegisteredEquivalence)
+class Equivalence(object):
+ _skip_add = False
+
+class NumberDensityEquivalence(Equivalence):
+ _type_name = "number_density"
+ dims = (density,number_density,)
+
+ def convert(self, x, new_dims, mu=0.6):
+ if new_dims == number_density:
+ return x/(mu*pc.mh)
+ elif new_dims == density:
+ return x*mu*pc.mh
+
+ def __str__(self):
+ return "number density: density <-> number density"
+
+class ThermalEquivalence(Equivalence):
+ _type_name = "thermal"
+ dims = (temperature,energy,)
+
+ def convert(self, x, new_dims):
+ if new_dims == energy:
+ return pc.kboltz*x
+ elif new_dims == temperature:
+ return x/pc.kboltz
+
+ def __str__(self):
+ return "thermal: temperature <-> energy"
+
+class MassEnergyEquivalence(Equivalence):
+ _type_name = "mass_energy"
+ dims = (mass,energy,)
+
+ def convert(self, x, new_dims):
+ if new_dims == energy:
+ return x*pc.clight*pc.clight
+ elif new_dims == mass:
+ return x/(pc.clight*pc.clight)
+
+ def __str__(self):
+ return "mass_energy: mass <-> energy"
+
+class SpectralEquivalence(Equivalence):
+ _type_name = "spectral"
+ dims = (length,rate,energy,)
+
+ def convert(self, x, new_dims):
+ if new_dims == energy:
+ if x.units.dimensions == length:
+ nu = pc.clight/x
+ elif x.units.dimensions == rate:
+ nu = x
+ return pc.hcgs*nu
+ elif new_dims == length:
+ if x.units.dimensions == rate:
+ return pc.clight/x
+ elif x.units.dimensions == energy:
+ return pc.hcgs*pc.clight/x
+ elif new_dims == rate:
+ if x.units.dimensions == length:
+ return pc.clight/x
+ elif x.units.dimensions == energy:
+ return x/pc.hcgs
+
+ def __str__(self):
+ return "spectral: length <-> rate <-> energy"
+
+class SoundSpeedEquivalence(Equivalence):
+ _type_name = "sound_speed"
+ dims = (velocity,temperature,energy,)
+
+ def convert(self, x, new_dims, mu=0.6, gamma=5./3.):
+ if new_dims == velocity:
+ if x.units.dimensions == temperature:
+ kT = pc.kboltz*x
+ elif x.units.dimensions == energy:
+ kT = x
+ return np.sqrt(gamma*kT/(mu*pc.mh))
+ else:
+ kT = x*x*mu*pc.mh/gamma
+ if new_dims == temperature:
+ return kT/pc.kboltz
+ else:
+ return kT
+
+ def __str__(self):
+ return "sound_speed (ideal gas): velocity <-> temperature <-> energy"
+
+class LorentzEquivalence(Equivalence):
+ _type_name = "lorentz"
+ dims = (dimensionless,velocity,)
+
+ def convert(self, x, new_dims):
+ if new_dims == dimensionless:
+ beta = x.in_cgs()/pc.clight
+ return 1./np.sqrt(1.-beta**2)
+ elif new_dims == velocity:
+ return pc.clight*np.sqrt(1.-1./(x*x))
+
+ def __str__(self):
+ return "lorentz: velocity <-> dimensionless"
+
+class SchwarzschildEquivalence(Equivalence):
+ _type_name = "schwarzschild"
+ dims = (mass,length,)
+
+ def convert(self, x, new_dims):
+ if new_dims == length:
+ return 2.*pc.G*x/(pc.clight*pc.clight)
+ elif new_dims == mass:
+ return 0.5*x*pc.clight*pc.clight/pc.G
+
+ def __str__(self):
+ return "schwarzschild: mass <-> length"
+
+class ComptonEquivalence(Equivalence):
+ _type_name = "compton"
+ dims = (mass,length,)
+
+ def convert(self, x, new_dims):
+ return pc.hcgs/(x*pc.clight)
+
+ def __str__(self):
+ return "compton: mass <-> length"
+
+class EffectiveTemperature(Equivalence):
+ _type_name = "effective_temperature"
+ dims = (flux,temperature,)
+
+ def convert(self, x, new_dims):
+ if new_dims == flux:
+ return pc.stefan_boltzmann_constant_cgs*x**4
+ elif new_dims == temperature:
+ return (x/pc.stefan_boltzmann_constant_cgs)**0.25
+
+ def __str__(self):
+ return "effective_temperature: flux <-> temperature"
+
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/tests/test_ytarray.py
--- a/yt/units/tests/test_ytarray.py
+++ b/yt/units/tests/test_ytarray.py
@@ -28,7 +28,8 @@
assert_array_equal, \
assert_equal, assert_raises, \
assert_array_almost_equal_nulp, \
- assert_array_almost_equal
+ assert_array_almost_equal, \
+ assert_allclose
from numpy import array
from yt.units.yt_array import \
YTArray, YTQuantity, \
@@ -471,21 +472,21 @@
km = YTQuantity(1, 'km')
balmy = YTQuantity(300, 'K')
- balmy_F = YTQuantity(80.33, 'F')
- balmy_C = YTQuantity(26.85, 'C')
+ balmy_F = YTQuantity(80.33, 'degF')
+ balmy_C = YTQuantity(26.85, 'degC')
balmy_R = YTQuantity(540, 'R')
- assert_array_almost_equal(balmy.in_units('F'), balmy_F)
- assert_array_almost_equal(balmy.in_units('C'), balmy_C)
+ assert_array_almost_equal(balmy.in_units('degF'), balmy_F)
+ assert_array_almost_equal(balmy.in_units('degC'), balmy_C)
assert_array_almost_equal(balmy.in_units('R'), balmy_R)
balmy_view = balmy.ndarray_view()
- balmy.convert_to_units('F')
+ balmy.convert_to_units('degF')
yield assert_true, balmy_view.base is balmy.base
yield assert_array_almost_equal, np.array(balmy), np.array(balmy_F)
- balmy.convert_to_units('C')
+ balmy.convert_to_units('degC')
yield assert_true, balmy_view.base is balmy.base
yield assert_array_almost_equal, np.array(balmy), np.array(balmy_C)
@@ -493,13 +494,13 @@
yield assert_true, balmy_view.base is balmy.base
yield assert_array_almost_equal, np.array(balmy), np.array(balmy_R)
- balmy.convert_to_units('F')
+ balmy.convert_to_units('degF')
yield assert_true, balmy_view.base is balmy.base
yield assert_array_almost_equal, np.array(balmy), np.array(balmy_F)
yield assert_raises, InvalidUnitOperation, np.multiply, balmy, km
- # Does CGS convergion from F to K work?
+ # Does CGS conversion from F to K work?
yield assert_array_almost_equal, balmy.in_cgs(), YTQuantity(300, 'K')
@@ -859,6 +860,114 @@
os.chdir(curdir)
shutil.rmtree(tmpdir)
+def test_cgs_conversions():
+ from yt.utilities.physical_constants import qp, eps_0, clight
+ from yt.units.dimensions import current_mks, time
+
+ qp_mks = qp.in_units("C")
+ yield assert_equal, qp_mks.units.dimensions, current_mks*time
+ yield assert_array_almost_equal, qp_mks.in_cgs(), qp
+ qp_mks_k = qp_mks.in_units("kC")
+ yield assert_equal, qp_mks.units.dimensions, current_mks*time
+ yield assert_array_almost_equal, qp_mks_k.in_cgs(), qp
+ yield assert_array_almost_equal, qp_mks_k.in_units("kesu"), qp.in_units("kesu")
+
+ K = 1.0/(4*np.pi*eps_0)
+ yield assert_array_almost_equal, K.in_cgs(), 1.0
+
+ B = YTQuantity(1.0, "T")
+ yield assert_array_almost_equal, B.in_units("gauss"), YTQuantity(1.0e4, "gauss")
+
+ I = YTQuantity(1.0, "A")
+ yield assert_array_almost_equal, I.in_units("statA"), YTQuantity(0.1*clight, "statA")
+
+def test_equivalencies():
+ from yt.utilities.physical_constants import clight, mp, kboltz, hcgs, mh, me, \
+ mass_sun_cgs, G, stefan_boltzmann_constant_cgs
+ import yt.units as u
+
+ # Mass-energy
+
+ E = mp.to_equivalent("keV","mass_energy")
+ yield assert_equal, E, mp*clight*clight
+ yield assert_allclose, mp, E.to_equivalent("g", "mass_energy")
+
+ # Thermal
+
+ T = YTQuantity(1.0e8,"K")
+ E = T.to_equivalent("W*hr","thermal")
+ yield assert_equal, E, (kboltz*T).in_units("W*hr")
+ yield assert_allclose, T, E.to_equivalent("K", "thermal")
+
+ # Spectral
+
+ l = YTQuantity(4000.,"angstrom")
+ nu = l.to_equivalent("Hz","spectral")
+ yield assert_equal, nu, clight/l
+ E = hcgs*nu
+ l2 = E.to_equivalent("angstrom", "spectral")
+ yield assert_allclose, l, l2
+ nu2 = clight/l2.in_units("cm")
+ yield assert_allclose, nu, nu2
+ E2 = nu2.to_equivalent("keV", "spectral")
+ yield assert_allclose, E2, E.in_units("keV")
+
+ # Sound-speed
+
+ mu = 0.6
+ gg = 5./3.
+ c_s = T.to_equivalent("km/s","sound_speed")
+ yield assert_equal, c_s, np.sqrt(gg*kboltz*T/(mu*mh))
+ yield assert_allclose, T, c_s.to_equivalent("K","sound_speed")
+
+ mu = 0.5
+ gg = 4./3.
+ c_s = T.to_equivalent("km/s","sound_speed", mu=mu, gamma=gg)
+ yield assert_equal, c_s, np.sqrt(gg*kboltz*T/(mu*mh))
+ yield assert_allclose, T, c_s.to_equivalent("K","sound_speed",
+ mu=mu, gamma=gg)
+
+ # Lorentz
+
+ v = 0.8*clight
+ g = v.to_equivalent("dimensionless","lorentz")
+ g2 = YTQuantity(1./np.sqrt(1.-0.8*0.8), "dimensionless")
+ yield assert_allclose, g, g2
+ v2 = g2.to_equivalent("mile/hr", "lorentz")
+ yield assert_allclose, v2, v.in_units("mile/hr")
+
+ # Schwarzschild
+
+ R = mass_sun_cgs.to_equivalent("kpc","schwarzschild")
+ yield assert_equal, R.in_cgs(), 2*G*mass_sun_cgs/(clight*clight)
+ yield assert_allclose, mass_sun_cgs, R.to_equivalent("g", "schwarzschild")
+
+ # Compton
+
+ l = me.to_equivalent("angstrom","compton")
+ yield assert_equal, l, hcgs/(me*clight)
+ yield assert_allclose, me, l.to_equivalent("g", "compton")
+
+ # Number density
+
+ rho = mp/u.cm**3
+
+ n = rho.to_equivalent("cm**-3","number_density")
+ yield assert_equal, n, rho/(mh*0.6)
+ yield assert_allclose, rho, n.to_equivalent("g/cm**3","number_density")
+
+ n = rho.to_equivalent("cm**-3","number_density", mu=0.75)
+ yield assert_equal, n, rho/(mh*0.75)
+ yield assert_allclose, rho, n.to_equivalent("g/cm**3","number_density", mu=0.75)
+
+ # Effective temperature
+
+ T = YTQuantity(1.0e4, "K")
+ F = T.to_equivalent("erg/s/cm**2","effective_temperature")
+ yield assert_equal, F, stefan_boltzmann_constant_cgs*T**4
+ yield assert_allclose, T, F.to_equivalent("K", "effective_temperature")
+
+
def test_numpy_wrappers():
a1 = YTArray([1, 2, 3], 'cm')
a2 = YTArray([2, 3, 4, 5, 6], 'cm')
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/unit_lookup_table.py
--- a/yt/units/unit_lookup_table.py
+++ b/yt/units/unit_lookup_table.py
@@ -20,7 +20,8 @@
metallicity_sun, erg_per_eV, amu_grams, mass_electron_grams, \
cm_per_ang, jansky_cgs, mass_jupiter_grams, mass_earth_grams, \
boltzmann_constant_erg_per_K, kelvin_per_rankine, \
- speed_of_light_cm_per_s
+ speed_of_light_cm_per_s, planck_length_cm, planck_charge_esu, \
+ planck_energy_erg, planck_mass_grams, planck_temperature_K, planck_time_s
import numpy as np
# Lookup a unit symbol with the symbol string, and provide a tuple with the
@@ -37,20 +38,25 @@
# other cgs
"dyne": (1.0, dimensions.force),
"erg": (1.0, dimensions.energy),
- "esu": (1.0, dimensions.charge),
- "gauss": (1.0, dimensions.magnetic_field),
- "C" : (1.0, dimensions.temperature, -273.15),
+ "esu": (1.0, dimensions.charge_cgs),
+ "gauss": (1.0, dimensions.magnetic_field_cgs),
+ "degC": (1.0, dimensions.temperature, -273.15),
+ "statA": (1.0, dimensions.current_cgs),
# some SI
"m": (1.0e2, dimensions.length),
"J": (1.0e7, dimensions.energy),
"W": (1.0e7, dimensions.power),
"Hz": (1.0, dimensions.rate),
+ "N": (1.0e5, dimensions.force),
+ "C": (0.1*speed_of_light_cm_per_s, dimensions.charge_mks),
+ "A": (0.1*speed_of_light_cm_per_s, dimensions.current_mks),
+ "T": (1.0e4, dimensions.magnetic_field_mks),
# Imperial units
"ft": (30.48, dimensions.length),
"mile": (160934, dimensions.length),
- "F": (kelvin_per_rankine, dimensions.temperature, -459.67),
+ "degF": (kelvin_per_rankine, dimensions.temperature, -459.67),
"R": (kelvin_per_rankine, dimensions.temperature),
# dimensionless stuff
@@ -93,13 +99,11 @@
# misc
"eV": (erg_per_eV, dimensions.energy),
"amu": (amu_grams, dimensions.mass),
- "me": (mass_electron_grams, dimensions.mass),
"angstrom": (cm_per_ang, dimensions.length),
"Jy": (jansky_cgs, dimensions.specific_flux),
"counts": (1.0, dimensions.dimensionless),
- "kB": (boltzmann_constant_erg_per_K,
- dimensions.energy/dimensions.temperature),
"photons": (1.0, dimensions.dimensionless),
+ "me": (mass_electron_grams, dimensions.mass),
# for AstroPy compatibility
"solMass": (mass_sun_grams, dimensions.mass),
@@ -109,11 +113,19 @@
"sr": (1.0, dimensions.solid_angle),
"rad": (1.0, dimensions.solid_angle),
"deg": (np.pi/180., dimensions.angle),
- "Fr": (1.0, dimensions.charge),
- "G": (1.0, dimensions.magnetic_field),
+ "Fr": (1.0, dimensions.charge_cgs),
+ "G": (1.0, dimensions.magnetic_field_cgs),
"d": (1.0, dimensions.time),
"Angstrom": (cm_per_ang, dimensions.length),
+ # Planck units
+ "m_pl": (planck_mass_grams, dimensions.mass),
+ "l_pl": (planck_length_cm, dimensions.length),
+ "t_pl": (planck_time_s, dimensions.time),
+ "T_pl": (planck_temperature_K, dimensions.temperature),
+ "q_pl": (planck_charge_esu, dimensions.charge_cgs),
+ "E_pl": (planck_energy_erg, dimensions.energy),
+
}
# Add LaTeX representations for units with trivial representations.
@@ -177,7 +189,12 @@
"Hz",
"W",
"gauss",
+ "G",
"Jy",
+ "N",
+ "T",
+ "A",
+ "C",
)
cgs_base_units = {
@@ -194,4 +211,16 @@
dimensions.time:'s',
dimensions.temperature:'K',
dimensions.angle:'radian',
+ dimensions.current_mks:'A',
}
+
+cgs_conversions = {
+ "C":"esu",
+ "T":"gauss",
+ "A":"statA",
+}
+
+for conv in cgs_conversions.keys():
+ if conv in prefixable_units:
+ for p in unit_prefixes:
+ cgs_conversions[p+conv] = p+cgs_conversions[conv]
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/unit_object.py
--- a/yt/units/unit_object.py
+++ b/yt/units/unit_object.py
@@ -17,7 +17,7 @@
Pow, Symbol, Integer, \
Float, Basic, Rational, sqrt
from sympy.core.numbers import One
-from sympy import sympify, latex
+from sympy import sympify, latex, symbols
from sympy.parsing.sympy_parser import \
parse_expr, auto_number, rationalize
from keyword import iskeyword
@@ -26,7 +26,7 @@
from yt.units.unit_lookup_table import \
latex_symbol_lut, unit_prefixes, \
prefixable_units, cgs_base_units, \
- mks_base_units
+ mks_base_units, cgs_conversions
from yt.units.unit_registry import UnitRegistry
import copy
@@ -116,7 +116,7 @@
# Extra attributes
__slots__ = ["expr", "is_atomic", "cgs_value", "cgs_offset", "dimensions",
- "registry"]
+ "registry", "cgs_conversion", "is_mks"]
def __new__(cls, unit_expr=sympy_one, cgs_value=None, cgs_offset=0.0,
dimensions=None, registry=None, **assumptions):
@@ -191,12 +191,11 @@
validate_dimensions(dimensions)
else:
# lookup the unit symbols
- try:
- cgs_value, dimensions = \
- _get_unit_data_from_expr(unit_expr, registry.lut)
- except ValueError:
- cgs_value, dimensions, cgs_offset = \
- _get_unit_data_from_expr(unit_expr, registry.lut)
+ unit_data = _get_unit_data_from_expr(unit_expr, registry.lut)
+ cgs_value = unit_data[0]
+ dimensions = unit_data[1]
+ if len(unit_data) == 3:
+ cgs_offset = unit_data[2]
# Create obj with superclass construct.
obj = Expr.__new__(cls, **assumptions)
@@ -209,6 +208,22 @@
obj.dimensions = dimensions
obj.registry = registry
+ check_atoms = [atom for atom in unit_expr.free_symbols
+ if str(atom) in cgs_conversions]
+ if len(check_atoms) > 0:
+ conversions = []
+ for atom in check_atoms:
+ conversions.append((atom,symbols(cgs_conversions[str(atom)])))
+ conversion = unit_expr.subs(conversions)
+ conversion = Unit(unit_expr=conversion, cgs_value=1.0,
+ dimensions=None, registry=registry)
+ is_mks = True
+ else:
+ conversion = None
+ is_mks = False
+ obj.cgs_conversion = conversion
+ obj.is_mks = is_mks
+
if unit_key:
registry.unit_objs[unit_key] = obj
@@ -264,7 +279,7 @@
cgs_offset = self.cgs_offset
else:
raise InvalidUnitOperation("Quantities with units of Fahrenheit "
- "and Celcius cannot be multiplied.")
+ "and Celsius cannot be multiplied.")
return Unit(self.expr * u.expr,
cgs_value=(self.cgs_value * u.cgs_value),
@@ -287,7 +302,7 @@
cgs_offset = self.cgs_offset
else:
raise InvalidUnitOperation("Quantities with units of Farhenheit "
- "and Celcius cannot be multiplied.")
+ "and Celsius cannot be multiplied.")
return Unit(self.expr / u.expr,
cgs_value=(self.cgs_value / u.cgs_value),
@@ -307,7 +322,8 @@
"it to a float." % (p, type(p)) )
return Unit(self.expr**p, cgs_value=(self.cgs_value**p),
- dimensions=(self.dimensions**p), registry=self.registry)
+ dimensions=(self.dimensions**p),
+ registry=self.registry)
def __eq__(self, u):
""" Test unit equality. """
@@ -340,6 +356,14 @@
def same_dimensions_as(self, other_unit):
""" Test if dimensions are the same. """
+ first_check = False
+ second_check = False
+ if self.cgs_conversion:
+ first_check = self.cgs_conversion.dimensions / other_unit.dimensions == sympy_one
+ if other_unit.cgs_conversion:
+ second_check = other_unit.cgs_conversion.dimensions / self.dimensions == sympy_one
+ if first_check or second_check:
+ return True
return (self.dimensions / other_unit.dimensions) == sympy_one
@property
@@ -372,20 +396,28 @@
Create and return dimensionally-equivalent cgs units.
"""
- units_string = self._get_system_unit_string(cgs_base_units)
+ if self.cgs_conversion:
+ units = self.cgs_conversion
+ else:
+ units = self
+ units_string = units._get_system_unit_string(cgs_base_units)
return Unit(units_string, cgs_value=1.0,
- dimensions=self.dimensions, registry=self.registry)
+ dimensions=units.dimensions, registry=self.registry)
def get_mks_equivalent(self):
"""
Create and return dimensionally-equivalent mks units.
"""
- units_string = self._get_system_unit_string(mks_base_units)
- cgs_value = (get_conversion_factor(self, self.get_cgs_equivalent())[0] /
- get_conversion_factor(self, Unit(units_string))[0])
+ if self.cgs_conversion and not self.is_mks:
+ units = self.cgs_conversion
+ else:
+ units = self
+ units_string = units._get_system_unit_string(mks_base_units)
+ cgs_value = (get_conversion_factor(units, units.get_cgs_equivalent())[0] /
+ get_conversion_factor(units, Unit(units_string))[0])
return Unit(units_string, cgs_value=cgs_value,
- dimensions=self.dimensions, registry=self.registry)
+ dimensions=units.dimensions, registry=self.registry)
def get_conversion_factor(self, other_units):
return get_conversion_factor(self, other_units)
@@ -430,7 +462,7 @@
return ratio, ratio*old_units.cgs_offset - new_units.cgs_offset
else:
raise InvalidUnitOperation(
- "Fahrenheit and Celsius are not absoulte temperature scales "
+ "Fahrenheit and Celsius are not absolute temperature scales "
"and cannot be used in compound unit symbols.")
#
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/unit_symbols.py
--- a/yt/units/unit_symbols.py
+++ b/yt/units/unit_symbols.py
@@ -140,8 +140,8 @@
MeV = mega_electron_volt = quan(1.0, "MeV")
GeV = giga_electron_volt = quan(1.0, "GeV")
amu = atomic_mass_unit = quan(1.0, "amu")
+angstrom = quan(1.0, "angstrom")
me = electron_mass = quan(1.0, "me")
-angstrom = quan(1.0, "angstrom")
#
# Angle units
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/units/yt_array.py
--- a/yt/units/yt_array.py
+++ b/yt/units/yt_array.py
@@ -33,10 +33,12 @@
from yt.units.dimensions import dimensionless
from yt.utilities.exceptions import \
YTUnitOperationError, YTUnitConversionError, \
- YTUfuncUnitError, YTIterableUnitCoercionError
+ YTUfuncUnitError, YTIterableUnitCoercionError, \
+ YTInvalidUnitEquivalence
from numbers import Number as numeric_type
from yt.utilities.on_demand_imports import _astropy
from sympy import Rational
+from yt.units.unit_lookup_table import unit_prefixes, prefixable_units
NULL_UNIT = Unit()
@@ -458,6 +460,44 @@
"""
return self.in_units(self.units.get_mks_equivalent())
+ def to_equivalent(self, unit, equiv, **kwargs):
+ """
+ Convert a YTArray or YTQuantity to an equivalent, e.g., something that is
+ related by only a constant factor but not in the same units.
+
+ Parameters
+ ----------
+ unit : string
+ 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 :method:`list_equivalencies`
+ method.
+
+ Examples
+ --------
+ >>> a = yt.YTArray(1.0e7,"K")
+ >>> a.to_equivalent("keV", "thermal")
+ """
+ from equivalencies import equivalence_registry
+ this_equiv = equivalence_registry[equiv]()
+ old_dims = self.units.dimensions
+ new_dims = YTQuantity(1.0, unit, registry=self.units.registry).units.dimensions
+ if old_dims in this_equiv.dims and new_dims in this_equiv.dims:
+ return this_equiv.convert(self, new_dims, **kwargs).in_units(unit)
+ else:
+ raise YTInvalidUnitEquivalence(equiv, self.units, unit)
+
+ def list_equivalencies(self):
+ """
+ Lists the possible equivalencies associated with this YTArray or
+ YTQuantity.
+ """
+ from equivalencies import equivalence_registry
+ for k,v in equivalence_registry.items():
+ if self.units.dimensions in v.dims:
+ print v()
+
def ndarray_view(self):
"""
Returns a view into the array, but as an ndarray rather than ytarray.
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/utilities/exceptions.py
--- a/yt/utilities/exceptions.py
+++ b/yt/utilities/exceptions.py
@@ -426,7 +426,6 @@
def __str__(self):
return "A file already exists at %s and clobber=False." % self.filename
-
class YTGDFUnknownGeometry(Exception):
def __init__(self, geometry):
self.geometry = geometry
@@ -435,4 +434,14 @@
return '''Unknown geometry %i. Please refer to GDF standard
for more information''' % self.geometry
+class YTInvalidUnitEquivalence(Exception):
+ def __init__(self, equiv, unit1, unit2):
+ self.equiv = equiv
+ self.unit1 = unit1
+ self.unit2 = unit2
+ def __str__(self):
+ return "The unit equivalence '%s' does not exist for the units '%s' and '%s.'" % (self.equiv,
+ self.unit1,
+ self.unit2)
+
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/utilities/physical_constants.py
--- a/yt/utilities/physical_constants.py
+++ b/yt/utilities/physical_constants.py
@@ -1,5 +1,6 @@
from yt.utilities.physical_ratios import *
from yt.units.yt_array import YTQuantity
+from math import pi
mass_electron_cgs = YTQuantity(mass_electron_grams, 'g')
amu_cgs = YTQuantity(amu_grams, 'g')
@@ -14,6 +15,7 @@
# Charge
charge_proton_cgs = YTQuantity(4.8032056e-10, 'esu')
+elementary_charge = charge_proton_cgs
# Physical Constants
boltzmann_constant_cgs = YTQuantity(boltzmann_constant_erg_per_K, 'erg/K')
@@ -35,7 +37,7 @@
mass_mars_cgs = YTQuantity(mass_mars_grams, 'g')
mass_saturn_cgs = YTQuantity(mass_saturn_grams, 'g')
mass_uranus_cgs = YTQuantity(mass_uranus_grams, 'g')
-mass_neptun_cgs = YTQuantity(mass_neptun_grams, 'g')
+mass_neptune_cgs = YTQuantity(mass_neptune_grams, 'g')
#Short cuts
G = gravitational_constant_cgs
@@ -48,5 +50,17 @@
kboltz = boltzmann_constant_cgs
kb = kboltz
hcgs = planck_constant_cgs
+hbar = 0.5*hcgs/pi
sigma_thompson = cross_section_thompson_cgs
Na = 1 / amu_cgs
+
+#Planck units
+m_pl = planck_mass = YTQuantity(planck_mass_grams, "g")
+l_pl = planck_length = YTQuantity(planck_length_cm, "cm")
+t_pl = planck_time = YTQuantity(planck_time_s, "s")
+E_pl = planck_energy = YTQuantity(planck_energy_erg, "erg")
+q_pl = planck_charge = YTQuantity(planck_charge_esu, "esu")
+T_pl = planck_temperature = YTQuantity(planck_temperature_K, "K")
+
+mu_0 = YTQuantity(4.0e-7*pi, "N/A**2")
+eps_0 = (1.0/(clight.in_mks()**2*mu_0)).in_units("C**2/N/m**2")
diff -r c91388a959906a04a5db1e773537e94bcb636a62 -r 1e0faa5cb0fc8a94a358a7d90fd672176db42eeb yt/utilities/physical_ratios.py
--- a/yt/utilities/physical_ratios.py
+++ b/yt/utilities/physical_ratios.py
@@ -92,7 +92,7 @@
mass_mars_grams = mass_sun_grams / 3098708.0
mass_saturn_grams = mass_sun_grams / 3497.898
mass_uranus_grams = mass_sun_grams / 22902.98
-mass_neptun_grams = mass_sun_grams / 19412.24
+mass_neptune_grams = mass_sun_grams / 19412.24
# flux
jansky_cgs = 1.0e-23
@@ -109,3 +109,11 @@
# Miscellaneous
HUGE = 1.0e90
TINY = 1.0e-40
+
+# Planck units
+planck_mass_grams = 2.17650925245e-05
+planck_length_cm = 1.6161992557e-33
+planck_time_s = planck_length_cm / speed_of_light_cm_per_s
+planck_energy_erg = planck_mass_grams * speed_of_light_cm_per_s * speed_of_light_cm_per_s
+planck_temperature_K = planck_energy_erg / boltzmann_constant_erg_per_K
+planck_charge_esu = 5.62274532302e-09
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