Add recombination current params to all bishop88 functions

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

@@ -14,6 +14,8 @@ Enhancements
 ~~~~~~~~~~~~
 * Created two new incidence angle modifier functions: :py:func:`pvlib.pvsystem.iam_martin_ruiz`
   and :py:func:`pvlib.pvsystem.iam_interp`. (:issue:`751`)
+* Added recombination current parameters to bishop88 single-diode functions and also
+  to :py:func:`pvlib.pvsystem.max_power_point`. (:issue:`762`)
 
 Bug fixes
 ~~~~~~~~~

pvlib/pvsystem.py

@@ -2473,7 +2473,8 @@ def singlediode(photocurrent, saturation_current, resistance_series,
 
 
 def max_power_point(photocurrent, saturation_current, resistance_series,
-                    resistance_shunt, nNsVth, method='brentq'):
+                    resistance_shunt, nNsVth, d2mutau=0, NsVbi=np.Inf,
+                    method='brentq'):
     """
     Given the single diode equation coefficients, calculates the maximum power
     point (MPP).
@@ -2491,6 +2492,14 @@ def max_power_point(photocurrent, saturation_current, resistance_series,
     nNsVth : numeric
         product of thermal voltage ``Vth`` [V], diode ideality factor ``n``,
         and number of serices cells ``Ns``
+    d2mutau : numeric
+        PVSyst thin-film recombination parameter that is the ratio of thickness
+        of the intrinsic layer squared :math:`d^2` and the diffusion length of
+        charge carriers :math:`\\mu \\tau`, in volts [V], defaults to 0[V]
+    NsVbi : numeric
+        PVSyst thin-film recombination parameter that is the product of the PV
+        module number of series cells ``Ns`` and the builtin voltage ``Vbi`` of
+        the intrinsic layer, in volts [V], defaults to ``np.inf``
     method : str
         either ``'newton'`` or ``'brentq'``
 
@@ -2508,7 +2517,8 @@ def max_power_point(photocurrent, saturation_current, resistance_series,
     """
     i_mp, v_mp, p_mp = _singlediode.bishop88_mpp(
         photocurrent, saturation_current, resistance_series,
-        resistance_shunt, nNsVth, method=method.lower()
+        resistance_shunt, nNsVth, d2mutau=0, NsVbi=np.Inf,
+        method=method.lower()
     )
     if isinstance(photocurrent, pd.Series):
         ivp = {'i_mp': i_mp, 'v_mp': v_mp, 'p_mp': p_mp}

pvlib/singlediode.py

@@ -173,7 +173,7 @@ def bishop88(diode_voltage, photocurrent, saturation_current,
 
 def bishop88_i_from_v(voltage, photocurrent, saturation_current,
                       resistance_series, resistance_shunt, nNsVth,
-                      method='newton'):
+                      d2mutau=0, NsVbi=np.Inf, method='newton'):
     """
     Find current given any voltage.
 
@@ -192,6 +192,14 @@ def bishop88_i_from_v(voltage, photocurrent, saturation_current,
     nNsVth : numeric
         product of diode ideality factor (n), number of series cells (Ns), and
         thermal voltage (Vth = k_b * T / q_e) in volts [V]
+    d2mutau : numeric
+        PVSyst thin-film recombination parameter that is the ratio of thickness
+        of the intrinsic layer squared :math:`d^2` and the diffusion length of
+        charge carriers :math:`\\mu \\tau`, in volts [V], defaults to 0[V]
+    NsVbi : numeric
+        PVSyst thin-film recombination parameter that is the product of the PV
+        module number of series cells ``Ns`` and the builtin voltage ``Vbi`` of
+        the intrinsic layer, in volts [V], defaults to ``np.inf``
     method : str
         one of two optional search methods: either ``'brentq'``, a reliable and
         bounded method or ``'newton'`` which is the default.
@@ -203,7 +211,7 @@ def bishop88_i_from_v(voltage, photocurrent, saturation_current,
     """
     # collect args
     args = (photocurrent, saturation_current, resistance_series,
-            resistance_shunt, nNsVth)
+            resistance_shunt, nNsVth, d2mutau, NsVbi)
 
     def fv(x, v, *a):
         # calculate voltage residual given diode voltage "x"
@@ -216,8 +224,9 @@ def fv(x, v, *a):
         # 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
-        def vd_from_brent(voc, v, iph, isat, rs, rsh, gamma):
-            return brentq(fv, 0.0, voc, args=(v, iph, isat, rs, rsh, gamma))
+        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))
 
         vd_from_brent_vectorized = np.vectorize(vd_from_brent)
         vd = vd_from_brent_vectorized(voc_est, voltage, *args)
@@ -235,7 +244,7 @@ def vd_from_brent(voc, v, iph, isat, rs, rsh, gamma):
 
 def bishop88_v_from_i(current, photocurrent, saturation_current,
                       resistance_series, resistance_shunt, nNsVth,
-                      method='newton'):
+                      d2mutau=0, NsVbi=np.Inf, method='newton'):
     """
     Find voltage given any current.
 
@@ -254,6 +263,14 @@ def bishop88_v_from_i(current, photocurrent, saturation_current,
     nNsVth : numeric
         product of diode ideality factor (n), number of series cells (Ns), and
         thermal voltage (Vth = k_b * T / q_e) in volts [V]
+    d2mutau : numeric
+        PVSyst thin-film recombination parameter that is the ratio of thickness
+        of the intrinsic layer squared :math:`d^2` and the diffusion length of
+        charge carriers :math:`\\mu \\tau`, in volts [V], defaults to 0[V]
+    NsVbi : numeric
+        PVSyst thin-film recombination parameter that is the product of the PV
+        module number of series cells ``Ns`` and the builtin voltage ``Vbi`` of
+        the intrinsic layer, in volts [V], defaults to ``np.inf``
     method : str
         one of two optional search methods: either ``'brentq'``, a reliable and
         bounded method or ``'newton'`` which is the default.
@@ -265,7 +282,7 @@ def bishop88_v_from_i(current, photocurrent, saturation_current,
     """
     # collect args
     args = (photocurrent, saturation_current, resistance_series,
-            resistance_shunt, nNsVth)
+            resistance_shunt, nNsVth, d2mutau, NsVbi)
     # first bound the search using voc
     voc_est = estimate_voc(photocurrent, saturation_current, nNsVth)
 
@@ -277,8 +294,9 @@ def fi(x, i, *a):
         # 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
-        def vd_from_brent(voc, i, iph, isat, rs, rsh, gamma):
-            return brentq(fi, 0.0, voc, args=(i, iph, isat, rs, rsh, gamma))
+        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))
 
         vd_from_brent_vectorized = np.vectorize(vd_from_brent)
         vd = vd_from_brent_vectorized(voc_est, current, *args)
@@ -295,7 +313,8 @@ def vd_from_brent(voc, i, iph, isat, rs, rsh, gamma):
 
 
 def bishop88_mpp(photocurrent, saturation_current, resistance_series,
-                 resistance_shunt, nNsVth, method='newton'):
+                 resistance_shunt, nNsVth, d2mutau=0, NsVbi=np.Inf,
+                 method='newton'):
     """
     Find max power point.
 
@@ -312,6 +331,14 @@ def bishop88_mpp(photocurrent, saturation_current, resistance_series,
     nNsVth : numeric
         product of diode ideality factor (n), number of series cells (Ns), and
         thermal voltage (Vth = k_b * T / q_e) in volts [V]
+    d2mutau : numeric
+        PVSyst thin-film recombination parameter that is the ratio of thickness
+        of the intrinsic layer squared :math:`d^2` and the diffusion length of
+        charge carriers :math:`\\mu \\tau`, in volts [V], defaults to 0[V]
+    NsVbi : numeric
+        PVSyst thin-film recombination parameter that is the product of the PV
+        module number of series cells ``Ns`` and the builtin voltage ``Vbi`` of
+        the intrinsic layer, in volts [V], defaults to ``np.inf``
     method : str
         one of two optional search methods: either ``'brentq'``, a reliable and
         bounded method or ``'newton'`` which is the default.
@@ -324,7 +351,7 @@ def bishop88_mpp(photocurrent, saturation_current, resistance_series,
     """
     # collect args
     args = (photocurrent, saturation_current, resistance_series,
-            resistance_shunt, nNsVth)
+            resistance_shunt, nNsVth, d2mutau, NsVbi)
     # first bound the search using voc
     voc_est = estimate_voc(photocurrent, saturation_current, nNsVth)
 
@@ -334,8 +361,9 @@ def fmpp(x, *a):
     if method.lower() == 'brentq':
         # break out arguments for numpy.vectorize to handle broadcasting
         vec_fun = np.vectorize(
-            lambda voc, iph, isat, rs, rsh, gamma:
-                brentq(fmpp, 0.0, voc, args=(iph, isat, rs, rsh, gamma))
+            lambda voc, iph, isat, rs, rsh, gamma, d2mutau, NsVbi:
+                brentq(fmpp, 0.0, voc,
+                       args=(iph, isat, rs, rsh, gamma, d2mutau, NsVbi))
         )
         vd = vec_fun(voc_est, *args)
     elif method.lower() == 'newton':

pvlib/test/test_singlediode.py

@@ -4,7 +4,8 @@
 
 import numpy as np
 from pvlib import pvsystem
-from pvlib.singlediode import bishop88, estimate_voc, VOLTAGE_BUILTIN
+from pvlib.singlediode import (bishop88_mpp, estimate_voc, VOLTAGE_BUILTIN,
+                               bishop88, bishop88_i_from_v, bishop88_v_from_i)
 import pytest
 from conftest import requires_scipy
 
@@ -153,9 +154,13 @@ def get_pvsyst_fs_495():
         'temp_ref': 25, 'irrad_ref': 1000, 'I_L_ref': 1.5743233463848496
     }
 
+# DeSoto @(888[W/m**2], 55[degC]) = {Pmp: 72.71, Isc: 1.402, Voc: 75.42)
 
+
+@requires_scipy
 @pytest.mark.parametrize(
     'poa, temp_cell, expected, tol', [
+        # reference conditions
         (
             get_pvsyst_fs_495()['irrad_ref'],
             get_pvsyst_fs_495()['temp_ref'],
@@ -167,8 +172,19 @@ def get_pvsyst_fs_495():
             },
             (5e-4, 0.04)
         ),
-        (POA, TCELL, {'pmp': 76.26, 'isc': 1.387, 'voc': 79.29}, (1e-3, 1e-3))]
-)  # DeSoto @(888[W/m**2], 55[degC]) = {Pmp: 72.71, Isc: 1.402, Voc: 75.42)
+        # other conditions
+        (
+            POA,
+            TCELL,
+            {
+                'pmp': 76.262,
+                'isc': 1.3868,
+                'voc': 79.292
+            },
+            (1e-4, 1e-4)
+        )
+    ]
+)
 def test_pvsyst_recombination_loss(poa, temp_cell, expected, tol):
     """test PVSst recombination loss"""
     pvsyst_fs_495 = get_pvsyst_fs_495()
@@ -199,9 +215,46 @@ def test_pvsyst_recombination_loss(poa, temp_cell, expected, tol):
     )
     # test max power
     assert np.isclose(max(pvsyst[2]), expected['pmp'], *tol)
+
     # test short circuit current
     isc_pvsyst = np.interp(0, pvsyst[1], pvsyst[0])
     assert np.isclose(isc_pvsyst, expected['isc'], *tol)
-    # test open circuit current
+
+    # test open circuit voltage
     voc_pvsyst = np.interp(0, pvsyst[0][::-1], pvsyst[1][::-1])
     assert np.isclose(voc_pvsyst, expected['voc'], *tol)
+
+    # repeat tests as above with specialized bishop88 functions
+    y = dict(d2mutau=pvsyst_fs_495['d2mutau'],
+             NsVbi=VOLTAGE_BUILTIN*pvsyst_fs_495['cells_in_series'])
+
+    mpp_88 = bishop88_mpp(*x, **y, method='newton')
+    assert np.isclose(mpp_88[2], expected['pmp'], *tol)
+
+    isc_88 = bishop88_i_from_v(0, *x, **y, method='newton')
+    assert np.isclose(isc_88, expected['isc'], *tol)
+
+    voc_88 = bishop88_v_from_i(0, *x, **y, method='newton')
+    assert np.isclose(voc_88, expected['voc'], *tol)
+
+    ioc_88 = bishop88_i_from_v(voc_88, *x, **y, method='newton')
+    assert np.isclose(ioc_88, 0.0, *tol)
+
+    vsc_88 = bishop88_v_from_i(isc_88, *x, **y, method='newton')
+    assert np.isclose(vsc_88, 0.0, *tol)
+
+    # now repeat with brentq
+    mpp_88 = bishop88_mpp(*x, **y, method='brentq')
+    assert np.isclose(mpp_88[2], expected['pmp'], *tol)
+
+    isc_88 = bishop88_i_from_v(0, *x, **y, method='brentq')
+    assert np.isclose(isc_88, expected['isc'], *tol)
+
+    voc_88 = bishop88_v_from_i(0, *x, **y, method='brentq')
+    assert np.isclose(voc_88, expected['voc'], *tol)
+
+    ioc_88 = bishop88_i_from_v(voc_88, *x, **y, method='brentq')
+    assert np.isclose(ioc_88, 0.0, *tol)
+
+    vsc_88 = bishop88_v_from_i(isc_88, *x, **y, method='brentq')
+    assert np.isclose(vsc_88, 0.0, *tol)