diff --git a/docs/sphinx/source/whatsnew/v0.10.0.rst b/docs/sphinx/source/whatsnew/v0.10.0.rst index 87b62c8a05..ea8ba7df87 100644 --- a/docs/sphinx/source/whatsnew/v0.10.0.rst +++ b/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`) - diff --git a/pvlib/singlediode.py b/pvlib/singlediode.py index 81d6ce3761..a4a760f1b9 100644 --- a/pvlib/singlediode.py +++ b/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,16 +384,54 @@ 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) @@ -352,7 +439,7 @@ 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, diff --git a/pvlib/tests/test_singlediode.py b/pvlib/tests/test_singlediode.py index 9fddae07f4..1223e73e8e 100644 --- a/pvlib/tests/test_singlediode.py +++ b/pvlib/tests/test_singlediode.py @@ -412,3 +412,135 @@ 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) + + +@pytest.fixture +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, method_kwargs', [ + ('newton', { + 'tol': 1e-8, + 'rtol': 1e-8, + 'maxiter': 30, + }), + ('brentq', { + 'xtol': 1e-8, + 'rtol': 1e-8, + 'maxiter': 30, + }) +]) +def test_bishop88_kwargs_transfer(method, method_kwargs, mocker, + bishop88_arguments): + """test method_kwargs modifying optimizer does not break anything""" + # patch method namespace at singlediode module namespace + optimizer_mock = mocker.patch('pvlib.singlediode.' + method) + + # check kwargs passed to bishop_.* are a subset of the call args + # since they are called with more keyword arguments + + bishop88_i_from_v(0, **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + _, kwargs = optimizer_mock.call_args + assert method_kwargs.items() <= kwargs.items() + + bishop88_v_from_i(0, **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + _, kwargs = optimizer_mock.call_args + assert method_kwargs.items() <= kwargs.items() + + bishop88_mpp(**bishop88_arguments, method=method, + method_kwargs=method_kwargs) + _, kwargs = optimizer_mock.call_args + assert method_kwargs.items() <= kwargs.items() + + +@pytest.mark.parametrize('method, 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, method_kwargs, bishop88_arguments): + """test invalid method_kwargs passed onto the optimizer fail""" + + pytest.raises(TypeError, bishop88_i_from_v, + 0, **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + + pytest.raises(TypeError, bishop88_v_from_i, + 0, **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + + pytest.raises(TypeError, bishop88_mpp, + **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + + +@pytest.mark.parametrize('method', ['newton', 'brentq']) +def test_bishop88_full_output_kwarg(method, bishop88_arguments): + """test call to bishop88_.* with full_output=True return values are ok""" + method_kwargs = {'full_output': True} + + ret_val = bishop88_i_from_v(0, **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + assert isinstance(ret_val, tuple) # ret_val must be a tuple + assert len(ret_val) == 2 # of two elements + assert isinstance(ret_val[0], float) # first one has bishop88 result + assert isinstance(ret_val[1], tuple) # second is output from optimizer + # any root finder returns at least 2 elements with full_output=True + assert len(ret_val[1]) >= 2 + + ret_val = bishop88_v_from_i(0, **bishop88_arguments, method=method, + method_kwargs=method_kwargs) + assert isinstance(ret_val, tuple) # ret_val must be a tuple + assert len(ret_val) == 2 # of two elements + assert isinstance(ret_val[0], float) # first one has bishop88 result + assert isinstance(ret_val[1], tuple) # second is output from optimizer + # any root finder returns at least 2 elements with full_output=True + assert len(ret_val[1]) >= 2 + + ret_val = bishop88_mpp(**bishop88_arguments, method=method, + method_kwargs=method_kwargs) + assert isinstance(ret_val, tuple) # ret_val must be a tuple + assert len(ret_val) == 2 # of two elements + assert isinstance(ret_val[0], tuple) # first one has bishop88 result + assert len(ret_val[0]) == 3 # of three elements (I,V,P) + assert isinstance(ret_val[1], tuple) # second is output from optimizer + # any root finder returns at least 2 elements with full_output=True + assert len(ret_val[1]) >= 2