Skip to content

Add scipy.optimize.elementwise.find_root (method='chandrupatla') to bishop88 functions #2498

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 11 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 11 additions & 1 deletion docs/sphinx/source/whatsnew/v0.13.1.rst
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,16 @@ Bug fixes

Enhancements
~~~~~~~~~~~~

* Add ``method='chandrupatla'`` (faster than ``brentq`` and slower than ``newton``,
but convergence is guaranteed) as an option for
:py:func:`pvlib.pvsystem.singlediode`,
:py:func:`~pvlib.pvsystem.i_from_v`,
:py:func:`~pvlib.pvsystem.v_from_i`,
:py:func:`~pvlib.pvsystem.max_power_point`,
:py:func:`~pvlib.singlediode.bishop88_mpp`,
:py:func:`~pvlib.singlediode.bishop88_v_from_i`, and
:py:func:`~pvlib.singlediode.bishop88_i_from_v`. (:issue:`2497`, :pull:`2498`)


Documentation
~~~~~~~~~~~~~
Expand All @@ -44,3 +53,4 @@ Maintenance
Contributors
~~~~~~~~~~~~
* Elijah Passmore (:ghuser:`eljpsm`)
* Kevin Anderson (:ghuser:`kandersolar`)
28 changes: 22 additions & 6 deletions pvlib/pvsystem.py
Original file line number Diff line number Diff line change
Expand Up @@ -2498,7 +2498,11 @@

method : str, default 'lambertw'
Determines the method used to calculate points on the IV curve. The
options are ``'lambertw'``, ``'newton'``, or ``'brentq'``.
options are ``'lambertw'``, ``'newton'``, ``'brentq'``, or
``'chandrupatla'``.

.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.

Returns
-------
Expand Down Expand Up @@ -2630,7 +2634,11 @@
cells ``Ns`` and the builtin voltage ``Vbi`` of the intrinsic layer.
[V].
method : str
either ``'newton'`` or ``'brentq'``
either ``'newton'``, ``'brentq'``, or ``'chandrupatla'``.

.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.


Returns
-------
Expand Down Expand Up @@ -2713,8 +2721,12 @@
0 < nNsVth

method : str
Method to use: ``'lambertw'``, ``'newton'``, or ``'brentq'``. *Note*:
``'brentq'`` is limited to 1st quadrant only.
Method to use: ``'lambertw'``, ``'newton'``, ``'brentq'``, or
``'chandrupatla'``. *Note*: ``'brentq'`` is limited to non-negative current.

Check failure on line 2725 in pvlib/pvsystem.py

View workflow job for this annotation

GitHub Actions / flake8-linter

E501 line too long (84 > 79 characters)

.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.


Returns
-------
Expand Down Expand Up @@ -2795,8 +2807,12 @@
0 < nNsVth

method : str
Method to use: ``'lambertw'``, ``'newton'``, or ``'brentq'``. *Note*:
``'brentq'`` is limited to 1st quadrant only.
Method to use: ``'lambertw'``, ``'newton'``, ``'brentq'``, or
``'chandrupatla'``. *Note*: ``'brentq'`` is limited to non-negative current.

Check failure on line 2811 in pvlib/pvsystem.py

View workflow job for this annotation

GitHub Actions / flake8-linter

E501 line too long (84 > 79 characters)

.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.


Returns
-------
Expand Down
166 changes: 127 additions & 39 deletions pvlib/singlediode.py
Original file line number Diff line number Diff line change
Expand Up @@ -109,13 +109,13 @@
(a-Si) modules that is the product of the PV module number of series
cells :math:`N_{s}` and the builtin voltage :math:`V_{bi}` of the
intrinsic layer. [V].
breakdown_factor : float, default 0
breakdown_factor : numeric, default 0
fraction of ohmic current involved in avalanche breakdown :math:`a`.
Default of 0 excludes the reverse bias term from the model. [unitless]
breakdown_voltage : float, default -5.5
breakdown_voltage : numeric, default -5.5
reverse breakdown voltage of the photovoltaic junction :math:`V_{br}`
[V]
breakdown_exp : float, default 3.28
breakdown_exp : numeric, default 3.28
avalanche breakdown exponent :math:`m` [unitless]
gradients : bool
False returns only I, V, and P. True also returns gradients
Expand Down Expand Up @@ -162,12 +162,11 @@
# calculate temporary values to simplify calculations
v_star = diode_voltage / nNsVth # non-dimensional diode voltage
g_sh = 1.0 / resistance_shunt # conductance
if breakdown_factor > 0: # reverse bias is considered
Copy link
Member Author

@kandersolar kandersolar Jul 8, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Had to get rid of the if for vectorized evaluation in find_root. The math works out the same either way, it's just a bit of unnecessary computation when breakdown_factor is zero.

(see also #1821)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not sure now why that was there, maybe thinking to save few microseconds.

brk_term = 1 - diode_voltage / breakdown_voltage
brk_pwr = np.power(brk_term, -breakdown_exp)
i_breakdown = breakdown_factor * diode_voltage * g_sh * brk_pwr
else:
i_breakdown = 0.

brk_term = 1 - diode_voltage / breakdown_voltage
brk_pwr = np.power(brk_term, -breakdown_exp)
i_breakdown = breakdown_factor * diode_voltage * g_sh * brk_pwr

i = (photocurrent - saturation_current * np.expm1(v_star) # noqa: W503
- diode_voltage * g_sh - i_recomb - i_breakdown) # noqa: W503
v = diode_voltage - i * resistance_series
Expand All @@ -177,18 +176,14 @@
grad_i_recomb = np.where(is_recomb, i_recomb / v_recomb, 0)
grad_2i_recomb = np.where(is_recomb, 2 * grad_i_recomb / v_recomb, 0)
g_diode = saturation_current * np.exp(v_star) / nNsVth # conductance
if breakdown_factor > 0: # reverse bias is considered
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same here

brk_pwr_1 = np.power(brk_term, -breakdown_exp - 1)
brk_pwr_2 = np.power(brk_term, -breakdown_exp - 2)
brk_fctr = breakdown_factor * g_sh
grad_i_brk = brk_fctr * (brk_pwr + diode_voltage *
-breakdown_exp * brk_pwr_1)
grad2i_brk = (brk_fctr * -breakdown_exp # noqa: W503
* (2 * brk_pwr_1 + diode_voltage # noqa: W503
* (-breakdown_exp - 1) * brk_pwr_2)) # noqa: W503
else:
grad_i_brk = 0.
grad2i_brk = 0.
brk_pwr_1 = np.power(brk_term, -breakdown_exp - 1)
brk_pwr_2 = np.power(brk_term, -breakdown_exp - 2)
brk_fctr = breakdown_factor * g_sh
grad_i_brk = brk_fctr * (brk_pwr + diode_voltage *
-breakdown_exp * brk_pwr_1)
grad2i_brk = (brk_fctr * -breakdown_exp # noqa: W503
* (2 * brk_pwr_1 + diode_voltage # noqa: W503
* (-breakdown_exp - 1) * brk_pwr_2)) # noqa: W503
grad_i = -g_diode - g_sh - grad_i_recomb - grad_i_brk # di/dvd
grad_v = 1.0 - grad_i * resistance_series # dv/dvd
# dp/dv = d(iv)/dv = v * di/dv + i
Expand Down Expand Up @@ -247,12 +242,19 @@
breakdown_exp : float, default 3.28
avalanche breakdown exponent :math:`m` [unitless]
method : str, default 'newton'
Either ``'newton'`` or ``'brentq'``. ''method'' must be ``'newton'``
if ``breakdown_factor`` is not 0.
Either ``'newton'``, ``'brentq'``, or ``'chandrupatla'``.
''method'' must be ``'newton'`` if ``breakdown_factor`` is not 0.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why this requirement for 'newton'? The arguments passed when method='brentq' or method='chandrupatla' include all the reverse bias parameters.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I do not know the answer to this. I see in #948 you wrote this:

The Lambert W technique cannot solve the single diode equation with the term for current at reverse bias. The Bishop '88 methods can solve the equation. Tested with method='newton', not clear if method='brentq' is reliable for this application.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

My apologies for asking you to explain something I added 5 years ago :). I didn't pick up that all you did was reformat a bit of documentation.

I think that comment was added because the tests didn't confirm that 'brentq' works. I can't think of a reason why it wouldn't work. I suspect there was a vectorization issue that I didn't see how to fix at the time. It may be gone with removal of the if breakdown > 0: statement.


.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.

method_kwargs : dict, optional
Keyword arguments passed to root finder method. See
:py:func:`scipy:scipy.optimize.brentq` and
:py:func:`scipy:scipy.optimize.newton` parameters.
Keyword arguments passed to the root finder. For options, see:

* ``method='brentq'``: :py:func:`scipy:scipy.optimize.brentq`
* ``method='newton'``: :py:func:`scipy:scipy.optimize.newton`
* ``method='chandrupatla'``: :py:func:`scipy:scipy.optimize.elementwise.find_root`

``'full_output': True`` is allowed, and ``optimizer_output`` would be
returned. See examples section.

Expand Down Expand Up @@ -291,7 +293,7 @@
.. [1] "Computer simulation of the effects of electrical mismatches in
photovoltaic cell interconnection circuits" JW Bishop, Solar Cell (1988)
:doi:`10.1016/0379-6787(88)90059-2`
"""
""" # noqa: E501
# collect args
args = (photocurrent, saturation_current,
resistance_series, resistance_shunt, nNsVth, d2mutau, NsVbi,
Expand Down Expand Up @@ -333,6 +335,30 @@
vd = newton(func=lambda x, *a: fv(x, voltage, *a), x0=x0,
fprime=lambda x, *a: bishop88(x, *a, gradients=True)[4],
args=args, **method_kwargs)
elif method == 'chandrupatla':
try:
from scipy.optimize.elementwise import find_root
except ModuleNotFoundError as e:

Check warning on line 341 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L338-L341

Added lines #L338 - L341 were not covered by tests
# TODO remove this when our minimum scipy version is >=1.15
msg = (

Check warning on line 343 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L343

Added line #L343 was not covered by tests
"method='chandrupatla' requires scipy v1.15 or greater "
"(available for Python 3.10+). "
"Select another method, or update your version of scipy."
)
raise ImportError(msg) from e

Check warning on line 348 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L348

Added line #L348 was not covered by tests

voc_est = estimate_voc(photocurrent, saturation_current, nNsVth)
shape = _shape_of_max_size(voltage, voc_est)
vlo = np.zeros(shape)
vhi = np.full(shape, voc_est)
bounds = (vlo, vhi)
kwargs_trimmed = method_kwargs.copy()
kwargs_trimmed.pop("full_output", None) # not valid for find_root

Check warning on line 356 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L350-L356

Added lines #L350 - L356 were not covered by tests

result = find_root(fv, bounds, args=(voltage, *args), **kwargs_trimmed)
vd = result.x
if method_kwargs.get('full_output'):
vd = (vd, result) # mimic the other methods

Check warning on line 361 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L358-L361

Added lines #L358 - L361 were not covered by tests
else:
raise NotImplementedError("Method '%s' isn't implemented" % method)

Expand Down Expand Up @@ -388,12 +414,19 @@
breakdown_exp : float, default 3.28
avalanche breakdown exponent :math:`m` [unitless]
method : str, default 'newton'
Either ``'newton'`` or ``'brentq'``. ''method'' must be ``'newton'``
if ``breakdown_factor`` is not 0.
Either ``'newton'``, ``'brentq'``, or ``'chandrupatla'``.
''method'' must be ``'newton'`` if ``breakdown_factor`` is not 0.

.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.

method_kwargs : dict, optional
Keyword arguments passed to root finder method. See
:py:func:`scipy:scipy.optimize.brentq` and
:py:func:`scipy:scipy.optimize.newton` parameters.
Keyword arguments passed to the root finder. For options, see:

* ``method='brentq'``: :py:func:`scipy:scipy.optimize.brentq`
* ``method='newton'``: :py:func:`scipy:scipy.optimize.newton`
* ``method='chandrupatla'``: :py:func:`scipy:scipy.optimize.elementwise.find_root`

``'full_output': True`` is allowed, and ``optimizer_output`` would be
returned. See examples section.

Expand Down Expand Up @@ -432,7 +465,7 @@
.. [1] "Computer simulation of the effects of electrical mismatches in
photovoltaic cell interconnection circuits" JW Bishop, Solar Cell (1988)
:doi:`10.1016/0379-6787(88)90059-2`
"""
""" # noqa: E501
# collect args
args = (photocurrent, saturation_current,
resistance_series, resistance_shunt, nNsVth, d2mutau, NsVbi,
Expand Down Expand Up @@ -474,6 +507,29 @@
vd = newton(func=lambda x, *a: fi(x, current, *a), x0=x0,
fprime=lambda x, *a: bishop88(x, *a, gradients=True)[3],
args=args, **method_kwargs)
elif method == 'chandrupatla':
try:
from scipy.optimize.elementwise import find_root
except ModuleNotFoundError as e:

Check warning on line 513 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L510-L513

Added lines #L510 - L513 were not covered by tests
# TODO remove this when our minimum scipy version is >=1.15
msg = (

Check warning on line 515 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L515

Added line #L515 was not covered by tests
"method='chandrupatla' requires scipy v1.15 or greater "
"(available for Python 3.10+). "
"Select another method, or update your version of scipy."
)
raise ImportError(msg) from e

Check warning on line 520 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L520

Added line #L520 was not covered by tests

shape = _shape_of_max_size(current, voc_est)
vlo = np.zeros(shape)
vhi = np.full(shape, voc_est)
bounds = (vlo, vhi)
kwargs_trimmed = method_kwargs.copy()
kwargs_trimmed.pop("full_output", None) # not valid for find_root

Check warning on line 527 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L522-L527

Added lines #L522 - L527 were not covered by tests

result = find_root(fi, bounds, args=(current, *args), **kwargs_trimmed)
vd = result.x
if method_kwargs.get('full_output'):
vd = (vd, result) # mimic the other methods

Check warning on line 532 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L529-L532

Added lines #L529 - L532 were not covered by tests
else:
raise NotImplementedError("Method '%s' isn't implemented" % method)

Expand Down Expand Up @@ -526,12 +582,19 @@
breakdown_exp : numeric, default 3.28
avalanche breakdown exponent :math:`m` [unitless]
method : str, default 'newton'
Either ``'newton'`` or ``'brentq'``. ''method'' must be ``'newton'``
if ``breakdown_factor`` is not 0.
Either ``'newton'``, ``'brentq'``, or ``'chandrupatla'``.
''method'' must be ``'newton'`` if ``breakdown_factor`` is not 0.

.. note::
``'chandrupatla'`` requires scipy 1.15 or greater.

method_kwargs : dict, optional
Keyword arguments passed to root finder method. See
:py:func:`scipy:scipy.optimize.brentq` and
:py:func:`scipy:scipy.optimize.newton` parameters.
Keyword arguments passed to the root finder. For options, see:

* ``method='brentq'``: :py:func:`scipy:scipy.optimize.brentq`
* ``method='newton'``: :py:func:`scipy:scipy.optimize.newton`
* ``method='chandrupatla'``: :py:func:`scipy:scipy.optimize.elementwise.find_root`

``'full_output': True`` is allowed, and ``optimizer_output`` would be
returned. See examples section.

Expand Down Expand Up @@ -571,7 +634,7 @@
.. [1] "Computer simulation of the effects of electrical mismatches in
photovoltaic cell interconnection circuits" JW Bishop, Solar Cell (1988)
:doi:`10.1016/0379-6787(88)90059-2`
"""
""" # noqa: E501
# collect args
args = (photocurrent, saturation_current,
resistance_series, resistance_shunt, nNsVth, d2mutau, NsVbi,
Expand Down Expand Up @@ -611,6 +674,31 @@
vd = newton(func=fmpp, x0=x0,
fprime=lambda x, *a: bishop88(x, *a, gradients=True)[7],
args=args, **method_kwargs)
elif method == 'chandrupatla':
try:
from scipy.optimize.elementwise import find_root
except ModuleNotFoundError as e:

Check warning on line 680 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L677-L680

Added lines #L677 - L680 were not covered by tests
# TODO remove this when our minimum scipy version is >=1.15
msg = (

Check warning on line 682 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L682

Added line #L682 was not covered by tests
"method='chandrupatla' requires scipy v1.15 or greater "
"(available for Python 3.10+). "
"Select another method, or update your version of scipy."
)
raise ImportError(msg) from e

Check warning on line 687 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L687

Added line #L687 was not covered by tests

vlo = np.zeros_like(photocurrent)
vhi = np.full_like(photocurrent, voc_est)
kwargs_trimmed = method_kwargs.copy()
kwargs_trimmed.pop("full_output", None) # not valid for find_root

Check warning on line 692 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L689-L692

Added lines #L689 - L692 were not covered by tests

result = find_root(fmpp,

Check warning on line 694 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L694

Added line #L694 was not covered by tests
(vlo, vhi),
args=args,
**kwargs_trimmed)
vd = result.x
if method_kwargs.get('full_output'):
vd = (vd, result) # mimic the other methods

Check warning on line 700 in pvlib/singlediode.py

View check run for this annotation

Codecov / codecov/patch

pvlib/singlediode.py#L698-L700

Added lines #L698 - L700 were not covered by tests

else:
raise NotImplementedError("Method '%s' isn't implemented" % method)

Expand Down
13 changes: 13 additions & 0 deletions tests/conftest.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@

import pandas as pd
import os
import sys
from packaging.version import Version
import pytest
from functools import wraps
Expand Down Expand Up @@ -194,6 +195,18 @@ def has_spa_c():
reason="requires pandas>=2.0.0")


# single-diode equation functions have method=='chandrupatla', which relies
# on scipy.optimize.elementwise.find_root, which is only available in
# scipy>=1.15. That is only available for python 3.10 and above, so
# we need to skip those tests on python 3.9.
# TODO remove this when we drop support for python 3.9.
chandrupatla_available = sys.version_info >= (3, 10)
chandrupatla = pytest.param(
"chandrupatla", marks=pytest.mark.skipif(not chandrupatla_available,
reason="needs scipy 1.15")
)


@pytest.fixture()
def golden():
return Location(39.742476, -105.1786, 'America/Denver', 1830.14)
Expand Down
Loading
Loading