diff --git a/Cargo.toml b/Cargo.toml
index 94faf5c..4fa159c 100644
--- a/Cargo.toml
+++ b/Cargo.toml
@@ -17,9 +17,9 @@ exclude = ["/.github/*", "*.ipynb", "/docs"]
 rustdoc-args = [ "--html-in-header", "./docs-header.html" ]
 
 [dependencies]
-quantity = "0.4"
+quantity = "0.5"
 approx = "0.4"
-num-dual = { version = "0.4", features = ["ndarray"] }
+num-dual = { version = "0.5", features = ["linalg"] }
 ndarray = { version = "0.15", features = ["serde"] }
 num-traits = "0.2"
 thiserror = "1.0"
@@ -27,8 +27,8 @@ serde = { version = "1.0", features = ["derive"] }
 serde_json = "1.0"
 indexmap = "1.7"
 either = "1.6"
-numpy = { version = "0.15", optional = true }
-pyo3 = { version = "0.15", optional = true }
+numpy = { version = "0.16", optional = true }
+pyo3 = { version = "0.16", optional = true }
 
 [features]
 default = []
diff --git a/build_wheel/Cargo.toml b/build_wheel/Cargo.toml
index eaf3852..895c05c 100644
--- a/build_wheel/Cargo.toml
+++ b/build_wheel/Cargo.toml
@@ -10,5 +10,5 @@ crate-type = ["cdylib"]
 
 [dependencies]
 feos-core = { path = "..", features = ["python"] }
-pyo3 = { version = "0.15", features = ["extension-module", "abi3", "abi3-py36"] }
+pyo3 = { version = "0.16", features = ["extension-module", "abi3", "abi3-py37"] }
 
diff --git a/src/phase_equilibria/mod.rs b/src/phase_equilibria/mod.rs
index 27ae334..400eb42 100644
--- a/src/phase_equilibria/mod.rs
+++ b/src/phase_equilibria/mod.rs
@@ -18,6 +18,7 @@ pub use phase_diagram_pure::PhaseDiagramPure;
 
 /// Level of detail in the iteration output.
 #[derive(Copy, Clone, PartialOrd, PartialEq)]
+#[cfg_attr(feature = "python", pyo3::pyclass)]
 pub enum Verbosity {
     /// Do not print output.
     None,
diff --git a/src/phase_equilibria/vle_pure.rs b/src/phase_equilibria/vle_pure.rs
index 5624df8..ec886a3 100644
--- a/src/phase_equilibria/vle_pure.rs
+++ b/src/phase_equilibria/vle_pure.rs
@@ -260,7 +260,7 @@ impl<U: EosUnit, E: EquationOfState> PhaseEquilibrium<U, E, 2> {
         let density = 0.75 * eos.max_density(None)?;
         let liquid = State::new_nvt(eos, temperature, U::reference_moles() / density, &m)?;
         let z = liquid.compressibility(Contributions::Total);
-        let mu = liquid.chemical_potential(Contributions::Residual);
+        let mu = liquid.chemical_potential(Contributions::ResidualNvt);
         let p = temperature
             * density
             * U::gas_constant()
@@ -427,7 +427,8 @@ impl<U: EosUnit, E: EquationOfState> PhaseEquilibrium<U, E, 2> {
         (0..eos.components())
             .map(|i| {
                 let pure_eos = Rc::new(eos.subset(&[i]));
-                PhaseEquilibrium::pure_t(&pure_eos, temperature, None, SolverOptions::default()).ok()
+                PhaseEquilibrium::pure_t(&pure_eos, temperature, None, SolverOptions::default())
+                    .ok()
             })
             .collect()
     }
diff --git a/src/python/cubic.rs b/src/python/cubic.rs
index fac555b..89bd9f4 100644
--- a/src/python/cubic.rs
+++ b/src/python/cubic.rs
@@ -3,7 +3,6 @@ use crate::joback::JobackRecord;
 use crate::parameter::{IdentifierOption, Parameter, ParameterError, PureRecord};
 use crate::python::joback::PyJobackRecord;
 use crate::python::parameter::{PyBinaryRecord, PyChemicalRecord, PyIdentifier};
-use crate::python::{PyContributions, PyVerbosity};
 use crate::*;
 use numpy::convert::ToPyArray;
 use numpy::{PyArray1, PyArray2};
@@ -25,10 +24,7 @@ impl PyPengRobinsonRecord {
     fn new(tc: f64, pc: f64, acentric_factor: f64) -> Self {
         Self(PengRobinsonRecord::new(tc, pc, acentric_factor))
     }
-}
 
-#[pyproto]
-impl pyo3::class::basic::PyObjectProtocol for PyPengRobinsonRecord {
     fn __repr__(&self) -> PyResult<String> {
         Ok(self.0.to_string())
     }
@@ -113,8 +109,8 @@ impl_vle_state!(PengRobinson, PyPengRobinson);
 #[pymodule]
 pub fn cubic(_py: Python<'_>, m: &PyModule) -> PyResult<()> {
     m.add_class::<PyIdentifier>()?;
-    m.add_class::<PyVerbosity>()?;
-    m.add_class::<PyContributions>()?;
+    m.add_class::<Verbosity>()?;
+    m.add_class::<Contributions>()?;
     m.add_class::<PyChemicalRecord>()?;
     m.add_class::<PyJobackRecord>()?;
 
diff --git a/src/python/joback.rs b/src/python/joback.rs
index 23c1f9d..d2302b2 100644
--- a/src/python/joback.rs
+++ b/src/python/joback.rs
@@ -38,10 +38,7 @@ impl PyJobackRecord {
     fn new(a: f64, b: f64, c: f64, d: f64, e: f64) -> Self {
         Self(JobackRecord::new(a, b, c, d, e))
     }
-}
 
-#[pyproto]
-impl pyo3::class::basic::PyObjectProtocol for PyJobackRecord {
     fn __repr__(&self) -> PyResult<String> {
         Ok(self.0.to_string())
     }
diff --git a/src/python/mod.rs b/src/python/mod.rs
index 4a1cb76..d8b510c 100644
--- a/src/python/mod.rs
+++ b/src/python/mod.rs
@@ -1,8 +1,8 @@
-use crate::{Contributions, EosError, Verbosity};
+use crate::EosError;
 use pyo3::exceptions::PyRuntimeError;
 use pyo3::prelude::*;
 use pyo3::{wrap_pymodule, PyErr};
-use quantity::python::PyInit_quantity;
+use quantity::python::__PYO3_PYMODULE_DEF_QUANTITY;
 
 mod cubic;
 mod equation_of_state;
@@ -14,83 +14,8 @@ mod statehd;
 mod user_defined;
 mod utils;
 
-pub use cubic::PyInit_cubic;
-pub use user_defined::PyInit_user_defined;
-
-/// Helmholtz energy contributions to consider
-/// when computing a property.
-#[pyclass(name = "Contributions")]
-#[derive(Copy, Clone)]
-pub struct PyContributions(pub Contributions);
-
-#[pymethods]
-impl PyContributions {
-    /// Only compute ideal gas contribution.
-    #[classattr]
-    #[allow(non_snake_case)]
-    pub fn IdealGas() -> Self {
-        Self(Contributions::IdealGas)
-    }
-
-    /// Only compute residual contribution with respect
-    /// to an ideal gas contribution which is defined at
-    /// T, V, {n}.
-    ///
-    /// See also
-    /// --------
-    /// ResidualP: to use an ideal gas reference defined at T, p, {n}
-    #[classattr]
-    #[allow(non_snake_case)]
-    pub fn Residual() -> Self {
-        Self(Contributions::Residual)
-    }
-
-    /// Only compute residual contribution with respect
-    /// to an ideal gas contribution which is defined at
-    /// T, p, {n}.
-    ///
-    /// See also
-    /// --------
-    /// Residual: to use an ideal gas reference defined at T, V, {n}
-    #[classattr]
-    #[allow(non_snake_case)]
-    pub fn ResidualP() -> Self {
-        Self(Contributions::ResidualP)
-    }
-
-    /// Compute all contributions
-    ///
-    /// Note
-    /// ----
-    /// This is the default for most properties.
-    #[classattr]
-    #[allow(non_snake_case)]
-    pub fn Total() -> Self {
-        Self(Contributions::Total)
-    }
-}
-
-/// Verbosity levels for iterative solvers.
-#[pyclass(name = "Verbosity")]
-#[derive(Copy, Clone)]
-pub struct PyVerbosity(pub Verbosity);
-
-#[pymethods]
-impl PyVerbosity {
-    /// Print a status message at the end of the iteration.
-    #[classattr]
-    #[allow(non_snake_case)]
-    pub fn Result() -> Self {
-        Self(Verbosity::Result)
-    }
-
-    /// Print a detailed progress of the iteration.
-    #[classattr]
-    #[allow(non_snake_case)]
-    pub fn Iter() -> Self {
-        Self(Verbosity::Iter)
-    }
-}
+pub use cubic::__PYO3_PYMODULE_DEF_CUBIC;
+pub use user_defined::__PYO3_PYMODULE_DEF_USER_DEFINED;
 
 impl From<EosError> for PyErr {
     fn from(e: EosError) -> PyErr {
diff --git a/src/python/parameter.rs b/src/python/parameter.rs
index 810a84b..d894e04 100644
--- a/src/python/parameter.rs
+++ b/src/python/parameter.rs
@@ -111,10 +111,7 @@ impl PyIdentifier {
     fn set_formula(&mut self, formula: &str) {
         self.0.formula = Some(formula.to_string());
     }
-}
 
-#[pyproto]
-impl pyo3::class::basic::PyObjectProtocol for PyIdentifier {
     fn __repr__(&self) -> PyResult<String> {
         Ok(self.0.to_string())
     }
@@ -205,10 +202,7 @@ impl PyChemicalRecord {
                 .to_object(py),
         }
     }
-}
 
-#[pyproto]
-impl pyo3::class::basic::PyObjectProtocol for PyChemicalRecord {
     fn __repr__(&self) -> PyResult<String> {
         Ok(self.0.to_string())
     }
@@ -271,10 +265,7 @@ impl PyBinaryRecord {
     fn set_model_record(&mut self, model_record: f64) {
         self.0.model_record = model_record;
     }
-}
 
-#[pyproto]
-impl pyo3::class::basic::PyObjectProtocol for PyBinaryRecord {
     fn __repr__(&self) -> PyResult<String> {
         Ok(self.0.to_string())
     }
@@ -337,10 +328,7 @@ impl PyBinarySegmentRecord {
     fn set_model_record(&mut self, model_record: f64) {
         self.0.model_record = model_record;
     }
-}
 
-#[pyproto]
-impl pyo3::class::basic::PyObjectProtocol for PyBinarySegmentRecord {
     fn __repr__(&self) -> PyResult<String> {
         Ok(self.0.to_string())
     }
@@ -428,10 +416,7 @@ macro_rules! impl_pure_record {
             fn set_ideal_gas_record(&mut self, ideal_gas_record: $py_ideal_gas_record) {
                 self.0.ideal_gas_record = Some(ideal_gas_record.0);
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PyPureRecord {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
@@ -540,10 +525,7 @@ macro_rules! impl_segment_record {
             fn set_ideal_gas_record(&mut self, ideal_gas_record: $py_ideal_gas_record) {
                 self.0.ideal_gas_record = Some(ideal_gas_record.0);
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PySegmentRecord {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
diff --git a/src/python/phase_equilibria.rs b/src/python/phase_equilibria.rs
index bfd3972..1081273 100644
--- a/src/python/phase_equilibria.rs
+++ b/src/python/phase_equilibria.rs
@@ -43,13 +43,13 @@ macro_rules! impl_vle_state {
                 initial_state: Option<&PyPhaseEquilibrium>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 Ok(Self(PhaseEquilibrium::pure_t(
                     &eos.0,
                     temperature.into(),
                     initial_state.and_then(|s| Some(&s.0)),
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                 )?))
             }
 
@@ -88,13 +88,13 @@ macro_rules! impl_vle_state {
                 initial_state: Option<&PyPhaseEquilibrium>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 Ok(Self(PhaseEquilibrium::pure_p(
                     &eos.0,
                     pressure.into(),
                     initial_state.and_then(|s| Some(&s.0)),
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                 )?))
             }
 
@@ -141,7 +141,7 @@ macro_rules! impl_vle_state {
                 init_vle_state: Option<&PyPhaseEquilibrium>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
                 non_volatile_components: Option<Vec<usize>>,
             ) -> PyResult<Self> {
                 Ok(Self(PhaseEquilibrium::tp_flash(
@@ -150,7 +150,7 @@ macro_rules! impl_vle_state {
                     pressure.into(),
                     feed,
                     init_vle_state.and_then(|s| Some(&s.0)),
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(), non_volatile_components
+                    (max_iter, tol, verbosity).into(), non_volatile_components
                 )?))
             }
 
@@ -195,7 +195,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let x = vapor_molefracs.and_then(|m| Some(m.to_owned_array()));
                 Ok(Self(PhaseEquilibrium::bubble_point_tx(
@@ -205,8 +205,8 @@ macro_rules! impl_vle_state {
                     &liquid_molefracs.to_owned_array(),
                     x.as_ref(),
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into()
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into()
                     )
                 )?))
             }
@@ -252,7 +252,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let x = vapor_molefracs.and_then(|m| Some(m.to_owned_array()));
                 Ok(Self(PhaseEquilibrium::bubble_point_px(
@@ -262,8 +262,8 @@ macro_rules! impl_vle_state {
                     &liquid_molefracs.to_owned_array(),
                     x.as_ref(),
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into()
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into()
                     )
                 )?))
             }
@@ -309,7 +309,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let x = liquid_molefracs.and_then(|m| Some(m.to_owned_array()));
                 Ok(Self(PhaseEquilibrium::dew_point_tx(
@@ -319,8 +319,8 @@ macro_rules! impl_vle_state {
                     &vapor_molefracs.to_owned_array(),
                     x.as_ref(),
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into()
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into()
                     )
                 )?))
             }
@@ -366,7 +366,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let x = liquid_molefracs.and_then(|m| Some(m.to_owned_array()));
                 Ok(Self(PhaseEquilibrium::dew_point_px(
@@ -376,8 +376,8 @@ macro_rules! impl_vle_state {
                     &vapor_molefracs.to_owned_array(),
                     x.as_ref(),
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into()
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into()
                     )
                 )?))
             }
@@ -498,10 +498,7 @@ macro_rules! impl_vle_state {
             fn _repr_markdown_(&self) -> String {
                 self.0._repr_markdown_()
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PyPhaseEquilibrium {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
@@ -549,21 +546,21 @@ macro_rules! impl_vle_state {
                 x_init: (f64, f64),
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
                 max_iter_bd_inner: Option<usize>,
                 max_iter_bd_outer: Option<usize>,
                 tol_bd_inner: Option<f64>,
                 tol_bd_outer: Option<f64>,
-                verbosity_bd: Option<PyVerbosity>,
+                verbosity_bd: Option<Verbosity>,
             ) -> PyResult<PyThreePhaseEquilibrium> {
                 Ok(PyThreePhaseEquilibrium(PhaseEquilibrium::heteroazeotrope_t(
                     &eos.0,
                     temperature.into(),
                     x_init,
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                     (
-                        (max_iter_bd_inner, tol_bd_inner, verbosity_bd.map(|v| v.0)).into(),
-                        (max_iter_bd_outer, tol_bd_outer, verbosity_bd.map(|v| v.0)).into(),
+                        (max_iter_bd_inner, tol_bd_inner, verbosity_bd).into(),
+                        (max_iter_bd_outer, tol_bd_outer, verbosity_bd).into(),
                     )
                 )?))
             }
@@ -603,21 +600,21 @@ macro_rules! impl_vle_state {
                 x_init: (f64, f64),
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
                 max_iter_bd_inner: Option<usize>,
                 max_iter_bd_outer: Option<usize>,
                 tol_bd_inner: Option<f64>,
                 tol_bd_outer: Option<f64>,
-                verbosity_bd: Option<PyVerbosity>,
+                verbosity_bd: Option<Verbosity>,
             ) -> PyResult<PyThreePhaseEquilibrium> {
                 Ok(PyThreePhaseEquilibrium(PhaseEquilibrium::heteroazeotrope_p(
                     &eos.0,
                     pressure.into(),
                     x_init,
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                     (
-                        (max_iter_bd_inner, tol_bd_inner, verbosity_bd.map(|v| v.0)).into(),
-                        (max_iter_bd_outer, tol_bd_outer, verbosity_bd.map(|v| v.0)).into(),
+                        (max_iter_bd_inner, tol_bd_inner, verbosity_bd).into(),
+                        (max_iter_bd_outer, tol_bd_outer, verbosity_bd).into(),
                     )
                 )?))
             }
@@ -643,10 +640,7 @@ macro_rules! impl_vle_state {
             fn _repr_markdown_(&self) -> String {
                 self.0._repr_markdown_()
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PyThreePhaseEquilibrium {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
@@ -682,12 +676,12 @@ macro_rules! impl_vle_state {
                 init_vle_state: Option<&PyPhaseEquilibrium>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
                 non_volatile_components: Option<Vec<usize>>,
             ) -> PyResult<PyPhaseEquilibrium> {
                 Ok(PyPhaseEquilibrium(self.0.tp_flash(
                     init_vle_state.and_then(|s| Some(&s.0)),
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                     non_volatile_components
                 )?))
             }
@@ -731,14 +725,14 @@ macro_rules! impl_vle_state {
                 critical_temperature: Option<PySINumber>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let dia = PhaseDiagramPure::new(
                     &eos.0,
                     min_temperature.into(),
                     npoints,
                     critical_temperature.map(|t| t.into()),
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                 )?;
                 Ok(Self(dia))
             }
@@ -866,7 +860,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let dia = PhaseDiagramBinary::new_txy(
                     &eos.0,
@@ -874,8 +868,8 @@ macro_rules! impl_vle_state {
                     npoints,
                     x_lle,
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into(),
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into(),
                     )
                 )?;
                 Ok(Self(dia))
@@ -919,7 +913,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let dia = PhaseDiagramBinary::new_pxy(
                     &eos.0,
@@ -927,8 +921,8 @@ macro_rules! impl_vle_state {
                     npoints,
                     x_lle,
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into(),
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into(),
                     )
                 )?;
                 Ok(Self(dia))
@@ -1097,7 +1091,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let dia = PhaseDiagramHetero::new_txy(
                     &eos.0,
@@ -1107,8 +1101,8 @@ macro_rules! impl_vle_state {
                     npoints_vle,
                     npoints_lle,
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into(),
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into(),
                     )
                 )?;
                 Ok(Self(dia))
@@ -1159,7 +1153,7 @@ macro_rules! impl_vle_state {
                 max_iter_outer: Option<usize>,
                 tol_inner: Option<f64>,
                 tol_outer: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 let dia = PhaseDiagramHetero::new_pxy(
                     &eos.0,
@@ -1169,8 +1163,8 @@ macro_rules! impl_vle_state {
                     npoints_vle,
                     npoints_lle,
                     (
-                        (max_iter_inner, tol_inner, verbosity.map(|v| v.0)).into(),
-                        (max_iter_outer, tol_outer, verbosity.map(|v| v.0)).into(),
+                        (max_iter_inner, tol_inner, verbosity).into(),
+                        (max_iter_outer, tol_outer, verbosity).into(),
                     )
                 )?;
                 Ok(Self(dia))
diff --git a/src/python/state.rs b/src/python/state.rs
index 3e16bc7..5028d37 100644
--- a/src/python/state.rs
+++ b/src/python/state.rs
@@ -135,10 +135,10 @@ macro_rules! impl_state {
                 initial_temperature: Option<PySINumber>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Vec<Self>> {
                 let t = initial_temperature.and_then(|t0| Some(t0.into()));
-                let cp = State::critical_point_pure(&eos.0, t, (max_iter, tol, verbosity.map(|v| v.0)).into())?;
+                let cp = State::critical_point_pure(&eos.0, t, (max_iter, tol, verbosity).into())?;
                 Ok(cp.into_iter().map(Self).collect())
             }
 
@@ -172,13 +172,13 @@ macro_rules! impl_state {
                 initial_temperature: Option<PySINumber>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 Ok(PyState(State::critical_point(
                     &eos.0,
                     moles.as_deref(),
                     initial_temperature.map(|t| t.into()),
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                 )?))
             }
 
@@ -211,13 +211,13 @@ macro_rules! impl_state {
                 initial_molefracs: Option<[f64; 2]>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 Ok(PyState(State::critical_point_binary_t(
                     &eos.0,
                     temperature.into(),
                     initial_molefracs,
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                 )?))
             }
 
@@ -253,14 +253,14 @@ macro_rules! impl_state {
                 initial_molefracs: Option<[f64; 2]>,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Self> {
                 Ok(PyState(State::critical_point_binary_p(
                     &eos.0,
                     pressure.into(),
                     initial_temperature.map(|t| t.into()),
                     initial_molefracs,
-                    (max_iter, tol, verbosity.map(|v| v.0)).into(),
+                    (max_iter, tol, verbosity).into(),
                 )?))
             }
 
@@ -283,11 +283,11 @@ macro_rules! impl_state {
             fn stability_analysis(&self,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<Vec<Self>> {
                 Ok(self
                     .0
-                    .stability_analysis((max_iter, tol, verbosity.map(|v| v.0)).into())?
+                    .stability_analysis((max_iter, tol, verbosity).into())?
                     .into_iter()
                     .map(Self)
                     .collect())
@@ -312,9 +312,9 @@ macro_rules! impl_state {
             fn is_stable(&self,
                 max_iter: Option<usize>,
                 tol: Option<f64>,
-                verbosity: Option<PyVerbosity>,
+                verbosity: Option<Verbosity>,
             ) -> PyResult<bool> {
-                Ok(self.0.is_stable((max_iter, tol, verbosity.map(|v| v.0)).into())?)
+                Ok(self.0.is_stable((max_iter, tol, verbosity).into())?)
             }
 
             /// Return pressure.
@@ -328,10 +328,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
-            #[pyo3(text_signature = "($self, contributions=PyContributions::Total())")]
-            fn pressure(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.pressure(contributions.0))
+            #[args(contributions = "Contributions::Total")]
+            #[pyo3(text_signature = "($self, contributions)")]
+            fn pressure(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.pressure(contributions))
             }
 
             /// Return pressure contributions.
@@ -359,10 +359,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// float
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn compressibility(&self, contributions: PyContributions) -> f64 {
-                self.0.compressibility(contributions.0)
+            fn compressibility(&self, contributions: Contributions) -> f64 {
+                self.0.compressibility(contributions)
             }
 
             /// Return partial derivative of pressure w.r.t. volume.
@@ -376,10 +376,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dp_dv(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.dp_dv(contributions.0))
+            fn dp_dv(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.dp_dv(contributions))
             }
 
             /// Return partial derivative of pressure w.r.t. density.
@@ -393,10 +393,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dp_drho(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.dp_drho(contributions.0))
+            fn dp_drho(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.dp_drho(contributions))
             }
 
             /// Return partial derivative of pressure w.r.t. temperature.
@@ -410,10 +410,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dp_dt(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.dp_dt(contributions.0))
+            fn dp_dt(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.dp_dt(contributions))
             }
 
             /// Return partial derivative of pressure w.r.t. amount of substance.
@@ -427,10 +427,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray1
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dp_dni(&self, contributions: PyContributions) -> PySIArray1 {
-                PySIArray1::from(self.0.dp_dni(contributions.0))
+            fn dp_dni(&self, contributions: Contributions) -> PySIArray1 {
+                PySIArray1::from(self.0.dp_dni(contributions))
             }
 
             /// Return second partial derivative of pressure w.r.t. volume.
@@ -444,10 +444,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn d2p_dv2(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.d2p_dv2(contributions.0))
+            fn d2p_dv2(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.d2p_dv2(contributions))
             }
 
             /// Return second partial derivative of pressure w.r.t. density.
@@ -461,10 +461,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn d2p_drho2(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.d2p_drho2(contributions.0))
+            fn d2p_drho2(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.d2p_drho2(contributions))
             }
 
             /// Return molar volume of each component.
@@ -478,10 +478,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray1
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn molar_volume(&self, contributions: PyContributions) -> PySIArray1 {
-                PySIArray1::from(self.0.molar_volume(contributions.0))
+            fn molar_volume(&self, contributions: Contributions) -> PySIArray1 {
+                PySIArray1::from(self.0.molar_volume(contributions))
             }
 
             /// Return chemical potential of each component.
@@ -495,10 +495,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray1
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn chemical_potential(&self, contributions: PyContributions) -> PySIArray1 {
-                PySIArray1::from(self.0.chemical_potential(contributions.0))
+            fn chemical_potential(&self, contributions: Contributions) -> PySIArray1 {
+                PySIArray1::from(self.0.chemical_potential(contributions))
             }
 
             /// Return chemical potential contributions.
@@ -532,10 +532,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray1
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dmu_dt(&self, contributions: PyContributions) -> PySIArray1 {
-                PySIArray1::from(self.0.dmu_dt(contributions.0))
+            fn dmu_dt(&self, contributions: Contributions) -> PySIArray1 {
+                PySIArray1::from(self.0.dmu_dt(contributions))
             }
 
             /// Return derivative of chemical potential w.r.t amount of substance.
@@ -549,10 +549,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray2
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dmu_dni(&self, contributions: PyContributions) -> PySIArray2 {
-                PySIArray2::from(self.0.dmu_dni(contributions.0))
+            fn dmu_dni(&self, contributions: Contributions) -> PySIArray2 {
+                PySIArray2::from(self.0.dmu_dni(contributions))
             }
 
             /// Return logarithmic fugacity coefficient.
@@ -616,10 +616,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn c_v(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.c_v(contributions.0))
+            fn c_v(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.c_v(contributions))
             }
 
             /// Return derivative of isochoric heat capacity w.r.t. temperature.
@@ -633,10 +633,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn dc_v_dt(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.dc_v_dt(contributions.0))
+            fn dc_v_dt(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.dc_v_dt(contributions))
             }
 
             /// Return isobaric heat capacity.
@@ -650,10 +650,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn c_p(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.c_p(contributions.0))
+            fn c_p(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.c_p(contributions))
             }
 
 	        /// Return entropy.
@@ -667,10 +667,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn entropy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.entropy(contributions.0))
+            fn entropy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.entropy(contributions))
             }
 
             /// Return derivative of entropy with respect to temperature.
@@ -684,10 +684,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn ds_dt(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.ds_dt(contributions.0))
+            fn ds_dt(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.ds_dt(contributions))
             }
 
             /// Return molar entropy.
@@ -701,10 +701,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn molar_entropy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.molar_entropy(contributions.0))
+            fn molar_entropy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.molar_entropy(contributions))
             }
 
 
@@ -719,10 +719,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray1
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn partial_molar_entropy(&self, contributions: PyContributions) -> PySIArray1 {
-                PySIArray1::from(self.0.partial_molar_entropy(contributions.0))
+            fn partial_molar_entropy(&self, contributions: Contributions) -> PySIArray1 {
+                PySIArray1::from(self.0.partial_molar_entropy(contributions))
             }
 
             /// Return enthalpy.
@@ -736,10 +736,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn enthalpy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.enthalpy(contributions.0))
+            fn enthalpy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.enthalpy(contributions))
             }
 
             /// Return molar enthalpy.
@@ -753,10 +753,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn molar_enthalpy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.molar_enthalpy(contributions.0))
+            fn molar_enthalpy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.molar_enthalpy(contributions))
             }
 
 
@@ -771,10 +771,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SIArray1
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn partial_molar_enthalpy(&self, contributions: PyContributions) -> PySIArray1 {
-                PySIArray1::from(self.0.partial_molar_enthalpy(contributions.0))
+            fn partial_molar_enthalpy(&self, contributions: Contributions) -> PySIArray1 {
+                PySIArray1::from(self.0.partial_molar_enthalpy(contributions))
             }
 
             /// Return helmholtz_energy.
@@ -788,10 +788,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn helmholtz_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.helmholtz_energy(contributions.0))
+            fn helmholtz_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.helmholtz_energy(contributions))
             }
 
             /// Return molar helmholtz_energy.
@@ -805,10 +805,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn molar_helmholtz_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.molar_helmholtz_energy(contributions.0))
+            fn molar_helmholtz_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.molar_helmholtz_energy(contributions))
             }
 
             /// Return helmholtz energy contributions.
@@ -836,10 +836,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn gibbs_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.gibbs_energy(contributions.0))
+            fn gibbs_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.gibbs_energy(contributions))
             }
 
             /// Return molar gibbs_energy.
@@ -853,10 +853,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn molar_gibbs_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.molar_gibbs_energy(contributions.0))
+            fn molar_gibbs_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.molar_gibbs_energy(contributions))
             }
 
 
@@ -871,10 +871,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn internal_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.internal_energy(contributions.0))
+            fn internal_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.internal_energy(contributions))
             }
 
             /// Return molar internal_energy.
@@ -888,10 +888,10 @@ macro_rules! impl_state {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn molar_internal_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.molar_internal_energy(contributions.0))
+            fn molar_internal_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.molar_internal_energy(contributions))
             }
 
             /// Return Joule Thomson coefficient.
@@ -982,10 +982,7 @@ macro_rules! impl_state {
                     )
                 }
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PyState {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
@@ -1069,10 +1066,10 @@ macro_rules! impl_state_molarweight {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn specific_helmholtz_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.specific_helmholtz_energy(contributions.0))
+            fn specific_helmholtz_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.specific_helmholtz_energy(contributions))
             }
 
             /// Return mass specific entropy.
@@ -1086,10 +1083,10 @@ macro_rules! impl_state_molarweight {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn specific_entropy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.specific_entropy(contributions.0))
+            fn specific_entropy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.specific_entropy(contributions))
             }
 
             /// Return mass specific internal_energy.
@@ -1103,10 +1100,10 @@ macro_rules! impl_state_molarweight {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn specific_internal_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.specific_internal_energy(contributions.0))
+            fn specific_internal_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.specific_internal_energy(contributions))
             }
 
             /// Return mass specific gibbs_energy.
@@ -1120,10 +1117,10 @@ macro_rules! impl_state_molarweight {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn specific_gibbs_energy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.specific_gibbs_energy(contributions.0))
+            fn specific_gibbs_energy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.specific_gibbs_energy(contributions))
             }
 
             /// Return mass specific enthalpy.
@@ -1137,10 +1134,10 @@ macro_rules! impl_state_molarweight {
             /// Returns
             /// -------
             /// SINumber
-            #[args(contributions = "PyContributions::Total()")]
+            #[args(contributions = "Contributions::Total")]
             #[pyo3(text_signature = "($self, contributions)")]
-            fn specific_enthalpy(&self, contributions: PyContributions) -> PySINumber {
-                PySINumber::from(self.0.specific_enthalpy(contributions.0))
+            fn specific_enthalpy(&self, contributions: Contributions) -> PySINumber {
+                PySINumber::from(self.0.specific_enthalpy(contributions))
             }
         }
     };
diff --git a/src/python/user_defined.rs b/src/python/user_defined.rs
index d773635..c7f8dd5 100644
--- a/src/python/user_defined.rs
+++ b/src/python/user_defined.rs
@@ -1,4 +1,4 @@
-use crate::python::{statehd::*, PyContributions, PyVerbosity};
+use crate::python::statehd::*;
 use crate::*;
 use ndarray::prelude::*;
 use num_dual::python::{PyDual3Dual64, PyDual3_64, PyDual64, PyHyperDual64, PyHyperDualDual64};
@@ -272,8 +272,8 @@ pub fn user_defined(_py: Python, m: &PyModule) -> PyResult<()> {
     m.add_class::<PyStateD3D>()?;
     m.add_class::<PyStateD3DV2>()?;
     m.add_class::<PyStateD3DV3>()?;
-    m.add_class::<PyVerbosity>()?;
-    m.add_class::<PyContributions>()?;
+    m.add_class::<Verbosity>()?;
+    m.add_class::<Contributions>()?;
     m.add_class::<PyUserDefinedEos>()?;
     m.add_class::<PyState>()?;
     m.add_class::<PyPhaseEquilibrium>()?;
diff --git a/src/python/utils.rs b/src/python/utils.rs
index e1b7811..6a5d975 100644
--- a/src/python/utils.rs
+++ b/src/python/utils.rs
@@ -236,10 +236,7 @@ macro_rules! impl_estimator {
             fn get_datapoints(&self) -> usize {
                 self.0.datapoints()
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PyDataSet {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
@@ -385,10 +382,7 @@ macro_rules! impl_estimator {
             fn _repr_markdown_(&self) -> String {
                 self.0._repr_markdownn_()
             }
-        }
 
-        #[pyproto]
-        impl pyo3::class::basic::PyObjectProtocol for PyEstimator {
             fn __repr__(&self) -> PyResult<String> {
                 Ok(self.0.to_string())
             }
diff --git a/src/state/properties.rs b/src/state/properties.rs
index e853f12..1da2814 100644
--- a/src/state/properties.rs
+++ b/src/state/properties.rs
@@ -17,13 +17,14 @@ pub(crate) enum Evaluate {
 
 /// Possible contributions that can be computed.
 #[derive(Clone, Copy)]
+#[cfg_attr(feature = "python", pyo3::pyclass)]
 pub enum Contributions {
     /// Only compute the ideal gas contribution
     IdealGas,
     /// Only compute the difference between the total and the ideal gas contribution
-    Residual,
+    ResidualNvt,
     /// Compute the differnce between the total and the ideal gas contribution for a (N,p,T) reference state
-    ResidualP,
+    ResidualNpt,
     /// Compute ideal gas and residual contributions
     Total,
 }
@@ -143,14 +144,14 @@ impl<U: EosUnit, E: EquationOfState> State<U, E> {
         match contributions {
             Contributions::IdealGas => f(self, Evaluate::IdealGas),
             Contributions::Total => f(self, Evaluate::Total),
-            Contributions::Residual => {
+            Contributions::ResidualNvt => {
                 if additive {
                     f(self, Evaluate::Residual)
                 } else {
                     f(self, Evaluate::Total) - f(self, Evaluate::IdealGas)
                 }
             }
-            Contributions::ResidualP => {
+            Contributions::ResidualNpt => {
                 let p = self.pressure_(Evaluate::Total);
                 let state_p = Self::new_nvt(
                     &self.eos,
@@ -292,7 +293,8 @@ impl<U: EosUnit, E: EquationOfState> State<U, E> {
 
     /// Logarithm of the fugacity coefficient: $\ln\varphi_i=\beta\mu_i^\mathrm{res}\left(T,p,\lbrace N_i\rbrace\right)$
     pub fn ln_phi(&self) -> Array1<f64> {
-        (self.chemical_potential(Contributions::ResidualP) / (U::gas_constant() * self.temperature))
+        (self.chemical_potential(Contributions::ResidualNpt)
+            / (U::gas_constant() * self.temperature))
             .into_value()
             .unwrap()
     }
@@ -304,18 +306,18 @@ impl<U: EosUnit, E: EquationOfState> State<U, E> {
                 - s.chemical_potential_(evaluate) / self.temperature)
                 / (U::gas_constant() * self.temperature)
         };
-        self.evaluate_property(func, Contributions::ResidualP, false)
+        self.evaluate_property(func, Contributions::ResidualNpt, false)
     }
 
     /// Partial derivative of the logarithm of the fugacity coefficient w.r.t. pressure: $\left(\frac{\partial\ln\varphi_i}{\partial p}\right)_{T,N_i}$
     pub fn dln_phi_dp(&self) -> QuantityArray1<U> {
-        self.molar_volume(Contributions::ResidualP) / (U::gas_constant() * self.temperature)
+        self.molar_volume(Contributions::ResidualNpt) / (U::gas_constant() * self.temperature)
     }
 
     /// Partial derivative of the logarithm of the fugacity coefficient w.r.t. moles: $\left(\frac{\partial\ln\varphi_i}{\partial N_j}\right)_{T,p,N_k}$
     pub fn dln_phi_dnj(&self) -> QuantityArray2<U> {
         let n = self.eos.components();
-        let dmu_dni = self.dmu_dni(Contributions::Residual);
+        let dmu_dni = self.dmu_dni(Contributions::ResidualNvt);
         let dp_dni = self.dp_dni(Contributions::Total);
         let dp_dv = self.dp_dv(Contributions::Total);
         let dp_dn_2 = QuantityArray::from_shape_fn((n, n), |(i, j)| dp_dni.get(i) * dp_dni.get(j));
@@ -619,7 +621,7 @@ impl<U: EosUnit, E: EquationOfState + EntropyScaling<U>> State<U, E> {
     /// Return the viscosity via entropy scaling.
     pub fn viscosity(&self) -> EosResult<QuantityScalar<U>> {
         let s = self
-            .molar_entropy(Contributions::Residual)
+            .molar_entropy(Contributions::ResidualNvt)
             .to_reduced(U::reference_molar_entropy())?;
         Ok(self
             .eos
@@ -633,7 +635,7 @@ impl<U: EosUnit, E: EquationOfState + EntropyScaling<U>> State<U, E> {
     /// that is used for entropy scaling.
     pub fn ln_viscosity_reduced(&self) -> EosResult<f64> {
         let s = self
-            .molar_entropy(Contributions::Residual)
+            .molar_entropy(Contributions::ResidualNvt)
             .to_reduced(U::reference_molar_entropy())?;
         self.eos.viscosity_correlation(s, &self.molefracs)
     }
@@ -647,7 +649,7 @@ impl<U: EosUnit, E: EquationOfState + EntropyScaling<U>> State<U, E> {
     /// Return the diffusion via entropy scaling.
     pub fn diffusion(&self) -> EosResult<QuantityScalar<U>> {
         let s = self
-            .molar_entropy(Contributions::Residual)
+            .molar_entropy(Contributions::ResidualNvt)
             .to_reduced(U::reference_molar_entropy())?;
         Ok(self
             .eos
@@ -661,7 +663,7 @@ impl<U: EosUnit, E: EquationOfState + EntropyScaling<U>> State<U, E> {
     /// that is used for entropy scaling.
     pub fn ln_diffusion_reduced(&self) -> EosResult<f64> {
         let s = self
-            .molar_entropy(Contributions::Residual)
+            .molar_entropy(Contributions::ResidualNvt)
             .to_reduced(U::reference_molar_entropy())?;
         self.eos.diffusion_correlation(s, &self.molefracs)
     }
@@ -675,7 +677,7 @@ impl<U: EosUnit, E: EquationOfState + EntropyScaling<U>> State<U, E> {
     /// Return the thermal conductivity via entropy scaling.
     pub fn thermal_conductivity(&self) -> EosResult<QuantityScalar<U>> {
         let s = self
-            .molar_entropy(Contributions::Residual)
+            .molar_entropy(Contributions::ResidualNvt)
             .to_reduced(U::reference_molar_entropy())?;
         Ok(self
             .eos
@@ -692,7 +694,7 @@ impl<U: EosUnit, E: EquationOfState + EntropyScaling<U>> State<U, E> {
     /// that is used for entropy scaling.
     pub fn ln_thermal_conductivity_reduced(&self) -> EosResult<f64> {
         let s = self
-            .molar_entropy(Contributions::Residual)
+            .molar_entropy(Contributions::ResidualNvt)
             .to_reduced(U::reference_molar_entropy())?;
         self.eos
             .thermal_conductivity_correlation(s, &self.molefracs)