diff --git a/docs/sphinx/source/reference/pv_modeling/sdm.rst b/docs/sphinx/source/reference/pv_modeling/sdm.rst index bfd5103ebe..2aff477883 100644 --- a/docs/sphinx/source/reference/pv_modeling/sdm.rst +++ b/docs/sphinx/source/reference/pv_modeling/sdm.rst @@ -17,6 +17,7 @@ Functions relevant for single diode models. pvsystem.v_from_i pvsystem.max_power_point ivtools.sdm.pvsyst_temperature_coeff + singlediode.batzelis_keypoints Low-level functions for solving the single diode equation. @@ -37,3 +38,4 @@ Functions for fitting diode models ivtools.sde.fit_sandia_simple ivtools.sdm.fit_cec_sam ivtools.sdm.fit_desoto + ivtools.sdm.fit_desoto_batzelis diff --git a/docs/sphinx/source/reference/pv_modeling/system_models.rst b/docs/sphinx/source/reference/pv_modeling/system_models.rst index fb637ee8ed..3b8f29d0b2 100644 --- a/docs/sphinx/source/reference/pv_modeling/system_models.rst +++ b/docs/sphinx/source/reference/pv_modeling/system_models.rst @@ -55,3 +55,11 @@ PVGIS model :toctree: ../generated/ pvarray.huld + +Other +^^^^^ + +.. autosummary:: + :toctree: ../generated/ + + pvarray.batzelis diff --git a/docs/sphinx/source/whatsnew/v0.13.2.rst b/docs/sphinx/source/whatsnew/v0.13.2.rst index 3177f7a022..5c132831f4 100644 --- a/docs/sphinx/source/whatsnew/v0.13.2.rst +++ b/docs/sphinx/source/whatsnew/v0.13.2.rst @@ -18,6 +18,13 @@ Bug fixes Enhancements ~~~~~~~~~~~~ +* Add :py:func:`~pvlib.ivtools.sdm.fit_desoto_batzelis`, a function to estimate + parameters for the De Soto single-diode model from datasheet values. (:pull:`2563`) +* Add :py:func:`~pvlib.singlediode.batzelis_keypoints`, a function to estimate + maximum power, open circuit, and short circuit points using parameters for + the single-diode equation. (:pull:`2563`) +* Add :py:func:`~pvlib.pvarray.batzelis`, a function to estimate maximum power + open circuit, and short circuit points from datasheet values. (:pull:`2563`) * Add ``method='chandrupatla'`` (faster than ``brentq`` and slower than ``newton``, but convergence is guaranteed) as an option for :py:func:`pvlib.pvsystem.singlediode`, diff --git a/pvlib/ivtools/sdm/__init__.py b/pvlib/ivtools/sdm/__init__.py index 8535f1f5e6..2bb8b9876b 100644 --- a/pvlib/ivtools/sdm/__init__.py +++ b/pvlib/ivtools/sdm/__init__.py @@ -10,7 +10,8 @@ from pvlib.ivtools.sdm.desoto import ( # noqa: F401 fit_desoto, - fit_desoto_sandia + fit_desoto_batzelis, + fit_desoto_sandia, ) from pvlib.ivtools.sdm.pvsyst import ( # noqa: F401 diff --git a/pvlib/ivtools/sdm/desoto.py b/pvlib/ivtools/sdm/desoto.py index 6b7fa619ca..389e0f0ef0 100644 --- a/pvlib/ivtools/sdm/desoto.py +++ b/pvlib/ivtools/sdm/desoto.py @@ -2,6 +2,7 @@ from scipy import constants from scipy import optimize +from scipy.special import lambertw from pvlib.ivtools.utils import rectify_iv_curve from pvlib.ivtools.sde import _fit_sandia_cocontent @@ -399,3 +400,74 @@ def _fit_desoto_sandia_diode(ee, voc, vth, tc, specs, const): new_x = sm.add_constant(x) res = sm.RLM(y, new_x).fit() return np.array(res.params)[1] + + +def fit_desoto_batzelis(v_mp, i_mp, v_oc, i_sc, alpha_sc, beta_voc): + """ + Determine De Soto single-diode model parameters from datasheet values + using Batzelis's method. + + This method is described in Section II.C of [1]_ and fully documented + in [2]_. + + Parameters + ---------- + v_mp : float + Maximum power point voltage at STC. [V] + i_mp : float + Maximum power point current at STC. [A] + v_oc : float + Open-circuit voltage at STC. [V] + i_sc : float + Short-circuit current at STC. [A] + alpha_sc : float + Short-circuit current temperature coefficient at STC. [A/K] + beta_voc : float + Open-circuit voltage temperature coefficient at STC. [V/K] + + Returns + ------- + dict + The returned dict contains the keys: + + * ``alpha_sc`` [A/K] + * ``a_ref`` [V] + * ``I_L_ref`` [A] + * ``I_o_ref`` [A] + * ``R_sh_ref`` [Ohm] + * ``R_s`` [Ohm] + + References + ---------- + .. [1] E. I. Batzelis, "Simple PV Performance Equations Theoretically Well + Founded on the Single-Diode Model," Journal of Photovoltaics vol. 7, + no. 5, pp. 1400-1409, Sep 2017, :doi:`10.1109/JPHOTOV.2017.2711431` + .. [2] E. I. Batzelis and S. A. Papathanassiou, "A method for the + analytical extraction of the single-diode PV model parameters," + IEEE Trans. Sustain. Energy, vol. 7, no. 2, pp. 504-512, Apr 2016. + :doi:`10.1109/TSTE.2015.2503435` + """ + # convert temp coeffs from A/K and V/K to 1/K + alpha_sc = alpha_sc / i_sc + beta_voc = beta_voc / v_oc + + # Equation numbers refer to [1] + t0 = 298.15 # K + del0 = (1 - beta_voc * t0) / (50.1 - alpha_sc * t0) # Eq 9 + w0 = np.real(lambertw(np.exp(1/del0 + 1))) + + # Eqs 11-15 + a0 = del0 * v_oc + Rs0 = (a0 * (w0 - 1) - v_mp) / i_mp + Rsh0 = a0 * (w0 - 1) / (i_sc * (1 - 1/w0) - i_mp) + Iph0 = (1 + Rs0 / Rsh0) * i_sc + Isat0 = Iph0 * np.exp(-1/del0) + + return { + 'alpha_sc': alpha_sc * i_sc, # convert 1/K to A/K + 'a_ref': a0, + 'I_L_ref': Iph0, + 'I_o_ref': Isat0, + 'R_sh_ref': Rsh0, + 'R_s': Rs0, + } diff --git a/pvlib/pvarray.py b/pvlib/pvarray.py index b156baa757..b43d8d349a 100644 --- a/pvlib/pvarray.py +++ b/pvlib/pvarray.py @@ -9,8 +9,9 @@ """ import numpy as np +import pandas as pd from scipy.optimize import curve_fit -from scipy.special import exp10 +from scipy.special import exp10, lambertw def pvefficiency_adr(effective_irradiance, temp_cell, @@ -394,3 +395,131 @@ def huld(effective_irradiance, temp_mod, pdc0, k=None, cell_type=None, k[5] * tprime**2 ) return pdc + + +def batzelis(effective_irradiance, temp_cell, + v_mp, i_mp, v_oc, i_sc, alpha_sc, beta_voc): + """ + Compute maximum power point, open circuit, and short circuit + values using Batzelis's method. + + Batzelis's method (described in Section III of [1]_) is a fast method + of computing the maximum power current and voltage. The calculations + are rooted in the De Soto single-diode model, but require only typical + datasheet information. + + Parameters + ---------- + effective_irradiance : numeric, non-negative + Effective irradiance incident on the PV module. [Wm⁻²] + temp_cell : numeric + PV module operating temperature. [°C] + v_mp : float + Maximum power point voltage at STC. [V] + i_mp : float + Maximum power point current at STC. [A] + v_oc : float + Open-circuit voltage at STC. [V] + i_sc : float + Short-circuit current at STC. [A] + alpha_sc : float + Short-circuit current temperature coefficient at STC. [A/K] + beta_voc : float + Open-circuit voltage temperature coefficient at STC. [V/K] + + Returns + ------- + dict + The returned dict-like object contains the keys/columns: + + * ``p_mp`` - power at maximum power point. [W] + * ``i_mp`` - current at maximum power point. [A] + * ``v_mp`` - voltage at maximum power point. [V] + * ``i_sc`` - short circuit current. [A] + * ``v_oc`` - open circuit voltage. [V] + + Notes + ----- + This method is the combination of three sub-methods for: + + 1. estimating single-diode model parameters from datasheet information + 2. translating SDM parameters from STC to operating conditions + (taken from the De Soto model) + 3. estimating the MPP, OC, and SC points on the resulting I-V curve. + + At extremely low irradiance (e.g. 1e-10 Wm⁻²), this model can produce + negative voltages. This function clips any negative voltages to zero. + + References + ---------- + .. [1] E. I. Batzelis, "Simple PV Performance Equations Theoretically Well + Founded on the Single-Diode Model," Journal of Photovoltaics vol. 7, + no. 5, pp. 1400-1409, Sep 2017, :doi:`10.1109/JPHOTOV.2017.2711431` + + Examples + -------- + >>> params = {'i_sc': 15.98, 'v_oc': 50.26, 'i_mp': 15.27, 'v_mp': 42.57, + ... 'alpha_sc': 0.007351, 'beta_voc': -0.120624} + >>> batzelis(np.array([1000, 800]), np.array([25, 30]), **params) + {'p_mp': array([650.0439 , 512.99199048]), + 'i_mp': array([15.27 , 12.23049303]), + 'v_mp': array([42.57 , 41.94368856]), + 'i_sc': array([15.98 , 12.813404]), + 'v_oc': array([50.26 , 49.26532902])} + """ + # convert temp coeffs from A/K and V/K to 1/K + alpha_sc = alpha_sc / i_sc + beta_voc = beta_voc / v_oc + + t0 = 298.15 + delT = temp_cell - (t0 - 273.15) + lamT = (temp_cell + 273.15) / t0 + g = effective_irradiance / 1000 + # for zero/negative irradiance, use lnG=large negative number so that + # computed voltages are negative and then clipped to zero + with np.errstate(divide='ignore'): # needed for pandas for some reason + lnG = np.log(g, out=np.full_like(g, -9e9), where=(g > 0)) + lnG = np.where(np.isfinite(g), lnG, np.nan) # also preserve nans + + # Eq 9-10 + del0 = (1 - beta_voc * t0) / (50.1 - alpha_sc * t0) + w0 = np.real(lambertw(np.exp(1/del0 + 1))) + + # Eqs 27-28 + alpha_imp = alpha_sc + (beta_voc - 1/t0) / (w0 - 1) + beta_vmp = (v_oc / v_mp) * ( + beta_voc / (1 + del0) + + (del0 * (w0 - 1) - 1/(1 + del0)) / t0 + ) + + # Eq 26 + eps0 = (del0 / (1 + del0)) * (v_oc / v_mp) + eps1 = del0 * (w0 - 1) * (v_oc / v_mp) - 1 + + # Eqs 22-25 + isc = g * i_sc * (1 + alpha_sc * delT) + voc = v_oc * (1 + del0 * lamT * lnG + beta_voc * delT) + imp = g * i_mp * (1 + alpha_imp * delT) + vmp = v_mp * (1 + eps0 * lamT * lnG + eps1 * (1 - g) + beta_vmp * delT) + + # handle negative voltages from zero and extremely small irradiance + vmp = np.clip(vmp, a_min=0, a_max=None) + voc = np.clip(voc, a_min=0, a_max=None) + + out = { + 'p_mp': vmp * imp, + 'i_mp': imp, + 'v_mp': vmp, + 'i_sc': isc, + 'v_oc': voc, + } + + # if pandas in, ensure pandas out + pandas_inputs = [ + x for x in [effective_irradiance, temp_cell] + if isinstance(x, pd.Series) + ] + if pandas_inputs: + out = pd.DataFrame(out, index=pandas_inputs[0].index) + + return out diff --git a/pvlib/singlediode.py b/pvlib/singlediode.py index 43be522437..ecc6c2d6cd 100644 --- a/pvlib/singlediode.py +++ b/pvlib/singlediode.py @@ -3,6 +3,7 @@ """ import numpy as np +import pandas as pd from pvlib.tools import _golden_sect_DataFrame from scipy.optimize import brentq, newton @@ -949,3 +950,85 @@ def _pwr_optfcn(df, loc): df['resistance_shunt'], df['nNsVth']) return current * df[loc] + + +def batzelis_keypoints(photocurrent, saturation_current, resistance_series, + resistance_shunt, nNsVth): + """ + Estimate maximum power, open-circuit, and short-circuit points from + single-diode equation parameters using Batzelis's method. + + This method is described in Section II.B of [1]_. + + Parameters + ---------- + photocurrent : numeric + Light-generated current. [A] + saturation_current : numeric + Diode saturation current. [A] + resistance_series : numeric + Series resistance. [Ohm] + resistance_shunt : numeric + Shunt resistance. [Ohm] + nNsVth : numeric + 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. [V] + + Returns + ------- + dict + The returned dict-like object contains the keys/columns: + + * ``p_mp`` - power at maximum power point. [W] + * ``i_mp`` - current at maximum power point. [A] + * ``v_mp`` - voltage at maximum power point. [V] + * ``i_sc`` - short circuit current. [A] + * ``v_oc`` - open circuit voltage. [V] + + References + ---------- + .. [1] E. I. Batzelis, "Simple PV Performance Equations Theoretically Well + Founded on the Single-Diode Model," Journal of Photovoltaics vol. 7, + no. 5, pp. 1400-1409, Sep 2017, :doi:`10.1109/JPHOTOV.2017.2711431` + """ + # convenience variables + Iph = photocurrent + Is = saturation_current + Rsh = resistance_shunt + Rs = resistance_series + a = nNsVth + + # Eqs 3-4 + isc = Iph / (Rs / Rsh + 1) # manipulated to handle Rsh=np.inf correctly + with np.errstate(divide='ignore'): # zero Iph + voc = a * np.log(Iph / Is) + + # Eqs 5-8 + w = np.real(lambertw(np.e * Iph / Is)) + # vmp = (1 + Rs/Rsh) * a * (w - 1) - Rs * Iph * (1 - 1/w) # not needed + with np.errstate(divide='ignore', invalid='ignore'): # zero Iph -> zero w + imp = Iph * (1 - 1/w) - a * (w - 1) / Rsh + + vmp = a * (w - 1) - Rs * imp + + vmp = np.where(Iph > 0, vmp, 0) + voc = np.where(Iph > 0, voc, 0) + imp = np.where(Iph > 0, imp, 0) + isc = np.where(Iph > 0, isc, 0) + + out = { + 'p_mp': imp * vmp, + 'i_mp': imp, + 'v_mp': vmp, + 'i_sc': isc, + 'v_oc': voc, + } + + # if pandas in, ensure pandas out + pandas_inputs = [ + x for x in [Iph, Is, Rsh, Rs, a] if isinstance(x, pd.Series)] + if pandas_inputs: + out = pd.DataFrame(out, index=pandas_inputs[0].index) + + return out diff --git a/tests/ivtools/sdm/test_desoto.py b/tests/ivtools/sdm/test_desoto.py index b861b26819..710bf15d51 100644 --- a/tests/ivtools/sdm/test_desoto.py +++ b/tests/ivtools/sdm/test_desoto.py @@ -92,3 +92,27 @@ def test_fit_desoto_sandia(cec_params_cansol_cs5p_220p): assert_allclose(result['dEgdT'], -0.0002677) assert_allclose(result['EgRef'], 1.3112547292120638) assert_allclose(result['cells_in_series'], specs['cells_in_series']) + + +def test_fit_desoto_batzelis(): + params = {'i_sc': 15.98, 'v_oc': 50.26, 'i_mp': 15.27, 'v_mp': 42.57, + 'alpha_sc': 0.007351, 'beta_voc': -0.120624} + expected = { # calculated with the function itself + 'alpha_sc': 0.007351, + 'a_ref': 1.7257632483733483, + 'I_L_ref': 15.985408866796396, + 'I_o_ref': 3.594308384705643e-12, + 'R_sh_ref': 389.4379947026243, + 'R_s': 0.13181590981241956, + } + out = sdm.fit_desoto_batzelis(**params) + for k in expected: + assert out[k] == pytest.approx(expected[k]) + + # ensure the STC values are reproduced + iv = pvsystem.singlediode(out['I_L_ref'], out['I_o_ref'], out['R_s'], + out['R_sh_ref'], out['a_ref']) + assert iv['i_sc'] == pytest.approx(params['i_sc']) + assert iv['i_mp'] == pytest.approx(params['i_mp'], rel=3e-3) + assert iv['v_oc'] == pytest.approx(params['v_oc'], rel=3e-4) + assert iv['v_mp'] == pytest.approx(params['v_mp'], rel=4e-3) diff --git a/tests/test_pvarray.py b/tests/test_pvarray.py index a1fc1c4ac3..e614f885db 100644 --- a/tests/test_pvarray.py +++ b/tests/test_pvarray.py @@ -1,4 +1,5 @@ import numpy as np +from numpy import nan import pandas as pd from numpy.testing import assert_allclose from .conftest import assert_series_equal @@ -114,3 +115,41 @@ def test_huld_errors(): pvarray.huld( eff_irr, temp_mod, pdc0, cell_type='csi', k_version='2021' ) + + +def test_batzelis(): + params = {'i_sc': 15.98, 'v_oc': 50.26, 'i_mp': 15.27, 'v_mp': 42.57, + 'alpha_sc': 0.007351, 'beta_voc': -0.120624} + g = np.array([1000, 500, 1200, 500, 1200, 0, nan, 1000]) + t = np.array([25, 20, 20, 50, 50, 25, 0, nan]) + expected = { # these values were computed using pvarray.batzelis itself + 'p_mp': [650.044, 328.599, 789.136, 300.079, 723.401, 0, nan, nan], + 'i_mp': [ 15.270, 7.626, 18.302, 7.680, 18.433, 0, nan, nan], + 'v_mp': [ 42.570, 43.090, 43.117, 39.071, 39.246, 0, nan, nan], + 'i_sc': [ 15.980, 7.972, 19.132, 8.082, 19.397, 0, nan, nan], + 'v_oc': [ 50.260, 49.687, 51.172, 45.948, 47.585, 0, nan, nan], + } + + # numpy array + actual = pvarray.batzelis(g, t, **params) + for key, exp in expected.items(): + np.testing.assert_allclose(actual[key], exp, atol=1e-3) + + # pandas series + actual = pvarray.batzelis(pd.Series(g), pd.Series(t), **params) + assert isinstance(actual, pd.DataFrame) + for key, exp in expected.items(): + np.testing.assert_allclose(actual[key], pd.Series(exp), atol=1e-3) + + # scalar + actual = pvarray.batzelis(g[1], t[1], **params) + for key, exp in expected.items(): + assert pytest.approx(exp[1], abs=1e-3) == actual[key] + + +def test_batzelis_negative_voltage(): + params = {'i_sc': 15.98, 'v_oc': 50.26, 'i_mp': 15.27, 'v_mp': 42.57, + 'alpha_sc': 0.007351, 'beta_voc': -0.120624} + actual = pvarray.batzelis(1e-10, 25, **params) + assert actual['v_mp'] == 0 + assert actual['v_oc'] == 0 diff --git a/tests/test_singlediode.py b/tests/test_singlediode.py index 8f3b4012c3..435aa91d84 100644 --- a/tests/test_singlediode.py +++ b/tests/test_singlediode.py @@ -7,7 +7,8 @@ import scipy from pvlib import pvsystem from pvlib.singlediode import (bishop88_mpp, estimate_voc, VOLTAGE_BUILTIN, - bishop88, bishop88_i_from_v, bishop88_v_from_i) + bishop88, bishop88_i_from_v, bishop88_v_from_i, + batzelis_keypoints) import pytest from numpy.testing import assert_array_equal from .conftest import TESTS_DATA_DIR @@ -617,9 +618,59 @@ def test_bishop88_init_cond(method): NsVbi=NsVbi)) bad_results = np.isnan(vmp2) | (vmp2 < 0) | (err > 0.00001) assert not bad_results.any() - # test v_from_i + # test i_from_v imp2 = bishop88_i_from_v(vmp, *sde_params, d2mutau=d2mutau, NsVbi=NsVbi) err = np.abs(_sde_check_solution(imp2, vmp, *sde_params, d2mutau=d2mutau, NsVbi=NsVbi)) bad_results = np.isnan(imp2) | (imp2 < 0) | (err > 0.00001) assert not bad_results.any() + + +def test_batzelis_keypoints(): + params = {'photocurrent': 10, 'saturation_current': 1e-10, + 'resistance_series': 0.2, 'resistance_shunt': 3000, + 'nNsVth': 1.7} + + exact_values = { # calculated using pvlib.pvsystem.singlediode + 'i_sc': 9.999333377550565, + 'v_oc': 43.05589965219406, + 'i_mp': 9.513255314772051, + 'v_mp': 35.97259289596944, + 'p_mp': 342.21646055371264, + } + rtol = 5e-3 # accurate to within half a percent in this case + + output = batzelis_keypoints(**params) + for key in exact_values: + assert output[key] == pytest.approx(exact_values[key], rel=rtol) + + # numpy arrays + params2 = {k: np.array([v] * 2) for k, v in params.items()} + output2 = batzelis_keypoints(**params2) + for key in exact_values: + exp = np.array([exact_values[key]] * 2) + np.testing.assert_allclose(output2[key], exp, rtol=rtol) + + # pandas + params3 = {k: pd.Series(v) for k, v in params2.items()} + output3 = batzelis_keypoints(**params3) + assert isinstance(output3, pd.DataFrame) + for key in exact_values: + exp = pd.Series([exact_values[key]] * 2) + np.testing.assert_allclose(output3[key], exp, rtol=rtol) + + +def test_batzelis_keypoints_night(): + # SDMs produce photocurrent=0 and resistance_shunt=inf at night + out = batzelis_keypoints(photocurrent=0, saturation_current=1e-10, + resistance_series=0.2, resistance_shunt=np.inf, + nNsVth=1.7) + for k, v in out.items(): + assert v == 0, k # ensure all outputs are zero (not nan, etc) + + # test also when Rsh=inf but Iph > 0 + out = batzelis_keypoints(photocurrent=0.1, saturation_current=1e-10, + resistance_series=0.2, resistance_shunt=np.inf, + nNsVth=1.7) + for k, v in out.items(): + assert v > 0, k # ensure all outputs >0 (not nan, etc)