{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Working with SI-units\n", "\n", "## Goal\n", "\n", "> Learn how to work with `SINumber` and `SIArray` objects which represent physical quantities, i.e. one or more floating point numbers with an associated unit.\n", "\n", "## The `feos.si` module\n", "\n", "Most interfaces in `FeOs` use dimensioned quantities as input. For example, to define a thermodynamic state at given temperature, pressure and amount of substance, all of these properties have to be multiplied by an apropriate unit before we can call the function that creates the state.\n", "\n", "`FeOs` uses the [quantity](https://itt-ustutt.github.io/quantity/) crate which generates the `si_units` Python module. You might find this module useful on it's own in which case you can also install it as a separate python package.\n", "For convenient use with `FeOs` however, it is re-exported as `feos.si`." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from feos.si import *\n", "import numpy as np" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `si` module contains units according to [the Standard Interational System of Units](https://en.wikipedia.org/wiki/International_System_of_Units) (SI), constants and prefixes.\n", "A **scalar** floating point number multiplied or divided by a unit has the **`SINumber`** data type." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "temperature : 300.15 K\n", "data type : \n" ] } ], "source": [ "temperature = 300.15 * KELVIN\n", "print(f\"temperature : {temperature}\")\n", "print(f\"data type : {type(temperature)}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The representation of a `SINumber` (e.g. using `print`) uses a prefix (e.g. `k` for `kilo`) so that the numerical value is convenient. For example 1000 g will be represented as 1 kg. Internally, all values are stored with respect to the basic SI unit, i.e. `METER`, `KILOGRAM`, `SECOND`, `AMPERE`, `MOL`, `KELVIN`, and `CANDELA`." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$1\\,\\mathrm{kg}$" ], "text/plain": [ "1 kg" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "1000 * GRAM" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "v = 0.07486757864516083 m³\n" ] } ], "source": [ "# volume of an ideal gas\n", "t = 300.15 * KELVIN\n", "p = 0.5 * BAR\n", "n = 1.5 * MOL\n", "v = n * RGAS * t / p\n", "print(f\"v = {v}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can use division to transform a `SINumber` to a `float` for the desired unit:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "t = 27.00 °C\n", "data type = \n" ] } ], "source": [ "# temperature in Celsius\n", "t_c = t / CELSIUS\n", "print(f\"t = {t_c:3.2f} °C\")\n", "print(f\"data type = {type(t_c)}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Mathematical operations and interaction with `numpy`\n", "\n", "The `SINumber` object supports some mathematical operations, i.e. we can\n", "\n", "- add or subtract two objects, **if they have the same unit**,\n", "- multiply or divide two objects, which returns a floating point number if they have the same unit,\n", "- raise an object to a power, and\n", "- take the square and cubic root, which only works if the units have the correct exponents to allow for the operation." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "pressure : 265 kPa\n", "temperature : 40 K\n", "velocity : 100 m/s\n", "velocity (cm/min) : 600000.0\n", "acceleration : 9.81 m s^-2\n" ] } ], "source": [ "# addition\n", "pressure = 2.5 * BAR + 15_000 * PASCAL\n", "print('pressure :', pressure)\n", "\n", "# subtraction\n", "temperature = 543.15 * KELVIN - 230.0 * CELSIUS\n", "print('temperature :', temperature)\n", "\n", "# division\n", "velocity = 360_000 * METER / HOUR\n", "print('velocity :', velocity)\n", "\n", "# division as transformation to target unit\n", "v_cm_minute = velocity / (CENTI * METER / MINUTE) # this is a floating point number\n", "print(f'velocity (cm/min) : {v_cm_minute:8.1f}')\n", "\n", "# power\n", "acceleration = 9.81 * METER / SECOND**2\n", "print('acceleration :', acceleration)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For the square and cubic root we can use the `numpy` functions, `numpy.sqrt` and `numpy.cbrt`, respectively, or use the methods provided by `si_units`." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "speed : 1 m/s\n", "length : 3 nm\n" ] } ], "source": [ "# square root\n", "speed = np.sqrt(METER**2 / SECOND**2)\n", "print('speed :', speed)\n", "\n", "# cubic root\n", "box_length = (27_000 * ANGSTROM**3).cbrt()\n", "print('length :', box_length)\n", "\n", "# both numpy and methods of SINumbers work\n", "assert(np.sqrt(METER**2 / SECOND**2) == (METER**2 / SECOND**2).sqrt())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The mathematical functions from the `math` module do not work because it is written in C and `SINumber`s are not compatible with the C function arguments.\n", "If you are having trouble with errors, consider transforming your properties to floats (by division with the proper unit) and then multiplying the correct unit to the result." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "ERROR: must be real number, not si_units.SINumber\n" ] } ], "source": [ "from math import sqrt\n", "\n", "try:\n", " sqrt(METER**2 / SECOND**2)\n", "except Exception as e:\n", " print(\"ERROR:\", e)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If you try to perform arithmetic operations with `SINumber`s that are not allowed (e.g. adding two properties with different units), an error is raised." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "ename": "PanicException", "evalue": "Inconsistent units Pa + K", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mPanicException\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mPASCAL\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mKELVIN\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mPanicException\u001b[0m: Inconsistent units Pa + K" ] } ], "source": [ "PASCAL + KELVIN" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can enforce correct units in our interfaces using the `has_unit` method:" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "p1 = 1.4550309581768168 MPa\n", "ERROR: Please provide the molar density, e.g. in units of mol/m³.\n" ] } ], "source": [ "MOL_M3 = MOL / METER**3 # we can create own units for future use\n", "\n", "def ideal_gas_pressure(density, temperature):\n", " if not density.has_unit(MOL_M3):\n", " raise ValueError(\"Please provide the molar density, e.g. in units of mol/m³.\")\n", " else:\n", " return density * temperature * RGAS\n", " \n", "try:\n", " p1 = ideal_gas_pressure(0.5 * KILO * MOL / METER**3, 350 * KELVIN)\n", " print('p1 = ', p1)\n", " p2 = ideal_gas_pressure(0.5 * KILOGRAM / METER**3, 350 * KELVIN)\n", " print('p2 = ', p2)\n", "except Exception as e:\n", " print(\"ERROR:\", e) " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Arrays of quantities\n", "\n", "A **numpy array** of floating point numbers multiplied or divided by a unit has the `SIArrayX` data type, where the `X` stands for the dimension of the numpy array." ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "pressures = [100000, 200000] Pa\n", "data type = \n" ] } ], "source": [ "ps = np.array([1.0, 2.0]) * BAR\n", "print(f\"pressures = {ps}\")\n", "print(f\"data type = {type(ps)}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Currently, indexing into an array only works for arrays of dimension one." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "200 kPa\n" ] } ], "source": [ "print(ps[1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are several useful numpy methods to create arrays. Most of them can be simply multiplied by units to convert them to `SIArray` objects.\n", "Some of these functions, such as `linspace` and `logspace` can directly be constructed using `SIArray1.linspace` and `SIArray1.logspace`, respectively." ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "pressures (numpy) = [100000, 125000, 150000, 175000, 200000] Pa\n", "pressures (si) = [100000, 125000, 150000, 175000, 200000] Pa\n" ] } ], "source": [ "ps_np = np.linspace(1, 2, 5) * BAR\n", "ps_si = SIArray1.linspace(1 * BAR, 2 * BAR, 5)\n", "print(f\"pressures (numpy) = {ps_np}\")\n", "print(f\"pressures (si) = {ps_si}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Just like `SINumber`, division by a unit that matches the property stored in an array yields a numpy ndarray." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "pressures = [0.1 0.125 0.15 0.175 0.2 ] MPa\n", "data type = \n" ] } ], "source": [ "ps_mpa = ps_np / (MEGA*PASCAL)\n", "print(f\"pressures = {ps_mpa} MPa\")\n", "print(f\"data type = {type(ps_mpa)}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Sometimes, we write functions that use arrays as input where it is possible that these arrays only contain one element.\n", "A numpy array with a single element qualifies for multiplication rules of a *scalar* and a `SINumber` which can cause problems.\n", "\n", "Consider the following example:" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Array \"[1, 2] mol\" contains 2 elements.\n", "\n", "ERROR: object of type 'si_units.SINumber' has no len()\n" ] } ], "source": [ "def print_length(array):\n", " print(f'Array \"{array}\" contains {len(array)} elements.\\n')\n", "\n", "n2 = np.array([1.0, 2.0]) * MOL # data type: SIArray1\n", "n1 = np.array([1.0]) * MOL # data type: SINumber\n", "\n", "try:\n", " print_length(n2)\n", " print_length(n1) # <- this fails, since n1 has data type SINumber, not SIArray1\n", "except Exception as e:\n", " print('ERROR:', e)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This example does not work, because `np.array([1.0]) * MOL` has data type `SINumber` instead of `SIArray1`. We can construct a `SIArray1` using the base constructor:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Array \"[1] mol\" contains 1 elements.\n", "\n" ] } ], "source": [ "# now it works\n", "n1 = SIArray1(1.0*MOL)\n", "print_length(n1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Derived units, constants and prefixes\n", "\n", "In conjunction to the base units, the `si` module also exports derived units (e.g. `HOUR = 3600 * SECOND`), constants and prefixes (such as `KILO` or `FEMTO`).\n", "Of course, we could multiply by a floating point number instead of using prefixes (e.g. `cm = 1e-2 * METER` vs. `cm = CENTI * METER`) but we think using prefixes make the code more readable.\n", "\n", "For a complete overview of exported units constants and prefixes, take a look at the [documentation of the `si_units` Python package](https://itt-ustutt.github.io/quantity/api.html#si-base-units-and-associated-constants)." ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Hyperfine transition frequency of Cs: 9.19263177 GHz\n", "Speed of light : 2.99792458e8 m/s\n", "Planck constant : 6.62607015e-34 Js\n", "Elementary charge : 1.602176634e-19 C\n", "Boltzmann constant : 1.380649e-23 J/K\n", "Avogradro constant : 6.02214076e23 mol^-1\n", "Luminous efficacy : 683 lm/W\n" ] } ], "source": [ "# The seven constants that inform the base units\n", "print('Hyperfine transition frequency of Cs:', DVCS)\n", "print('Speed of light :', CLIGHT)\n", "print('Planck constant :', PLANCK)\n", "print('Elementary charge :', QE)\n", "print('Boltzmann constant :', KB)\n", "print('Avogradro constant :', NAV)\n", "print('Luminous efficacy :', KCD)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Gravitational constant: 6.6743e-11 m³/kg/s²\n", "Ideal gas constant : 8.31446261815324 J/mol/K\n" ] } ], "source": [ "# Derived constants\n", "print('Gravitational constant:', G)\n", "print('Ideal gas constant :', RGAS)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Summary\n", "\n", "In `feos`, most interfaces use dimensioned units in form of `SINumber` and `SIArrayX` (where `X` stands for the dimnesion).\n", "This enables us to write functions that can check for proper units and thus circumvent unit errors that could occur, e.g. providing mass densities instead of molar densities.\n", "\n", "In this example, we learned how to create and use these objects in conjunction with methods provided by `numpy`.\n", "\n", "## Concluding remkars\n", "\n", "Hopefully you found this example helpful. If you have comments, critique or feedback, please let us know and consider [opening an issue on github](https://github.com/feos-org/feos/issues)." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.8.6" } }, "nbformat": 4, "nbformat_minor": 4 }