Merged commits from various branches

CHANGELOG.md

@@ -1,3 +1,36 @@
+# 5.1.0 (2024/MM/DD)
+
+- Added support for analyzing the peaks in DRT results by fitting skew normal distributions:
+  - Added an `analyze_peaks` method to the `DRTResult` class.
+  - Added a `DRTPeaks` class, which can be used to, e.g., calculate the area of a peak.
+- Added an implementation of the Loewner method for calculating the distribution of relaxation times:
+  - Added a `calculate_drt_lm` function.
+  - Added an `LMResult` class.
+  - Added `--model-order` and `--model-order-method` CLI arguments.
+  - Updated the `plot.mpl.plot_gamma` function to support plotting `LMResult` instances.
+- Added the ability to define constraints when fitting circuits:
+  - Added a `generate_fit_identifiers` function.
+  - Added a `FitIdentifiers` class.
+- Added support for plotting DRT results as gamma vs f.
+  - Added CLI argument for plotting DRT results as gamma vs f.
+- Added circuits for generating mock data (`CIRCUIT_13`, `CIRCUIT_13_INVALID`, `CIRCUIT_14`, and `CIRCUIT_14_INVALID`).
+- Updated the TR-RBF implementation to be based off of a newer version of pyDRTtools.
+- Updated plotting functions to support drawing smoother lines if the input result has a `circuit` property.
+- Updated an exception message to provide more information if a variable that is being fitted is out of bounds.
+- Updated documentation.
+- Updated the `generate_mock_data` function so that it attempts to cast arguments to the appropriate type if some other types of values (e.g., integers) are provided instead.
+- Updated the `circuit.registry.register_element` function to support a new `private` keyword argument.
+  - Updated the `circuit.registry.get_elements` function to not include by default `Element` classes that were registered with `private=True`.
+  - Updated the `KramersKronigRC` and `KramersKronigAdmittanceRC` classes to be registered with `private=True`.
+- Fixed a bug that caused methods such as `DRTResult.get_peaks` to miss peaks at the time constant extremes.
+- Fixed a bug that caused an element's parameters in `FitResult.to_parameters_dataframe` to not be in a sorted order.
+- Fixed the previously unimplemented `FitResult.get_parameters` method.
+- Fixed a bug that caused `FitResult.to_parameters_dataframe` to return negative standard errors when the fitted value was negative.
+- Fixed a bug that could cause `FitResult.to_parameters_dataframe` to divide by zero.
+- Fixed a bug where an exception would be raised when whitespace was included between keyword arguments when passing circuit identifiers or CDCs via the CLI (e.g., `<R(RC):noise=5e-2, log_max_f=4>` would previously raise an exception whereas `<R(RC):noise=5e-2,log_max_f=4>` would not).
+- Refactored some of the code.
+
+
 # 5.0.2 (2024/09/10)
 
 - Updated the documentation (e.g., various function and method docstrings have been updated). Notably, the functions related to handling progress updates are now included in the API documentation.

LICENSES/LICENSE-DRT-from-Loewner-framework.txt

@@ -0,0 +1,22 @@
+MIT License
+
+Copyright (c) 2023 projectsEECandDRI
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
+

LICENSES/LICENSE-pyDRTtools.txt

@@ -1,6 +1,6 @@
 MIT License
 
-Copyright (c) 2020 ciuccislab
+Copyright (c) 2023 ciuccislab
 
 Permission is hereby granted, free of charge, to any person obtaining a copy
 of this software and associated documentation files (the "Software"), to deal

LICENSES/README.md

@@ -8,6 +8,11 @@
 - License: GPLv3 or later
 - Dependency.
 
+# DRT-from-Loewner-framework
+- https://github.com/projectsEECandDRI/DRT-from-Loewner-framework
+- License: MIT
+- Ported to Python and modified.
+
 # DRT-python-code
 - https://github.com/akulikovsky/DRT-python-code
 - License: GPLv3 or later

docs/source/apidocs_drt.rst

@@ -8,7 +8,16 @@ A collection of functions and classes for calculating the distribution of relaxa
 Wrapper function and base class
 -------------------------------
 .. automodule:: pyimpspec
-   :members: calculate_drt, DRTResult
+   :members: calculate_drt
+
+.. automodule:: pyimpspec
+   :members: DRTResult
+
+
+Classes related to peak analysis
+--------------------------------
+.. automodule:: pyimpspec
+   :members: DRTPeaks, DRTPeak
 
 
 Method functions and classes
@@ -23,6 +32,15 @@ BHT method
    :members: BHTResult
 
 
+LM method
+~~~~~~~~~~
+.. automodule:: pyimpspec.analysis.drt
+   :members: calculate_drt_lm
+
+.. automodule:: pyimpspec.analysis.drt
+   :members: LMResult
+
+
 m(RQ)fit method
 ~~~~~~~~~~~~~~~
 .. automodule:: pyimpspec.analysis.drt

docs/source/apidocs_fitting.rst

@@ -9,12 +9,12 @@ A collection of functions and classes for fitting equivalent circuits to data se
 Function
 --------
 .. automodule:: pyimpspec
-   :members: fit_circuit
+   :members: fit_circuit, generate_fit_identifiers
 
 Classes
 -------
 .. automodule:: pyimpspec
-   :members: FitResult, FittedParameter
+   :members: FitResult, FittedParameter, FitIdentifiers
 
 .. raw:: latex
 

docs/source/guide_drt.rst

@@ -19,6 +19,7 @@ The supported methods
 Implementations based on the following approaches are included in pyimpspec:
 
 - **BHT**: The Bayesian Hilbert transform method (see `Liu et al. (2020) <https://doi.org/10.1016/j.electacta.2020.136864>`_) that was originally implemented in `DRTtools <https://github.com/ciuccislab/DRTtools>`_ and `pyDRTtools <https://github.com/ciuccislab/pyDRTtools>`_.
+- **LM**: The Loewner method (see `Mayo and Antoulas (2007) <https://doi.org/10.1016/j.laa.2007.03.008>`_, `Sorrentino et al. (2022) <http://doi.org/10.2139/ssrn.4217752>`_, `Rüther et al. (2023) <https://doi.org/10.3390/batteries9020132>`_, and `Sorrentino et al. (2023) <https://doi.org/10.1016/j.jpowsour.2023.233575>`_) that was implemented in `DRT-from-Loewner-framework <https://github.com/projectsEECandDRI/DRT-from-Loewner-framework>`_.
 - **m(RQ)fit**: The multi-(RQ)-fit method (see `Boukamp (2015) <https://doi.org/10.1016/j.electacta.2014.12.059>`_ and `Boukamp and Rolle (2017) <https://doi.org/10.1016/j.ssi.2016.10.009>`_).
 - **TR-NNLS**: Tikhonov regularization and non-negative least squares (see `Kulikovsky (2021) <https://doi.org/10.1149/1945-7111/abf508>`_) that was originally implemented in `DRT-python-code <https://github.com/akulikovsky/DRT-python-code>`_.
 - **TR-RBF**: Tikhonov regularization and radial basis function (or piecewise linear) discretization (see `Wan et al. (2015) <https://doi.org/10.1016/j.electacta.2015.09.097>`_, `Ciucci and Chen (2015) <https://doi.org/10.1016/j.electacta.2015.03.123>`_, `Effat and Ciucci (2017) <https://doi.org/10.1016/j.electacta.2017.07.050>`_, and `Maradesa et al. (2023) <https://doi.org/10.1149/1945-7111/acbca4>`_) that was originally implemented in DRTtools_ and pyDRTtools_.
@@ -51,10 +52,12 @@ Each method has its own function that can be used but there is also a wrapper fu
    ... )
    >>> from pyimpspec.analysis.drt import (
    ...   BHTResult,              # Result of the BHT method
+   ...   LMResult,               # Result of the Loewner method
    ...   MRQFitResult,           # Result of the m(RQ)fit method
    ...   TRNNLSResult,           # Result of the TR-NNLS method
    ...   TRRBFResult,            # Result of the TR-RBF method
    ...   calculate_drt_bht,      # BHT method
+   ...   calculate_drt_lm,       # Loewner method
    ...   calculate_drt_mrq_fit,  # m(RQ)fit method
    ...   calculate_drt_tr_nnls,  # TR-NNLS method
    ...   calculate_drt_tr_rbf,   # TR-RBF method
@@ -76,6 +79,7 @@ Below are some figures demonstrating the results of the methods listed above whe
    )
    from pyimpspec.analysis.drt import (
      calculate_drt_bht,
+     calculate_drt_lm,
      calculate_drt_mrq_fit,
      calculate_drt_tr_nnls,
      calculate_drt_tr_rbf,
@@ -87,11 +91,18 @@ Below are some figures demonstrating the results of the methods listed above whe
      ax.set_ylim(-100, 900)
 
    data = generate_mock_data("CIRCUIT_1", noise=5e-2, seed=42)[0]
+   figure, axes = mpl.plot_nyquist(data, colors=dict(impedance="black"), markers=dict(impedance="o"))
+   figure.tight_layout()
+
    drt = calculate_drt_bht(data)
    figure, axes = mpl.plot_gamma(drt)
    adjust_limits(axes[0])
    figure.tight_layout()
 
+   drt = calculate_drt_lm(data)
+   figure, axes = mpl.plot_gamma(drt)
+   figure.tight_layout()
+
    circuit = parse_cdc("R(RQ)(RQ)")
    fit = fit_circuit(circuit, data)
    drt = calculate_drt_mrq_fit(data, fit.circuit, fit=fit)
@@ -115,16 +126,216 @@ Below are some figures demonstrating the results of the methods listed above whe
     \clearpage
 
 
+.. note::
+
+   The Loewner method can provide an estimate for the high-frequency resistance, but the high-frequency (i.e., low time constant) resistive-capacitive peak may not always provide this value directly.
+
+
+.. plot::
+
+   from pyimpspec import (
+     generate_mock_data,
+     parse_cdc,
+   )
+   from pyimpspec.analysis.drt import calculate_drt_lm
+   from pyimpspec import mpl
+
+   cdc = "R{R=140}(R{R=230}C{C=1e-6})(R{R=576}C{C=1e-4})(R{R=150}L{L=4e1})"
+   circuit = parse_cdc(cdc)
+   drawing = circuit.to_drawing()
+   drawing.draw()
+
+   data = generate_mock_data(cdc, noise=0)[0]
+   drt = calculate_drt_lm(data)
+
+   figure, axes = mpl.plot_nyquist(data, colors=dict(impedance="black"), markers=dict(impedance="o"))
+   mpl.plot_nyquist(drt, colors=dict(impedance="red"), markers=dict(impedance="+"), figure=figure, axes=axes)
+   figure.tight_layout()
+
+   figure, axes = mpl.plot_gamma(drt)
+   figure.tight_layout()
+
+
+.. code::
+
+   |   tau, RC (s) |   gamma, RC (ohm) |   tau, RL (s) |   gamma, RL (ohm) |
+   |--------------:|------------------:|--------------:|------------------:|
+   |   8.00717e-12 |           289.999 |      0.266667 |           150.012 |
+   |   0.00023     |           230     |    nan        |           nan     |
+   |   0.0576      |           576.007 |    nan        |           nan     |
+
+
+This example shows that some additional processing of the obtained values may be necessary.
+In this case, the high-frequency resistive-capacitive peak does not directly tell us the value of :math:`R_1`.
+However, we can still obtain an estimated value: :math:`R_1 \approx 289.999\ \Omega - 150.012\ \Omega \approx 140\ \Omega`.
+Other estimated values such as capacitances and inductances can also be obtained:
+
+- :math:`R_2 \approx 230\ \Omega` and :math:`C_1 \approx 0.0023\ {\rm s} / 230\ \Omega \approx 10^{-6}\ {\rm F}`
+- :math:`R_3 \approx 576\ \Omega` and :math:`C_2 \approx 0.0576\ {\rm s} / 576.007\ \Omega \approx 10^{-4}\ {\rm F}`
+- :math:`R_4 \approx 150\ \Omega` and :math:`L_2 \approx 0.266667\ {\rm s} \times 150.012\ \Omega \approx 40\ {\rm H}`
+
+The above is a best-case scenario, but even with a bit of added noise (:math:`\sigma = 0.05\ \% \times |Z|`) one can still obtain decent estimates despite the additional peaks:
+
+- :math:`R1 \approx 289.925\ \Omega - 149.085\ \Omega \approx 141\ \Omega`.
+- :math:`R_2 \approx 230\ \Omega` and :math:`C_1 \approx 0.00230559\ {\rm s} / 230.233\ \Omega \approx 10^{-6}\ {\rm F}`
+- :math:`R_3 \approx 574\ \Omega` and :math:`C_2 \approx 0.0574022\ {\rm s} / 574.438\ \Omega \approx 10^{-4}\ {\rm F}`
+- :math:`R_4 \approx 149\ \Omega` and :math:`L_2 \approx 0.271608\ {\rm s} \times 149.085\ \Omega \approx 40\ {\rm H}`
+
+
+.. plot::
+
+   from pyimpspec import (
+     generate_mock_data,
+     parse_cdc,
+   )
+   from pyimpspec.analysis.drt import calculate_drt_lm
+   from pyimpspec import mpl
+
+   cdc = "R{R=140}(R{R=230}C{C=1e-6})(R{R=576}C{C=1e-4})(R{R=150}L{L=4e1})"
+   data = generate_mock_data(cdc, noise=5e-2, seed=42)[0]
+   drt = calculate_drt_lm(data)
+
+   figure, axes = mpl.plot_gamma(drt)
+   figure.tight_layout()
+
+
+.. code::
+
+   |   tau, RC (s) |   gamma, RC (ohm) |   tau, RL (s) |   gamma, RL (ohm) |
+   |--------------:|------------------:|--------------:|------------------:|
+   |   4.60645e-10 |      289.925      |   3.05543e-06 |        0.00891355 |
+   |   5.28003e-06 |        0.0198689  |   3.05543e-06 |        0.00891355 |
+   |   5.28003e-06 |        0.0198689  |   1.63833e-05 |        0.0126733  |
+   |   5.92308e-06 |        0.0270496  |   1.63833e-05 |        0.0126733  |
+   |   5.92308e-06 |        0.0270496  |   2.75326e-05 |        0.0099311  |
+   |   1.01729e-05 |        0.00475814 |   2.75326e-05 |        0.0099311  |
+   |   1.01729e-05 |        0.00475814 |   0.000109749 |        0.0314605  |
+   |   3.70669e-05 |        0.0147584  |   0.000109749 |        0.0314605  |
+   |   3.70669e-05 |        0.0147584  |   0.000166522 |        0.0195471  |
+   |   7.28352e-05 |        0.0196474  |   0.000166522 |        0.0195471  |
+   |   7.28352e-05 |        0.0196474  |   0.000272624 |        0.0393392  |
+   |   0.000230559 |      230.233      |   0.000272624 |        0.0393392  |
+   |   0.000771132 |        0.00893735 |   0.000512795 |        0.0399306  |
+   |   0.000771132 |        0.00893735 |   0.000512795 |        0.0399306  |
+   |   0.00105728  |        0.00459818 |   0.00789638  |        0.0116901  |
+   |   0.00105728  |        0.00459818 |   0.00789638  |        0.0116901  |
+   |   0.00185517  |        0.0331227  |   0.00956206  |        1.02824    |
+   |   0.00185517  |        0.0331227  |   0.271608    |      149.085      |
+   |   0.00370024  |        0.0500698  |   0.658958    |        0.0536769  |
+   |   0.00370024  |        0.0500698  |   0.658958    |        0.0536769  |
+   |   0.00546967  |        0.044606   |   1.45096     |        0.0433677  |
+   |   0.00546967  |        0.044606   |   1.45096     |        0.0433677  |
+   |   0.0126403   |        0.150534   | nan           |      nan          |
+   |   0.0126403   |        0.150534   | nan           |      nan          |
+   |   0.0285097   |        0.00708604 | nan           |      nan          |
+   |   0.0285097   |        0.00708604 | nan           |      nan          |
+   |   0.0455447   |        0.116258   | nan           |      nan          |
+   |   0.0455447   |        0.116258   | nan           |      nan          |
+   |   0.0574022   |      574.438      | nan           |      nan          |
+   |   0.0675327   |        0.0681987  | nan           |      nan          |
+   |   0.0675327   |        0.0681987  | nan           |      nan          |
+   |   0.131324    |        0.0287684  | nan           |      nan          |
+   |   0.131324    |        0.0287684  | nan           |      nan          |
+   |   0.228107    |        0.0382958  | nan           |      nan          |
+   |   0.228107    |        0.0382958  | nan           |      nan          |
+   |   0.276414    |        0.0018757  | nan           |      nan          |
+   |   0.276414    |        0.0018757  | nan           |      nan          |
+   |   1.00071     |        1.19876    | nan           |      nan          |
+
+
+.. raw:: latex
+
+    \clearpage
+
+
+Peak analysis
+-------------
+
+DRT results (aside from those obtained using the Loewner method) can be analyzed using skew normal distributions (see `Danzer (2019) <https://doi.org/10.3390/batteries5030053>`_ and `Plank et al. (2024) <https://doi.org/10.1016/j.jpowsour.2023.233845>`_).
+This facilitates the estimation of the polarization contribution of each peak based on the area of that peak.
+
+.. plot::
+
+   from pyimpspec import (
+     generate_mock_data,
+     parse_cdc,
+   )
+   from pyimpspec.analysis.drt import calculate_drt_tr_rbf
+   from pyimpspec import mpl
+
+   cdc = "R{R=100}(R{R=300}C{C=5e-6})(R{R=450}C{C=1e-5})"
+   circuit = parse_cdc(cdc)
+   drawing = circuit.to_drawing()
+   drawing.draw()
+
+   data = generate_mock_data(cdc, noise=5e-2, seed=42)[0]
+   drt = calculate_drt_tr_rbf(data)
+
+   figure, axes = mpl.plot_nyquist(data, colors=dict(impedance="black"), markers=dict(impedance="o"))
+   mpl.plot_nyquist(drt, colors=dict(impedance="red"), markers=dict(impedance="+"), figure=figure, axes=axes)
+   figure.tight_layout()
+
+   figure, axes = mpl.plot_gamma(drt)
+   figure.tight_layout()
+
+   peaks = drt.analyze_peaks()
+   figure, axes = mpl.plot_gamma(drt)
+   figure.tight_layout()
+
+   peaks.to_peaks_dataframe().to_markdown(index=False)
+
+
+.. code::
+
+   |     tau (s) |   gamma (ohm) |   R_peak (ohm) |
+   |------------:|--------------:|---------------:|
+   | 1.24844e-05 |   2.68972e-08 |    4.33529e-07 |
+   | 0.00155063  | 520.733       |  280.396       |
+   | 0.00419888  | 831.637       |  470.508       |
+   | 9.76195     |   0.364749    |    0.552689    |
+
+
+.. list-table:: True and estimated values of parallel RC elements.
+   :header-rows: 1
+
+   * - Component
+     - True value
+     - Estimated value
+   * - :math:`R_2`
+     - :math:`300\ \Omega`
+     - :math:`280\ \Omega`
+   * - :math:`C_1`
+     - :math:`5.0 \times 10^{-6}\ {\rm F}`
+     - :math:`5.5 \times 10^{-6}\ {\rm F}`
+   * - :math:`R_3`
+     - :math:`450\ \Omega`
+     - :math:`471\ \Omega`
+   * - :math:`C_2`
+     - :math:`1.0 \times 10^{-5}\ {\rm F}`
+     - :math:`8.9 \times 10^{-6}\ {\rm F}`
+
+
+.. raw:: latex
+
+    \clearpage
+
+
 References:
 
-- `Boukamp, B.A., 2015, Electrochim. Acta, 154, 35-46 <https://doi.org/10.1016/j.electacta.2014.12.059>`_
 - `Boukamp, B.A. and Rolle, A, 2017, Solid State Ionics, 302, 12-18 <https://doi.org/10.1016/j.ssi.2016.10.009>`_
+- `Boukamp, B.A., 2015, Electrochim. Acta, 154, 35-46 <https://doi.org/10.1016/j.electacta.2014.12.059>`_
 - `Ciucci, F. and Chen, C., 2015, Electrochim. Acta, 167, 439-454 <https://doi.org/10.1016/j.electacta.2015.03.123>`_
+- `Danzer, M.A., 2019, Batteries, 5, 53 <https://doi.org/10.3390/batteries5030053>`_
 - `Effat, M. B. and Ciucci, F., 2017, Electrochim. Acta, 247, 1117-1129 <https://doi.org/10.1016/j.electacta.2017.07.050>`_
 - `Kulikovsky, A., 2021, J. Electrochem. Soc., 168, 044512 <https://doi.org/10.1149/1945-7111/abf508>`_
 - `Liu, J., Wan, T. H., and Ciucci, F., 2020, Electrochim. Acta, 357, 136864 <https://doi.org/10.1016/j.electacta.2020.136864>`_
-- `Wan, T. H., Saccoccio, M., Chen, C., and Ciucci, F., 2015, Electrochim. Acta, 184, 483-499 <https://doi.org/10.1016/j.electacta.2015.09.097>`_
 - `Maradesa, A., Py, B., Wan, T.H., Effat, M.B., and Ciucci F., 2023, J. Electrochem. Soc, 170, 030502 <https://doi.org/10.1149/1945-7111/acbca4>`_
+- `Mayo, A.J. and Antoulas, A.C., 2007, Linear Algebra Its Appl. 425, 634–662 <https://doi.org/10.1016/j.laa.2007.03.008>`_,
+- `Plank, C., Rüther, T., Jahn, L., Schamel, M., Schmidt, J.P., Ciucci, F., and Danzer, M.A., 2024, Journal of Power Sources, 594, 233845 <https://doi.org/10.1016/j.jpowsour.2023.233845>`_
+- `Rüther, T., Gosea, I.V., Jahn, L., Antoulas, A.C., and Danzer, M.A., 2023, Batteries, 9, 132 <https://doi.org/10.3390/batteries9020132>`_
+- `Sorrentino, A., Patel, B., Gosea, I.V., Antoulas, A.C., and Vidaković-Koch, T., 2022, SSRN <http://doi.org/10.2139/ssrn.4217752>`_
+- `Sorrentino, A., Patel, B., Gosea, I.V., Antoulas, A.C., and Vidaković-Koch, T., 2023, J. Power Sources, 585, 223575 <https://doi.org/10.1016/j.jpowsour.2023.233575>`_
+- `Wan, T. H., Saccoccio, M., Chen, C., and Ciucci, F., 2015, Electrochim. Acta, 184, 483-499 <https://doi.org/10.1016/j.electacta.2015.09.097>`_
 
 .. raw:: latex