diff --git a/docs/sphinx/source/api.rst b/docs/sphinx/source/api.rst index c37bedf7df..26786d4f88 100644 --- a/docs/sphinx/source/api.rst +++ b/docs/sphinx/source/api.rst @@ -303,6 +303,7 @@ Functions for fitting PV models ivtools.fit_sde_sandia ivtools.fit_sdm_cec_sam + ivtools.fit_sdm_desoto Other ----- diff --git a/docs/sphinx/source/whatsnew/v0.7.0.rst b/docs/sphinx/source/whatsnew/v0.7.0.rst index 2b65b17a9e..fdc419ef6e 100644 --- a/docs/sphinx/source/whatsnew/v0.7.0.rst +++ b/docs/sphinx/source/whatsnew/v0.7.0.rst @@ -1,4 +1,4 @@ -.. _whatsnew_0700: +.. _whatsnew_0700: v0.7.0 (MONTH DAY, YEAR) ------------------------ @@ -111,6 +111,8 @@ Enhancements the single diode equation to an IV curve. * Add :py:func:`~pvlib.ivtools.fit_sdm_cec_sam`, a wrapper for the CEC single diode model fitting function '6parsolve' from NREL's System Advisor Model. +* Add :py:func:`~pvlib.ivtools.fit_sdm_desoto`, a method to fit the De Soto single + diode model to the typical specifications given in manufacturers datasheets. * Add `timeout` to :py:func:`pvlib.iotools.get_psm3`. Bug fixes @@ -162,4 +164,5 @@ Contributors * Anton Driesse (:ghuser:`adriesse`) * Alexander Morgan (:ghuser:`alexandermorgan`) * Miguel Sánchez de León Peque (:ghuser:`Peque`) +* Tanguy Lunel (:ghuser:`tylunel`) * Veronica Guo (:ghuser:`veronicaguo`) diff --git a/pvlib/ivtools.py b/pvlib/ivtools.py index 353a0b8a9a..4806b2ed16 100644 --- a/pvlib/ivtools.py +++ b/pvlib/ivtools.py @@ -262,6 +262,147 @@ def fit_sde_sandia(voltage, current, v_oc=None, i_sc=None, v_mp_i_mp=None, v_oc) +def fit_sdm_desoto(v_mp, i_mp, v_oc, i_sc, alpha_sc, beta_voc, + cells_in_series, EgRef=1.121, dEgdT=-0.0002677, + temp_ref=25, irrad_ref=1000, root_kwargs={}): + """ + Calculates the parameters for the De Soto single diode model using the + procedure described in [1]. This procedure has the advantage of + using common specifications given by manufacturers in the + datasheets of PV modules. + + The solution is found using the scipy.optimize.root() function, + with the corresponding default solver method 'hybr'. + No restriction is put on the fit variables, i.e. series + or shunt resistance could go negative. Nevertheless, if it happens, + check carefully the inputs and their units; alpha_sc and beta_voc are + often given in %/K in manufacturers datasheets and should be given + in A/K and V/K here. + + The parameters returned by this function can be used by + pvsystem.calcparams_desoto to calculate the values at different + irradiance and cell temperature. + + Parameters + ---------- + v_mp: float + Module voltage at the maximum-power point at reference conditions [V]. + i_mp: float + Module current at the maximum-power point at reference conditions [A]. + v_oc: float + Open-circuit voltage at reference conditions [V]. + i_sc: float + Short-circuit current at reference conditions [A]. + alpha_sc: float + The short-circuit current (i_sc) temperature coefficient of the + module [A/K]. + beta_voc: float + The open-circuit voltage (v_oc) temperature coefficient of the + module [V/K]. + cells_in_series: integer + Number of cell in the module. + EgRef: float, default 1.121 eV - value for silicon + Energy of bandgap of semi-conductor used [eV] + dEgdT: float, default -0.0002677 - value for silicon + Variation of bandgap according to temperature [eV/K] + temp_ref: float, default 25 + Reference temperature condition [C] + irrad_ref: float, default 1000 + Reference irradiance condition [W/m2] + root_kwargs: dictionary, default None + Dictionary of arguments to pass onto scipy.optimize.root() + + Returns + ------- + Tuple of the following elements: + + * Dictionary with the following elements: + I_L_ref: float + Light-generated current at reference conditions [A] + I_o_ref: float + Diode saturation current at reference conditions [A] + R_s: float + Series resistance [ohms] + R_sh_ref: float + Shunt resistance at reference conditions [ohms]. + a_ref: float + Modified ideality factor at reference conditions. + The product of the usual diode ideality factor (n, unitless), + number of cells in series (Ns), and cell thermal voltage at + specified effective irradiance and cell temperature. + alpha_sc: float + The short-circuit current (i_sc) temperature coefficient of the + module [A/K]. + EgRef: float + Energy of bandgap of semi-conductor used [eV] + dEgdT: float + Variation of bandgap according to temperature [eV/K] + irrad_ref: float + Reference irradiance condition [W/m2] + temp_ref: float + Reference temperature condition [C] + * scipy.optimize.OptimizeResult + Optimization result of scipy.optimize.root(). + See scipy.optimize.OptimizeResult for more details. + + References + ---------- + [1] W. De Soto et al., "Improvement and validation of a model for + photovoltaic array performance", Solar Energy, vol 80, pp. 78-88, + 2006. + + [2] John A Duffie, William A Beckman, "Solar Engineering of Thermal + Processes", Wiley, 2013 + """ + + try: + from scipy.optimize import root + from scipy import constants + except ImportError: + raise ImportError("The fit_sdm_desoto function requires scipy.") + + # Constants + k = constants.value('Boltzmann constant in eV/K') + Tref = temp_ref + 273.15 # [K] + + # initial guesses of variables for computing convergence: + # Values are taken from [2], p753 + Rsh_0 = 100.0 + a_0 = 1.5*k*Tref*cells_in_series + IL_0 = i_sc + Io_0 = i_sc * np.exp(-v_oc/a_0) + Rs_0 = (a_0*np.log1p((IL_0-i_mp)/Io_0) - v_mp)/i_mp + # params_i : initial values vector + params_i = np.array([IL_0, Io_0, a_0, Rsh_0, Rs_0]) + + # specs of module + specs = (i_sc, v_oc, i_mp, v_mp, beta_voc, alpha_sc, EgRef, dEgdT, + Tref, k) + + # computing with system of equations described in [1] + optimize_result = root(_system_of_equations_desoto, x0=params_i, + args=(specs,), **root_kwargs) + + if optimize_result.success: + sdm_params = optimize_result.x + else: + raise RuntimeError( + 'Parameter estimation failed:\n' + optimize_result.message) + + # results + return ({'I_L_ref': sdm_params[0], + 'I_o_ref': sdm_params[1], + 'a_ref': sdm_params[2], + 'R_sh_ref': sdm_params[3], + 'R_s': sdm_params[4], + 'alpha_sc': alpha_sc, + 'EgRef': EgRef, + 'dEgdT': dEgdT, + 'irrad_ref': irrad_ref, + 'temp_ref': temp_ref}, + optimize_result) + + def _find_mp(voltage, current): """ Finds voltage and current at maximum power point. @@ -348,3 +489,69 @@ def _calculate_sde_parameters(beta0, beta1, beta3, beta4, v_mp, i_mp, v_oc): else: # I0_voc > 0 I0 = I0_voc return (IL, I0, Rsh, Rs, nNsVth) + + +def _system_of_equations_desoto(params, specs): + """Evaluates the systems of equations used to solve for the single + diode equation parameters. Function designed to be used by + scipy.optimize.root() in fit_sdm_desoto(). + + Parameters + ---------- + params: ndarray + Array with parameters of the De Soto single diode model. Must be + given in the following order: IL, Io, a, Rsh, Rs + specs: tuple + Specifications of pv module given by manufacturer. Must be given + in the following order: Isc, Voc, Imp, Vmp, beta_oc, alpha_sc + + Returns + ------- + system of equations to solve with scipy.optimize.root(). + + + References + ---------- + [1] W. De Soto et al., "Improvement and validation of a model for + photovoltaic array performance", Solar Energy, vol 80, pp. 78-88, + 2006. + + [2] John A Duffie, William A Beckman, "Solar Engineering of Thermal + Processes", Wiley, 2013 + """ + + # six input known variables + Isc, Voc, Imp, Vmp, beta_oc, alpha_sc, EgRef, dEgdT, Tref, k = specs + + # five parameters vector to find + IL, Io, a, Rsh, Rs = params + + # five equation vector + y = [0, 0, 0, 0, 0] + + # 1st equation - short-circuit - eq(3) in [1] + y[0] = Isc - IL + Io * np.expm1(Isc * Rs / a) + Isc * Rs / Rsh + + # 2nd equation - open-circuit Tref - eq(4) in [1] + y[1] = -IL + Io * np.expm1(Voc / a) + Voc / Rsh + + # 3rd equation - Imp & Vmp - eq(5) in [1] + y[2] = Imp - IL + Io * np.expm1((Vmp + Imp * Rs) / a) \ + + (Vmp + Imp * Rs) / Rsh + + # 4th equation - Pmp derivated=0 - eq23.2.6 in [2] + # caution: eq(6) in [1] has a sign error + y[3] = Imp \ + - Vmp * ((Io / a) * np.exp((Vmp + Imp * Rs) / a) + 1.0 / Rsh) \ + / (1.0 + (Io * Rs / a) * np.exp((Vmp + Imp * Rs) / a) + Rs / Rsh) + + # 5th equation - open-circuit T2 - eq (4) at temperature T2 in [1] + T2 = Tref + 2 + Voc2 = (T2 - Tref) * beta_oc + Voc # eq (7) in [1] + a2 = a * T2 / Tref # eq (8) in [1] + IL2 = IL + alpha_sc * (T2 - Tref) # eq (11) in [1] + Eg2 = EgRef * (1 + dEgdT * (T2 - Tref)) # eq (10) in [1] + Io2 = Io * (T2 / Tref)**3 * np.exp(1 / k * (EgRef/Tref - Eg2/T2)) # eq (9) + y[4] = -IL2 + Io2 * np.expm1(Voc2 / a2) + Voc2 / Rsh # eq (4) at T2 + + return y diff --git a/pvlib/test/test_ivtools.py b/pvlib/test/test_ivtools.py index d4aca050c4..6f9fdd6fec 100644 --- a/pvlib/test/test_ivtools.py +++ b/pvlib/test/test_ivtools.py @@ -102,6 +102,35 @@ def test_fit_sdm_cec_sam(get_cec_params_cansol_cs5p_220p): cells_in_series=1, temp_ref=25) +@requires_scipy +def test_fit_sdm_desoto(): + result, _ = ivtools.fit_sdm_desoto(v_mp=31.0, i_mp=8.71, v_oc=38.3, + i_sc=9.43, alpha_sc=0.005658, + beta_voc=-0.13788, + cells_in_series=60) + result_expected = {'I_L_ref': 9.45232, + 'I_o_ref': 3.22460e-10, + 'a_ref': 1.59128, + 'R_sh_ref': 125.798, + 'R_s': 0.297814, + 'alpha_sc': 0.005658, + 'EgRef': 1.121, + 'dEgdT': -0.0002677, + 'irrad_ref': 1000, + 'temp_ref': 25} + assert np.allclose(pd.Series(result), pd.Series(result_expected), + rtol=1e-4) + + +@requires_scipy +def test_fit_sdm_desoto_failure(): + with pytest.raises(RuntimeError) as exc: + ivtools.fit_sdm_desoto(v_mp=31.0, i_mp=8.71, v_oc=38.3, i_sc=9.43, + alpha_sc=0.005658, beta_voc=-0.13788, + cells_in_series=10) + assert ('Parameter estimation failed') in str(exc.value) + + @pytest.fixture def get_bad_iv_curves(): # v1, i1 produces a bad value for I0_voc