Skip to main content
Donate to the Python Software Foundation or Purchase a PyCharm License to Benefit the PSF! Donate Now

ChemPy is a Python package useful for solving problems in chemistry.

Project description

ChemPy
======

.. image:: http://hera.physchem.kth.se:9090/api/badges/bjodah/chempy/status.svg
:target: http://hera.physchem.kth.se:9090/bjodah/chempy
:alt: Build status
.. image:: https://img.shields.io/pypi/v/chempy.svg
:target: https://pypi.python.org/pypi/chempy
:alt: PyPI version
.. image:: https://img.shields.io/badge/python-2.7,3.5,3.6-blue.svg
:target: https://www.python.org/
:alt: Python version
.. image:: https://zenodo.org/badge/8840/bjodah/chempy.svg
:target: https://zenodo.org/badge/latestdoi/8840/bjodah/chempy
.. image:: https://img.shields.io/pypi/l/chempy.svg
:target: https://github.com/bjodah/chempy/blob/master/LICENSE
:alt: License
.. image:: http://img.shields.io/badge/benchmarked%20by-asv-green.svg?style=flat
:target: http://hera.physchem.kth.se/~chempy/benchmarks
:alt: airspeedvelocity
.. image:: http://hera.physchem.kth.se/~chempy/branches/master/htmlcov/coverage.svg
:target: http://hera.physchem.kth.se/~chempy/branches/master/htmlcov
:alt: coverage


.. contents::


About ChemPy
------------
`ChemPy <https://github.com/bjodah/chempy>`_ is a `Python <https://www.python.org>`_ package useful for
chemistry (mainly physical/inorganic/analytical chemistry). Currently it includes:

- Numerical integration routines for chemical kinetics (ODE solver front-end)
- Integrated rate expressions (and convenience fitting routines)
- Solver for equilibria (including multiphase systems)
- Relations in physical chemistry:

- Debye-Hückel expressions
- Arrhenius & Eyring equation
- Einstein-Smoluchowski equation

- Properties (pure python implementations from the litterature)

- water density as function of temperature
- water permittivity as function of temperature and pressure
- water diffusivity as function of temperature
- water viscosity as function of temperature
- sulfuric acid density as function of temperature & weight fraction H₂SO₄
- More to come... (and contributions are most welcome!)


Documentation
-------------
The easiest way to get started is to have a look at the examples in this README,
and also the jupyter notebooks_. In addition there is auto-generated API documentation
for the latest stable release `here <https://bjodah.github.io/chempy/latest>`_
(and the API docs for the development version is found
`here <http://hera.physchem.kth.se/~chempy/branches/master/html>`_).

.. _notebooks: http://hera.physchem.kth.se/~chempy/branches/master/examples

Installation
------------
Simplest way to install ChemPy and its (optional) dependencies is to use the
`conda package manager <https://conda.pydata.org/docs/>`_::

$ conda install -c bjodah chempy pytest
$ pytest -rs --pyargs chempy

Optional dependencies
~~~~~~~~~~~~~~~~~~~~~
If you used ``conda`` to install ChemPy you can skip this section.
But if you use ``pip`` the default installation is achieved by writing::

$ python -m pip install --user --upgrade chempy pytest
$ python -m pytest -rs --pyargs chempy

you can skip the ``--user`` flag if you have got root permissions.
You may be interested in using additional backends (in addition to those provided by SciPy)
for solving ODE-systems and non-linear optimization problems::

$ pip install chempy[all]

(see `setup.py <setup.py>`_ for optional requirements.). Note that this option requires you
to have the following libraries installed (including their development headers):

- `pygslodeiv2 <https://github.com/bjodah/pygslodeiv2>`_ (requires GSL_ >=1.16)
- `pyodeint <https://github.com/bjodah/pyodeint>`_ (requires boost_ >=1.65.0)
- `pycvodes <https://github.com/bjodah/pycvodes>`_ (requires SUNDIALS_ >=2.7.0)

.. _GSL: https://www.gnu.org/software/gsl/
.. _boost: http://www.boost.org/
.. _SUNDIALS: https://computation.llnl.gov/projects/sundials

if you want to see what packages need to be installed on a Debian based system you may look at this
`Dockerfile <scripts/environment/Dockerfile>`_.

Examples
--------
See demonstration scripts in `examples/ <https://github.com/bjodah/chempy/tree/master/examples>`_,
and some rendered jupyter notebooks_.
You may also browse the documentation for more examples. Below you will find a few code snippets:

Parsing formulae
~~~~~~~~~~~~~~~~
.. code:: python

>>> from chempy import Substance
>>> ferricyanide = Substance.from_formula('Fe(CN)6-3')
>>> ferricyanide.composition == {0: -3, 26: 1, 6: 6, 7: 6} # 0 for charge
True
>>> print(ferricyanide.unicode_name)
Fe(CN)₆³⁻
>>> print(ferricyanide.latex_name + ", " + ferricyanide.html_name)
Fe(CN)_{6}^{3-}, Fe(CN)<sub>6</sub><sup>3-</sup>
>>> print('%.3f' % ferricyanide.mass)
211.955


as you see, in composition, the atomic numbers (and 0 for charge) is used as
keys and the count of each kind became respective value.

Balancing stoichiometry of a chemical reaction
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
.. code:: python

>>> from chempy import balance_stoichiometry # Main reaction in NASA's booster rockets:
>>> reac, prod = balance_stoichiometry({'NH4ClO4', 'Al'}, {'Al2O3', 'HCl', 'H2O', 'N2'})
>>> from pprint import pprint
>>> pprint(dict(reac))
{'Al': 10, 'NH4ClO4': 6}
>>> pprint(dict(prod))
{'Al2O3': 5, 'H2O': 9, 'HCl': 6, 'N2': 3}
>>> from chempy import mass_fractions
>>> for fractions in map(mass_fractions, [reac, prod]):
... pprint({k: '{0:.3g} wt%'.format(v*100) for k, v in fractions.items()})
...
{'Al': '27.7 wt%', 'NH4ClO4': '72.3 wt%'}
{'Al2O3': '52.3 wt%', 'H2O': '16.6 wt%', 'HCl': '22.4 wt%', 'N2': '8.62 wt%'}

ChemPy can even handle reactions with linear dependencies (underdetermined systems), e.g.:

.. code:: python

>>> pprint([dict(_) for _ in balance_stoichiometry({'C', 'O2'}, {'CO2', 'CO'})]) # doctest: +SKIP
[{'C': x1 + 2, 'O2': x1 + 1}, {'CO': 2, 'CO2': x1}]

that ``x1`` object is an instance of SymPy's Symbol_.


ChemPy can also balance reactions where the reacting species are more complex and
are better described in other terms than their molecular formula. A silly, yet
illustrative example would be how to make pancakes without any partially used packages:

.. code:: python

>>> substances = {s.name: s for s in [
... Substance('pancake', composition=dict(eggs=1, spoons_of_flour=2, cups_of_milk=1)),
... Substance('eggs_6pack', composition=dict(eggs=6)),
... Substance('milk_carton', composition=dict(cups_of_milk=4)),
... Substance('flour_bag', composition=dict(spoons_of_flour=60))
... ]}
>>> pprint([dict(_) for _ in balance_stoichiometry({'eggs_6pack', 'milk_carton', 'flour_bag'},
... {'pancake'}, substances=substances)])
[{'eggs_6pack': 10, 'flour_bag': 2, 'milk_carton': 15}, {'pancake': 60}]

.. _Symbol: http://docs.sympy.org/latest/modules/core.html#sympy.core.symbol.Symbol


Balancing reactions
~~~~~~~~~~~~~~~~~~~
.. code:: python

>>> from chempy import Equilibrium
>>> from sympy import symbols
>>> K1, K2, Kw = symbols('K1 K2 Kw')
>>> e1 = Equilibrium({'MnO4-': 1, 'H+': 8, 'e-': 5}, {'Mn+2': 1, 'H2O': 4}, K1)
>>> e2 = Equilibrium({'O2': 1, 'H2O': 2, 'e-': 4}, {'OH-': 4}, K2)
>>> coeff = Equilibrium.eliminate([e1, e2], 'e-')
>>> coeff
[4, -5]
>>> redox = e1*coeff[0] + e2*coeff[1]
>>> print(redox)
32 H+ + 4 MnO4- + 20 OH- = 26 H2O + 4 Mn+2 + 5 O2; K1**4/K2**5
>>> autoprot = Equilibrium({'H2O': 1}, {'H+': 1, 'OH-': 1}, Kw)
>>> n = redox.cancel(autoprot)
>>> n
20
>>> redox2 = redox + n*autoprot
>>> print(redox2)
12 H+ + 4 MnO4- = 6 H2O + 4 Mn+2 + 5 O2; K1**4*Kw**20/K2**5

Working with units
~~~~~~~~~~~~~~~~~~
Functions and objects useful
for working with units are available from the ``chempy.units`` module. Here is an
example of how ChemPy can check consistency of units:

.. code:: python

>>> from chempy import Reaction
>>> r = Reaction.from_string("H2O -> H+ + OH-; 1e-4/M/s")
Traceback (most recent call last):
...
ValueError: Check failed: 'consistent_units'
>>> r = Reaction.from_string("H2O -> H+ + OH-; 1e-4/s")
>>> from chempy.units import to_unitless, default_units as u
>>> to_unitless(r.param, 1/u.minute)
0.006

right now the ``.units`` module wraps the quantities_ package with some minor
additions and work-arounds. However, there is no guarantee that the underlying
package will not change in a future version of ChemPy (there are many packages
for dealing with units in the scientific Python ecosystem).

.. _quantities: http://python-quantities.readthedocs.io/en/latest/


Chemical equilibria
~~~~~~~~~~~~~~~~~~~
.. code:: python

>>> from chempy import Equilibrium
>>> from chempy.chemistry import Species
>>> water_autop = Equilibrium({'H2O'}, {'H+', 'OH-'}, 10**-14) # unit "molar" assumed
>>> ammonia_prot = Equilibrium({'NH4+'}, {'NH3', 'H+'}, 10**-9.24) # same here
>>> from chempy.equilibria import EqSystem
>>> substances = [Species.from_formula(f) for f in 'H2O OH- H+ NH3 NH4+'.split()]
>>> eqsys = EqSystem([water_autop, ammonia_prot], substances)
>>> print('\n'.join(map(str, eqsys.rxns))) # "rxns" short for "reactions"
H2O = H+ + OH-; 1e-14
NH4+ = H+ + NH3; 5.75e-10
>>> from collections import defaultdict
>>> init_conc = defaultdict(float, {'H2O': 1, 'NH3': 0.1})
>>> x, sol, sane = eqsys.root(init_conc)
>>> assert sol['success'] and sane
>>> print(', '.join('%.2g' % v for v in x))
1, 0.0013, 7.6e-12, 0.099, 0.0013


Concepts
~~~~~~~~~
ChemPy collects equations and utility functions for working with
concepts such as `ionic strength <https://en.wikipedia.org/wiki/Ionic_strength>`_:

.. code:: python

>>> from chempy.electrolytes import ionic_strength
>>> ionic_strength({'Fe+3': 0.050, 'ClO4-': 0.150}) == .3
True

note how ChemPy parsed the charges from the names of the substances. There are
also e.g. empirical equations and convenience classes for them available, e.g.:

.. code:: python

>>> from chempy.henry import Henry
>>> kH_O2 = Henry(1.2e-3, 1800, ref='carpenter_1966')
>>> print('%.1e' % kH_O2(298.15))
1.2e-03

to get more information about e.g. this class, you may can look at the API
`documentation <https://bjodah.github.io/chempy/latest/chempy.html#module-chempy.henry>`_ .


Chemical kinetics (system of ordinary differential equations)
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
A common task when modelling problems in chemistry is to investigate the time dependence
of a system. This branch of study is known as
`chemical kinetics <https://en.wikipedia.org/wiki/Chemical_kinetics>`_, and ChemPy has
some classes and functions for working with such problems:

.. code:: python

>>> from chempy import ReactionSystem # The rate constants below are arbitrary
>>> rsys = ReactionSystem.from_string("""2 Fe+2 + H2O2 -> 2 Fe+3 + 2 OH-; 42
... 2 Fe+3 + H2O2 -> 2 Fe+2 + O2 + 2 H+; 17
... H+ + OH- -> H2O; 1e10
... H2O -> H+ + OH-; 1e-4""") # "[H2O]" = 1.0 (actually 55.4 at RT)
>>> from chempy.kinetics.ode import get_odesys
>>> odesys, extra = get_odesys(rsys)
>>> from collections import defaultdict
>>> import numpy as np
>>> tout = sorted(np.concatenate((np.linspace(0, 23), np.logspace(-8, 1))))
>>> c0 = defaultdict(float, {'Fe+2': 0.05, 'H2O2': 0.1, 'H2O': 1.0, 'H+': 1e-2, 'OH-': 1e-12})
>>> result = odesys.integrate(tout, c0, atol=1e-12, rtol=1e-14)
>>> import matplotlib.pyplot as plt
>>> fig, axes = plt.subplots(1, 2, figsize=(12, 5))
>>> for ax in axes:
... _ = result.plot(names=[k for k in rsys.substances if k != 'H2O'], ax=ax)
... _ = ax.legend(loc='best', prop={'size': 9})
... _ = ax.set_xlabel('Time')
... _ = ax.set_ylabel('Concentration')
>>> _ = axes[1].set_ylim([1e-13, 1e-1])
>>> _ = axes[1].set_xscale('log')
>>> _ = axes[1].set_yscale('log')
>>> _ = fig.tight_layout()
>>> _ = fig.savefig('examples/kinetics.png', dpi=72)

.. image:: https://raw.githubusercontent.com/bjodah/chempy/master/examples/kinetics.png

Properties
~~~~~~~~~~
One of the fundamental tasks in science is the careful collection of data about the world
around us. ChemPy contains a growing collection of parametrizations from the scientific
litterature with relevance in chemistry. Here is how you use one of these formulations:

.. code:: python

>>> from chempy import Substance
>>> from chempy.properties.water_density_tanaka_2001 import water_density as rho
>>> from chempy.units import to_unitless, default_units as u
>>> water = Substance.from_formula('H2O')
>>> for T_C in (15, 25, 35):
... concentration_H2O = rho(T=(273.15 + T_C)*u.kelvin, units=u)/water.molar_mass(units=u)
... print('[H2O] = %.2f M (at %d °C)' % (to_unitless(concentration_H2O, u.molar), T_C))
...
[H2O] = 55.46 M (at 15 °C)
[H2O] = 55.35 M (at 25 °C)
[H2O] = 55.18 M (at 35 °C)


Run notebooks using binder
~~~~~~~~~~~~~~~~~~~~~~~~~~
Using only a web-browser (and an internet connection) it is possible to explore the
notebooks here: (by the courtesy of the people behind mybinder)

.. image:: http://mybinder.org/badge.svg
:target: https://mybinder.org/v2/gh/bjodah/chempy/0bd22289c734ab7ae908d1b03bdcaaa9463f5677?filepath=index.ipynb
:alt: Binder


Licensing
---------
The source code is Open Source and is released under the very permissive
`"simplified (2-clause) BSD license" <https://opensource.org/licenses/BSD-2-Clause>`_.
See `LICENSE <LICENSE>`_ for further details.

See also
--------
- `SymPy <https://github.com/sympy/sympy>`_
- `pyneqsys <https://github.com/bjodah/pyneqsys>`_
- `pyodesys <https://github.com/bjodah/pyodesys>`_
- `thermo <https://github.com/CalebBell/thermo>`_

Contributing
------------
Contributors are welcome to suggest improvements at https://github.com/bjodah/chempy
(see further details `here <CONTRIBUTING.rst>`_).


Author
------
Björn I. Dahlgren, contact:
- gmail address: bjodah
- kth.se address: bda

Project details


Download files

Download the file for your platform. If you're not sure which to choose, learn more about installing packages.

Filename, size & hash SHA256 hash help File type Python version Upload date
chempy-0.6.1.tar.gz (154.3 kB) Copy SHA256 hash SHA256 Source None

Supported by

Elastic Elastic Search Pingdom Pingdom Monitoring Google Google BigQuery Sentry Sentry Error logging AWS AWS Cloud computing DataDog DataDog Monitoring Fastly Fastly CDN SignalFx SignalFx Supporter DigiCert DigiCert EV certificate StatusPage StatusPage Status page