Allow users to interact with root finders parameters

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

@@ -21,6 +21,11 @@ Deprecations
 
 Enhancements
 ~~~~~~~~~~~~
+* Allow passing keyword arguments to :py:func:`scipy:scipy.optimize.brentq` and
+  :py:func:`scipy:scipy.optimize.newton` solvers in :py:func:`bishop88_mpp`,
+  :py:func:`bishop88_i_from_v` and :py:func:`bishop88_v_from_i`. Among others,
+  tolerance and number of iterations can be set.
+  (:issue:`1249`, :pull:`1764`)
 
 
 Bug fixes
@@ -45,3 +50,4 @@ Requirements
 Contributors
 ~~~~~~~~~~~~
 * Taos Transue (:ghuser:`reepoi`)
+* Echedey Luis (:ghuser:`echedey-ls`)

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', **kwargs):
     """
     Find current given any voltage.
 
@@ -247,6 +249,10 @@ 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.
+    kwargs : dict
+        Passed to root finder method. See
+        :py:func:`scipy:scipy.optimize.brentq` and
+        :py:func:`scipy:scipy.optimize.newton` parameters.
 
     Returns
     -------
@@ -257,12 +263,13 @@ def bishop88_i_from_v(voltage, photocurrent, saturation_current,
     args = (photocurrent, saturation_current, resistance_series,
             resistance_shunt, nNsVth, d2mutau, NsVbi,
             breakdown_factor, breakdown_voltage, breakdown_exp)
+    method = method.lower()
 
     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 +281,31 @@ 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),
+                          **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, kwargs = _prepare_newton_inputs((voltage,), args, voltage,
+                                                  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,
+                    **kwargs)
     else:
         raise NotImplementedError("Method '%s' isn't implemented" % method)
+
     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', **kwargs):
     """
     Find voltage given any current.
 
@@ -335,6 +346,10 @@ 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.
+    kwargs : dict
+        Passed to root finder method. See
+        :py:func:`scipy:scipy.optimize.brentq` and
+        :py:func:`scipy:scipy.optimize.newton` parameters.
 
     Returns
     -------
@@ -345,14 +360,15 @@ def bishop88_v_from_i(current, photocurrent, saturation_current,
     args = (photocurrent, saturation_current, resistance_series,
             resistance_shunt, nNsVth, d2mutau, NsVbi, breakdown_factor,
             breakdown_voltage, breakdown_exp)
+    method = method.lower()
     # 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,17 +377,20 @@ 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),
+                          **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, kwargs = _prepare_newton_inputs((current,), args, voc_est,
+                                                  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,
+                    **kwargs)
     else:
         raise NotImplementedError("Method '%s' isn't implemented" % method)
     return bishop88(vd, *args)[1]
@@ -380,7 +399,7 @@ def vd_from_brent(voc, i, iph, isat, rs, rsh, gamma, d2mutau, NsVbi,
 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', **kwargs):
     """
     Find max power point.
 
@@ -419,6 +438,10 @@ 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.
+    kwargs : dict
+        Passed to root finder method. See
+        :py:func:`scipy:scipy.optimize.brentq` and
+        :py:func:`scipy:scipy.optimize.newton` parameters.
 
     Returns
     -------
@@ -430,29 +453,31 @@ def bishop88_mpp(photocurrent, saturation_current, resistance_series,
     args = (photocurrent, saturation_current, resistance_series,
             resistance_shunt, nNsVth, d2mutau, NsVbi, breakdown_factor,
             breakdown_voltage, breakdown_exp)
+    method = method.lower()
     # 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),
+                            **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, kwargs = _prepare_newton_inputs((), args, voc_est, 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,
+            **kwargs)
     else:
         raise NotImplementedError("Method '%s' isn't implemented" % method)
     return bishop88(vd, *args)
@@ -482,7 +507,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, 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 +517,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 kwargs if not provided
+    kwargs['tol'] = kwargs.pop('tol', NEWTON_DEFAULT_PARAMS['tol'])
+    kwargs['maxiter'] = kwargs.pop('maxiter', NEWTON_DEFAULT_PARAMS['maxiter'])
+
+    return args, v0, kwargs
 
 
 def _lambertw_v_from_i(current, photocurrent, saturation_current,

pvlib/tests/test_singlediode.py

@@ -408,3 +408,112 @@ def test_pvsyst_breakdown(method, brk_params, recomb_params, poa, temp_cell,
 
     vsc_88 = bishop88_v_from_i(isc_88, *x, **y, method=method)
     assert np.isclose(vsc_88, 0.0, *tol)
+
+
+def bishop88_arguments():
+    pvsyst_fs_495 = get_pvsyst_fs_495()
+    # evaluate PVSyst model with thin-film recombination loss current
+    # at reference conditions
+    x = pvsystem.calcparams_pvsyst(
+        effective_irradiance=pvsyst_fs_495['irrad_ref'],
+        temp_cell=pvsyst_fs_495['temp_ref'],
+        alpha_sc=pvsyst_fs_495['alpha_sc'],
+        gamma_ref=pvsyst_fs_495['gamma_ref'],
+        mu_gamma=pvsyst_fs_495['mu_gamma'], I_L_ref=pvsyst_fs_495['I_L_ref'],
+        I_o_ref=pvsyst_fs_495['I_o_ref'], R_sh_ref=pvsyst_fs_495['R_sh_ref'],
+        R_sh_0=pvsyst_fs_495['R_sh_0'], R_sh_exp=pvsyst_fs_495['R_sh_exp'],
+        R_s=pvsyst_fs_495['R_s'],
+        cells_in_series=pvsyst_fs_495['cells_in_series'],
+        EgRef=pvsyst_fs_495['EgRef']
+    )
+    y = dict(d2mutau=pvsyst_fs_495['d2mutau'],
+             NsVbi=VOLTAGE_BUILTIN*pvsyst_fs_495['cells_in_series'])
+    # Convert (*x, **y) in a bishop88_.* call to dict of arguments
+    args_dict = {
+        'photocurrent': x[0],
+        'saturation_current': x[1],
+        'resistance_series': x[2],
+        'resistance_shunt': x[3],
+        'nNsVth': x[4],
+    }
+    args_dict.update(y)
+    return args_dict
+
+
+@pytest.mark.parametrize('method, kwargs', [
+    ('newton', {
+        'tol': 1e-3,
+        'rtol': 1e-3,
+        'maxiter': 30,
+    }),
+    ('brentq', {
+        'xtol': 1e-3,
+        'rtol': 1e-3,
+        'maxiter': 30,
+    })
+])
+def test_bishop88_kwargs_pass(method, kwargs):
+    """test kwargs modifying optimizer does not break anything"""
+    # build test tolerances from the kwargs passed to the method
+    tol = {
+        'atol': np.nanmax([kwargs.get('tol', np.nan),
+                           kwargs.get('xtol', np.nan)]),
+        'rtol': kwargs.get('rtol')
+    }
+    expected = {  # from reference conditions
+        'pmp': (get_pvsyst_fs_495()['I_mp_ref'] *
+                get_pvsyst_fs_495()['V_mp_ref']),
+        'isc': get_pvsyst_fs_495()['I_sc_ref'],
+        'voc': get_pvsyst_fs_495()['V_oc_ref']
+    }
+
+    # get arguments common to all bishop88_.* functions
+    bishop88_args = bishop88_arguments()
+
+    mpp_88 = bishop88_mpp(**bishop88_args, method=method, **kwargs)
+    assert np.isclose(mpp_88[2], expected['pmp'], **tol)
+
+    isc_88 = bishop88_i_from_v(0, **bishop88_args, method=method,
+                               **kwargs)
+    assert np.isclose(isc_88, expected['isc'], **tol)
+
+    voc_88 = bishop88_v_from_i(0, **bishop88_args, method=method,
+                               **kwargs)
+    assert np.isclose(voc_88, expected['voc'], **tol)
+
+    ioc_88 = bishop88_i_from_v(voc_88, **bishop88_args, method=method,
+                               **kwargs)
+    assert np.isclose(ioc_88, 0.0, **tol)
+
+    vsc_88 = bishop88_v_from_i(isc_88, **bishop88_args, method=method,
+                               **kwargs)
+    assert np.isclose(vsc_88, 0.0, **tol)
+
+
+@pytest.mark.parametrize('method, kwargs', [
+    ('newton', {
+        'tol': 1e-4,
+        'rtol': 1e-4,
+        'maxiter': 20,
+        '_inexistent_param': "0.01"
+    }),
+    ('brentq', {
+        'xtol': 1e-4,
+        'rtol': 1e-4,
+        'maxiter': 20,
+        '_inexistent_param': "0.01"
+    })
+])
+def test_bishop88_kwargs_fails(method, kwargs):
+    """test invalid kwargs passed onto the optimizer fail"""
+    # get arguments common to all bishop88_.* functions
+    bishop88_args = bishop88_arguments()
+
+    pytest.raises(TypeError, bishop88_mpp,
+                  **bishop88_args, method=method, **kwargs)
+
+    pytest.raises(TypeError, bishop88_i_from_v,
+                  0, **bishop88_args, method=method, **kwargs)
+
+    pytest.raises(TypeError, bishop88_v_from_i,
+                  0, **bishop88_args, method=method, **kwargs)