Allow users to interact with root finders parameters (#1764)

* Change of use in bishop88_i_from_v

* What the warnings test should look like

* Apply Adriesse's suggestion

Co-Authored-By: Anton Driesse <[email protected]>

* Update test_singlediode.py

* Minimun runnable test for pass conditions

* Cleanup import

* Add functional test for success and failing calls with kwargs

* Add documentation (minimal, only bishop88_i_from_v)

* Try to fix scipy links (I)

* It worked 🎉🎉🎉

* Move default newton values to _prepare_newton_inputs

* Modify bishop88_v_from_i, bishop88_mpp; update tests and docs

`maxiter` had to be changed so solutions converge near target value correctly.

* Typo in kwargs type error

Co-Authored-By: Cliff Hansen <[email protected]>

* Add what's new entry

* Delete fprime2 default newton param

Co-Authored-By: Kevin Anderson <[email protected]>

* Rename kwargs to method_kwargs

Co-Authored-By: Anton Driesse <[email protected]>

* Adriesse's suggestion to make full output available

Missing tests

Co-Authored-By: Anton Driesse <[email protected]>

* Update tests from Kevin suggestion

Co-Authored-By: Kevin Anderson <[email protected]>

* I won't miss Py3.7

* Uhm, yeah, don't like this solution but anyways...

* Tests should be fixed now

* full_output=True tests

* Update docs

* Now these are tests

* Bruh should have run locally

* Docs: add use cases and document as `**method_kwargs :`

* Tests: also test expected length of the first value

In case some unpacking went wrong

* Docs: doctest examples - make functional

* Behaviour & docs: Change how method_kwargs are passed, enhance docs

* omg who would have expected the mutable default parameter bug

* Apply suggestions from code review

Co-authored by @kandersolar

Co-authored-by: Kevin Anderson <[email protected]>

---------

Co-authored-by: Anton Driesse <[email protected]>
Co-authored-by: Cliff Hansen <[email protected]>
Co-authored-by: Kevin Anderson <[email protected]>
Co-authored-by: Kevin Anderson <[email protected]>

Author: 4 people

docs/sphinx/source/whatsnew/v0.10.0.rst

@@ -31,6 +31,13 @@ Enhancements
 ~~~~~~~~~~~~
 * Added `map_variables` parameter to :py:func:`pvlib.iotools.read_srml`
   and :py:func:`pvlib.iotools.read_srml_month_from_solardat` (:pull:`1773`)
+* Allow passing keyword arguments to :py:func:`scipy:scipy.optimize.brentq` and
+  :py:func:`scipy:scipy.optimize.newton` solvers in
+  :py:func:`~pvlib.singlediode.bishop88_mpp`,
+  :py:func:`~pvlib.singlediode.bishop88_i_from_v` and
+  :py:func:`~pvlib.singlediode.bishop88_v_from_i`. Among others,
+  tolerance and number of iterations can be set.
+  (:issue:`1249`, :pull:`1764`)
 * Improved `ModelChainResult.__repr__` (:pull:`1236`)
 
 
@@ -57,5 +64,5 @@ Contributors
 ~~~~~~~~~~~~
 * Taos Transue (:ghuser:`reepoi`)
 * Adam R. Jensen (:ghuser:`AdamRJensen`)
+* Echedey Luis (:ghuser:`echedey-ls`)
 * Cliff Hansen (:ghuser:`cwhanse`)
-

pvlib/singlediode.py

@@ -2,15 +2,17 @@
 Low-level functions for solving the single diode equation.
 """
 
-from functools import partial
 import numpy as np
 from pvlib.tools import _golden_sect_DataFrame
 
 from scipy.optimize import brentq, newton
 from scipy.special import lambertw
 
-# set keyword arguments for all uses of newton in this module
-newton = partial(newton, tol=1e-6, maxiter=100, fprime2=None)
+# newton method default parameters for this module
+NEWTON_DEFAULT_PARAMS = {
+    'tol': 1e-6,
+    'maxiter': 100
+}
 
 # intrinsic voltage per cell junction for a:Si, CdTe, Mertens et al.
 VOLTAGE_BUILTIN = 0.9  # [V]
@@ -206,7 +208,7 @@ def bishop88_i_from_v(voltage, photocurrent, saturation_current,
                       resistance_series, resistance_shunt, nNsVth,
                       d2mutau=0, NsVbi=np.Inf, breakdown_factor=0.,
                       breakdown_voltage=-5.5, breakdown_exp=3.28,
-                      method='newton'):
+                      method='newton', method_kwargs=None):
     """
     Find current given any voltage.
 
@@ -247,22 +249,59 @@ def bishop88_i_from_v(voltage, photocurrent, saturation_current,
     method : str, default 'newton'
        Either ``'newton'`` or ``'brentq'``. ''method'' must be ``'newton'``
        if ``breakdown_factor`` is not 0.
+    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.
+        ``'full_output': True`` is allowed, and ``optimizer_output`` would be
+        returned. See examples section.
 
     Returns
     -------
     current : numeric
         current (I) at the specified voltage (V). [A]
+    optimizer_output : tuple, optional, if specified in ``method_kwargs``
+        see root finder documentation for selected method.
+        Found root is diode voltage in [1]_.
+
+    Examples
+    --------
+    Using the following arguments that may come from any
+    `calcparams_.*` function in :py:mod:`pvlib.pvsystem`:
+
+    >>> args = {'photocurrent': 1., 'saturation_current': 9e-10, 'nNsVth': 4.,
+    ...         'resistance_series': 4., 'resistance_shunt': 5000.0}
+
+    Use default values:
+
+    >>> i = bishop88_i_from_v(0.0, **args)
+
+    Specify tolerances and maximum number of iterations:
+
+    >>> i = bishop88_i_from_v(0.0, **args, method='newton',
+    ...     method_kwargs={'tol': 1e-3, 'rtol': 1e-3, 'maxiter': 20})
+
+    Retrieve full output from the root finder:
+
+    >>> i, method_output = bishop88_i_from_v(0.0, **args, method='newton',
+    ...     method_kwargs={'full_output': True})
     """
     # collect args
     args = (photocurrent, saturation_current, resistance_series,
             resistance_shunt, nNsVth, d2mutau, NsVbi,
             breakdown_factor, breakdown_voltage, breakdown_exp)
+    method = method.lower()
+
+    # method_kwargs create dict if not provided
+    # this pattern avoids bugs with Mutable Default Parameters
+    if not method_kwargs:
+        method_kwargs = {}
 
     def fv(x, v, *a):
         # calculate voltage residual given diode voltage "x"
         return bishop88(x, *a)[1] - v
 
-    if method.lower() == 'brentq':
+    if method == 'brentq':
         # first bound the search using voc
         voc_est = estimate_voc(photocurrent, saturation_current, nNsVth)
 
@@ -274,27 +313,37 @@ def vd_from_brent(voc, v, iph, isat, rs, rsh, gamma, d2mutau, NsVbi,
             return brentq(fv, 0.0, voc,
                           args=(v, iph, isat, rs, rsh, gamma, d2mutau, NsVbi,
                                 breakdown_factor, breakdown_voltage,
-                                breakdown_exp))
+                                breakdown_exp),
+                          **method_kwargs)
 
         vd_from_brent_vectorized = np.vectorize(vd_from_brent)
         vd = vd_from_brent_vectorized(voc_est, voltage, *args)
-    elif method.lower() == 'newton':
+    elif method == 'newton':
         # make sure all args are numpy arrays if max size > 1
         # if voltage is an array, then make a copy to use for initial guess, v0
-        args, v0 = _prepare_newton_inputs((voltage,), args, voltage)
+        args, v0, method_kwargs = \
+            _prepare_newton_inputs((voltage,), args, voltage, method_kwargs)
         vd = newton(func=lambda x, *a: fv(x, voltage, *a), x0=v0,
                     fprime=lambda x, *a: bishop88(x, *a, gradients=True)[4],
-                    args=args)
+                    args=args,
+                    **method_kwargs)
     else:
         raise NotImplementedError("Method '%s' isn't implemented" % method)
-    return bishop88(vd, *args)[0]
+
+    # When 'full_output' parameter is specified, returned 'vd' is a tuple with
+    # many elements, where the root is the first one. So we use it to output
+    # the bishop88 result and return tuple(scalar, tuple with method results)
+    if method_kwargs.get('full_output') is True:
+        return (bishop88(vd[0], *args)[0], vd)
+    else:
+        return bishop88(vd, *args)[0]
 
 
 def bishop88_v_from_i(current, photocurrent, saturation_current,
                       resistance_series, resistance_shunt, nNsVth,
                       d2mutau=0, NsVbi=np.Inf, breakdown_factor=0.,
                       breakdown_voltage=-5.5, breakdown_exp=3.28,
-                      method='newton'):
+                      method='newton', method_kwargs=None):
     """
     Find voltage given any current.
 
@@ -335,24 +384,62 @@ def bishop88_v_from_i(current, photocurrent, saturation_current,
     method : str, default 'newton'
        Either ``'newton'`` or ``'brentq'``. ''method'' must be ``'newton'``
        if ``breakdown_factor`` is not 0.
+    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.
+        ``'full_output': True`` is allowed, and ``optimizer_output`` would be
+        returned. See examples section.
 
     Returns
     -------
     voltage : numeric
         voltage (V) at the specified current (I) in volts [V]
+    optimizer_output : tuple, optional, if specified in ``method_kwargs``
+        see root finder documentation for selected method.
+        Found root is diode voltage in [1]_.
+
+    Examples
+    --------
+    Using the following arguments that may come from any
+    `calcparams_.*` function in :py:mod:`pvlib.pvsystem`:
+
+    >>> args = {'photocurrent': 1., 'saturation_current': 9e-10, 'nNsVth': 4.,
+    ...         'resistance_series': 4., 'resistance_shunt': 5000.0}
+
+    Use default values:
+
+    >>> v = bishop88_v_from_i(0.0, **args)
+
+    Specify tolerances and maximum number of iterations:
+
+    >>> v = bishop88_v_from_i(0.0, **args, method='newton',
+    ...     method_kwargs={'tol': 1e-3, 'rtol': 1e-3, 'maxiter': 20})
+
+    Retrieve full output from the root finder:
+
+    >>> v, method_output = bishop88_v_from_i(0.0, **args, method='newton',
+    ...     method_kwargs={'full_output': True})
     """
     # collect args
     args = (photocurrent, saturation_current, resistance_series,
             resistance_shunt, nNsVth, d2mutau, NsVbi, breakdown_factor,
             breakdown_voltage, breakdown_exp)
+    method = method.lower()
+
+    # method_kwargs create dict if not provided
+    # this pattern avoids bugs with Mutable Default Parameters
+    if not method_kwargs:
+        method_kwargs = {}
+
     # first bound the search using voc
     voc_est = estimate_voc(photocurrent, saturation_current, nNsVth)
 
     def fi(x, i, *a):
         # calculate current residual given diode voltage "x"
         return bishop88(x, *a)[0] - i
 
-    if method.lower() == 'brentq':
+    if method == 'brentq':
         # brentq only works with scalar inputs, so we need a set up function
         # and np.vectorize to repeatedly call the optimizer with the right
         # arguments for possible array input
@@ -361,26 +448,36 @@ def vd_from_brent(voc, i, iph, isat, rs, rsh, gamma, d2mutau, NsVbi,
             return brentq(fi, 0.0, voc,
                           args=(i, iph, isat, rs, rsh, gamma, d2mutau, NsVbi,
                                 breakdown_factor, breakdown_voltage,
-                                breakdown_exp))
+                                breakdown_exp),
+                          **method_kwargs)
 
         vd_from_brent_vectorized = np.vectorize(vd_from_brent)
         vd = vd_from_brent_vectorized(voc_est, current, *args)
-    elif method.lower() == 'newton':
+    elif method == 'newton':
         # make sure all args are numpy arrays if max size > 1
         # if voc_est is an array, then make a copy to use for initial guess, v0
-        args, v0 = _prepare_newton_inputs((current,), args, voc_est)
+        args, v0, method_kwargs = \
+            _prepare_newton_inputs((current,), args, voc_est, method_kwargs)
         vd = newton(func=lambda x, *a: fi(x, current, *a), x0=v0,
                     fprime=lambda x, *a: bishop88(x, *a, gradients=True)[3],
-                    args=args)
+                    args=args,
+                    **method_kwargs)
     else:
         raise NotImplementedError("Method '%s' isn't implemented" % method)
-    return bishop88(vd, *args)[1]
+
+    # When 'full_output' parameter is specified, returned 'vd' is a tuple with
+    # many elements, where the root is the first one. So we use it to output
+    # the bishop88 result and return tuple(scalar, tuple with method results)
+    if method_kwargs.get('full_output') is True:
+        return (bishop88(vd[0], *args)[1], vd)
+    else:
+        return bishop88(vd, *args)[1]
 
 
 def bishop88_mpp(photocurrent, saturation_current, resistance_series,
                  resistance_shunt, nNsVth, d2mutau=0, NsVbi=np.Inf,
                  breakdown_factor=0., breakdown_voltage=-5.5,
-                 breakdown_exp=3.28, method='newton'):
+                 breakdown_exp=3.28, method='newton', method_kwargs=None):
     """
     Find max power point.
 
@@ -419,43 +516,91 @@ def bishop88_mpp(photocurrent, saturation_current, resistance_series,
     method : str, default 'newton'
        Either ``'newton'`` or ``'brentq'``. ''method'' must be ``'newton'``
        if ``breakdown_factor`` is not 0.
+    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.
+        ``'full_output': True`` is allowed, and ``optimizer_output`` would be
+        returned. See examples section.
 
     Returns
     -------
     tuple
         max power current ``i_mp`` [A], max power voltage ``v_mp`` [V], and
         max power ``p_mp`` [W]
+    optimizer_output : tuple, optional, if specified in ``method_kwargs``
+        see root finder documentation for selected method.
+        Found root is diode voltage in [1]_.
+
+    Examples
+    --------
+    Using the following arguments that may come from any
+    `calcparams_.*` function in :py:mod:`pvlib.pvsystem`:
+
+    >>> args = {'photocurrent': 1., 'saturation_current': 9e-10, 'nNsVth': 4.,
+    ...         'resistance_series': 4., 'resistance_shunt': 5000.0}
+
+    Use default values:
+
+    >>> i_mp, v_mp, p_mp = bishop88_mpp(**args)
+
+    Specify tolerances and maximum number of iterations:
+
+    >>> i_mp, v_mp, p_mp = bishop88_mpp(**args, method='newton',
+    ...     method_kwargs={'tol': 1e-3, 'rtol': 1e-3, 'maxiter': 20})
+
+    Retrieve full output from the root finder:
+
+    >>> (i_mp, v_mp, p_mp), method_output = bishop88_mpp(**args,
+    ...     method='newton', method_kwargs={'full_output': True})
     """
     # collect args
     args = (photocurrent, saturation_current, resistance_series,
             resistance_shunt, nNsVth, d2mutau, NsVbi, breakdown_factor,
             breakdown_voltage, breakdown_exp)
+    method = method.lower()
+
+    # method_kwargs create dict if not provided
+    # this pattern avoids bugs with Mutable Default Parameters
+    if not method_kwargs:
+        method_kwargs = {}
+
     # first bound the search using voc
     voc_est = estimate_voc(photocurrent, saturation_current, nNsVth)
 
     def fmpp(x, *a):
         return bishop88(x, *a, gradients=True)[6]
 
-    if method.lower() == 'brentq':
+    if method == 'brentq':
         # break out arguments for numpy.vectorize to handle broadcasting
         vec_fun = np.vectorize(
             lambda voc, iph, isat, rs, rsh, gamma, d2mutau, NsVbi, vbr_a, vbr,
             vbr_exp: brentq(fmpp, 0.0, voc,
                             args=(iph, isat, rs, rsh, gamma, d2mutau, NsVbi,
-                                  vbr_a, vbr, vbr_exp))
+                                  vbr_a, vbr, vbr_exp),
+                            **method_kwargs)
         )
         vd = vec_fun(voc_est, *args)
-    elif method.lower() == 'newton':
+    elif method == 'newton':
         # make sure all args are numpy arrays if max size > 1
         # if voc_est is an array, then make a copy to use for initial guess, v0
-        args, v0 = _prepare_newton_inputs((), args, voc_est)
+        args, v0, method_kwargs = \
+            _prepare_newton_inputs((), args, voc_est, method_kwargs)
         vd = newton(
             func=fmpp, x0=v0,
-            fprime=lambda x, *a: bishop88(x, *a, gradients=True)[7], args=args
-        )
+            fprime=lambda x, *a: bishop88(x, *a, gradients=True)[7], args=args,
+            **method_kwargs)
     else:
         raise NotImplementedError("Method '%s' isn't implemented" % method)
-    return bishop88(vd, *args)
+
+    # When 'full_output' parameter is specified, returned 'vd' is a tuple with
+    # many elements, where the root is the first one. So we use it to output
+    # the bishop88 result and return
+    # tuple(tuple with bishop88 solution, tuple with method results)
+    if method_kwargs.get('full_output') is True:
+        return (bishop88(vd[0], *args), vd)
+    else:
+        return bishop88(vd, *args)
 
 
 def _get_size_and_shape(args):
@@ -482,7 +627,7 @@ def _get_size_and_shape(args):
     return size, shape
 
 
-def _prepare_newton_inputs(i_or_v_tup, args, v0):
+def _prepare_newton_inputs(i_or_v_tup, args, v0, method_kwargs):
     # broadcast arguments for newton method
     # the first argument should be a tuple, eg: (i,), (v,) or ()
     size, shape = _get_size_and_shape(i_or_v_tup + args)
@@ -492,7 +637,12 @@ def _prepare_newton_inputs(i_or_v_tup, args, v0):
     # copy v0 to a new array and broadcast it to the shape of max size
     if shape is not None:
         v0 = np.broadcast_to(v0, shape).copy()
-    return args, v0
+
+    # set abs tolerance and maxiter from method_kwargs if not provided
+    # apply defaults, but giving priority to user-specified values
+    method_kwargs = {**NEWTON_DEFAULT_PARAMS, **method_kwargs}
+
+    return args, v0, method_kwargs
 
 
 def _lambertw_v_from_i(current, photocurrent, saturation_current,