ENH: implementing pvsyst recombination loss current for CdTe and a:Si (#504)

* implementing pvsyst recombination loss current for CdTe and a:Si

* closes #163
* add constant VOLTAGE_BUILTIN = 0.9 (V) to singlediode_methods.py
* add arguments d2mutau=0, voltage_builtin=VOLTAGE_BUILTIN, and
cells_in_series to calculate i_recomb and v_recomb
* add check for d2mutau to calc i_recomb, v_recomb, and gradients
* add terms to current and gradients

* ENH: need to set grad_i_recomb if d2mutau not given

* also set v_recomb to Inf, so that 1/v_recomb is always zero
* but precalculate grad_2i_recomb only if d2mutau > 0, otherwise, use
zero

* ENH: fix need to multiply vbi by Ns, add test for recombination

* WIP: ENH: TST: BUG: in test compare desoto vs. pvsyst

* ENH: respond to comments

* change order of bishop88 args: cells_in_series to be first before
d2mutau and voltage_builtin
* hardcode PVSYST_FS495 coefficients
* don't compare pvsyst to desoto, instead compare pvsyst to fixed values
at two conditions: reference and 888[W/m^2],55[degC]
* rename test to be agnostic to module SKU
* remove testing script if __name__ == '__main__': section, not for
release

* DOC: ENH: update bishop88 documentation for pvsyst thin-film recombination loss

* DOC: ENH: update what's new

* DOC: add utf-8 file encoding

* DOC: ENH: BUG: minor doc fixes and formatting

* need to indent bullets in warning
* only need cells in series for PVSyst thin-film recombination loss
* escape math latex twice inside docstrings without raw prefix, eg:
\\tau without extra escape is a tab character, ha ha ha
* in rst, italics is a single asterisk, not underscore, thats markdown

* ENH: improve code style, fix d2mutau is array

* expand module type abbrev. in warnings
* also expand EG to For example, too terse
* cells in series can be int _or_ `None`
* since d2mutau and vbi are arrays, then ambiguous to test for zero, so
use np.where to assign conditionally instead (in two places)

* ENH: TST: use fixture instead of module-level dictionary for test

* dangerous to use module-level dictionary inside test, b/c dictionary
values could be changed inadvertently, but not a fixture
* add fixture as argument in test
* change parametrize to call fixture instead of module-level dictionary

* ENH: replace voltage_builtin with NsVbi and remove cells_in_series

* related to #516, the use of cells_in_series here is problematic
* to calculate arrays of different cells or modules with different
parameters we use np.where(is_recomb, ..., ...) 3d0deb2 but what
should cells_in_series be if d2mutau is zero? doesn't make physical
sense for it to be zero, but None doesn't work with numbers
* remove cells_in_series arg from bishop88
* remove last warning about setting cells_in_series correctly, move it
to Notes, but worded differently in the context of setting NsVbi
* replace voltage_builtin with NsVbi, default to np.inf
* fix docstring parameter descriptions

Author: mikofski

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

@@ -58,6 +58,7 @@ Enhancements
   :func:`~pvlib.singlediode_methods.bishop88_i_from_v`,
   :func:`~pvlib.singlediode_methods.bishop88_v_from_i`, and
   :func:`~pvlib.singlediode_methods.bishop88_mpp`.
+* Add PVSyst thin-film recombination losses for CdTe and a:Si (:issue:`163`)
 * Python 3.7 officially supported. (:issue:`496`)
 
 

pvlib/singlediode_methods.py

@@ -1,3 +1,4 @@
+# -*- coding: utf-8 -*-
 """
 Low-level functions for solving the single diode equation.
 """
@@ -22,6 +23,9 @@
 # rename newton and set keyword arguments
 newton = partial(_array_newton, tol=1e-6, maxiter=100, fprime2=None)
 
+# intrinsic voltage per cell junction for a:Si, CdTe, Mertens et al.
+VOLTAGE_BUILTIN = 0.9  # [V]
+
 
 def estimate_voc(photocurrent, saturation_current, nNsVth):
     """
@@ -62,14 +66,16 @@ def estimate_voc(photocurrent, saturation_current, nNsVth):
 
 
 def bishop88(diode_voltage, photocurrent, saturation_current,
-             resistance_series, resistance_shunt, nNsVth, gradients=False):
+             resistance_series, resistance_shunt, nNsVth, d2mutau=0,
+             NsVbi=np.Inf, gradients=False):
     """
     Explicit calculation of points on the IV curve described by the single
-    diode equation [1].
+    diode equation [1]_.
 
-    [1] "Computer simulation of the effects of electrical mismatches in
-    photovoltaic cell interconnection circuits" JW Bishop, Solar Cell (1988)
-    https://doi.org/10.1016/0379-6787(88)90059-2
+    .. warning::
+       * Do not use ``d2mutau`` with CEC coefficients.
+       * Usage of ``d2mutau`` with PVSyst coefficients is required for cadmium-
+         telluride (CdTe) and amorphous-silicon (a:Si) PV modules only.
 
     Parameters
     ----------
@@ -86,6 +92,14 @@ def bishop88(diode_voltage, photocurrent, saturation_current,
     nNsVth : numeric
         product of thermal voltage ``Vth`` [V], diode ideality factor ``n``,
         and number of series 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``
     gradients : bool
         False returns only I, V, and P. True also returns gradients
 
@@ -96,22 +110,56 @@ def bishop88(diode_voltage, photocurrent, saturation_current,
         :math:`\\frac{dI}{dV_d}`, :math:`\\frac{dV}{dV_d}`,
         :math:`\\frac{dI}{dV}`, :math:`\\frac{dP}{dV}`, and
         :math:`\\frac{d^2 P}{dV dV_d}`
+
+    Notes
+    -----
+    The PVSyst thin-film recombination losses parameters ``d2mutau`` and
+    ``NsVbi`` are only applied to cadmium-telluride (CdTe) and amorphous-
+    silicon (a:Si) PV modules, [2]_, [3]_. The builtin voltage :math:`V_{bi}`
+    should account for all junctions. For example: tandem and triple junction
+    cells would have builtin voltages of 1.8[V] and 2.7[V] respectively, based
+    on the default of 0.9[V] for a single junction. The parameter ``NsVbi``
+    should only account for the number of series cells in a single parallel
+    sub-string if the module has cells in parallel greater than 1.
+
+    References
+    ----------
+    .. [1] "Computer simulation of the effects of electrical mismatches in
+       photovoltaic cell interconnection circuits" JW Bishop, Solar Cell (1988)
+       :doi:`10.1016/0379-6787(88)90059-2`
+
+    .. [2] "Improved equivalent circuit and Analytical Model for Amorphous
+       Silicon Solar Cells and Modules." J. Mertens, et al., IEEE Transactions
+       on Electron Devices, Vol 45, No 2, Feb 1998.
+       :doi:`10.1109/16.658676`
+
+    .. [3] "Performance assessment of a simulation model for PV modules of any
+       available technology", André Mermoud and Thibault Lejeune, 25th EUPVSEC,
+       2010
+       :doi:`10.4229/25thEUPVSEC2010-4BV.1.114`
     """
+    # calculate recombination loss current where d2mutau > 0
+    is_recomb = d2mutau > 0  # True where there is thin-film recombination loss
+    v_recomb = np.where(is_recomb, NsVbi - diode_voltage, np.inf)
+    i_recomb = np.where(is_recomb, photocurrent * d2mutau / v_recomb, 0)
     # calculate temporary values to simplify calculations
     v_star = diode_voltage / nNsVth  # non-dimensional diode voltage
     g_sh = 1.0 / resistance_shunt  # conductance
     i = (photocurrent - saturation_current * np.expm1(v_star)
-         - diode_voltage * g_sh)
+         - diode_voltage * g_sh - i_recomb)
     v = diode_voltage - i * resistance_series
     retval = (i, v, i*v)
     if gradients:
+        # calculate recombination loss current gradients where d2mutau > 0
+        grad_i_recomb = np.where(is_recomb, i_recomb / v_recomb, 0)
+        grad_2i_recomb = np.where(is_recomb, 2 * grad_i_recomb / v_recomb, 0)
         g_diode = saturation_current * np.exp(v_star) / nNsVth  # conductance
-        grad_i = -g_diode - g_sh  # di/dvd
+        grad_i = -g_diode - g_sh - grad_i_recomb # di/dvd
         grad_v = 1.0 - grad_i * resistance_series  # dv/dvd
         # dp/dv = d(iv)/dv = v * di/dv + i
         grad = grad_i / grad_v  # di/dv
         grad_p = v * grad + i  # dp/dv
-        grad2i = -g_diode / nNsVth  # d2i/dvd
+        grad2i = -g_diode / nNsVth - grad_2i_recomb  # d2i/dvd
         grad2v = -grad2i * resistance_series  # d2v/dvd
         grad2p = (
             grad_v * grad + v * (grad2i/grad_v - grad_i*grad2v/grad_v**2)

pvlib/test/test_singlediode_methods.py

@@ -4,6 +4,8 @@
 
 import numpy as np
 from pvlib import pvsystem
+from pvlib.singlediode_methods import bishop88, estimate_voc, VOLTAGE_BUILTIN
+import pytest
 from conftest import requires_scipy
 
 POA = 888
@@ -100,7 +102,7 @@ def test_brentq_spr_e20_327():
 
 @requires_scipy
 def test_brentq_fs_495():
-    """test pvsystem.singlediode with Brent method on SPR-E20-327"""
+    """test pvsystem.singlediode with Brent method on FS495"""
     fs_495 = CECMOD.First_Solar_FS_495
     x = pvsystem.calcparams_desoto(
         effective_irradiance=POA, temp_cell=TCELL,
@@ -125,3 +127,75 @@ def test_brentq_fs_495():
                                 method='lambertw')
     assert np.isclose(pvs_ixx, ixx)
     return isc, voc, imp, vmp, pmp, i, v, pvs
+
+
+@pytest.fixture
+def pvsyst_fs_495():
+    """
+    PVsyst parameters for First Solar FS-495 module from PVSyst-6.7.2 database.
+
+    I_L_ref derived from Isc_ref conditions::
+
+        I_L_ref = (I_sc_ref + Id + Ish) / (1 - d2mutau/(Vbi*N_s - Vd))
+
+    where::
+
+        Vd = I_sc_ref * R_s
+        Id = I_o_ref * (exp(Vd / nNsVt) - 1)
+        Ish = Vd / R_sh_ref
+
+    """
+    return {
+        'd2mutau': 1.31, 'alpha_sc': 0.00039, 'gamma_ref': 1.48,
+        'mu_gamma': 0.001, 'I_o_ref': 9.62e-10, 'R_sh_ref': 5000,
+        'R_sh_0': 12500, 'R_sh_exp': 3.1, 'R_s': 4.6, 'beta_oc': -0.2116,
+        'EgRef': 1.5, 'cells_in_series': 108, 'cells_in_parallel': 2,
+        'I_sc_ref': 1.55, 'V_oc_ref': 86.5, 'I_mp_ref': 1.4, 'V_mp_ref': 67.85,
+        'temp_ref': 25, 'irrad_ref': 1000, 'I_L_ref': 1.5743233463848496
+    }
+
+
+@pytest.mark.parametrize(
+    'poa, temp_cell, expected, tol',
+    [(pvsyst_fs_495()['irrad_ref'], pvsyst_fs_495()['temp_ref'],
+      {'pmp': pvsyst_fs_495()['I_mp_ref'] * pvsyst_fs_495()['V_mp_ref'],
+       'isc': pvsyst_fs_495()['I_sc_ref'], 'voc': pvsyst_fs_495()['V_oc_ref']},
+      (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)
+def test_pvsyst_recombination_loss(pvsyst_fs_495, poa, temp_cell, expected,
+                                   tol):
+    """test PVSst recombination loss"""
+    # first evaluate PVSyst model with thin-film recombination loss current
+    # at reference conditions
+    x = pvsystem.calcparams_pvsyst(
+        effective_irradiance=poa, temp_cell=temp_cell,
+        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']
+    )
+    il_pvsyst, io_pvsyst, rs_pvsyst, rsh_pvsyst, nnsvt_pvsyst = x
+    voc_est_pvsyst = estimate_voc(photocurrent=il_pvsyst,
+                                  saturation_current=io_pvsyst,
+                                  nNsVth=nnsvt_pvsyst)
+    vd_pvsyst = np.linspace(0, voc_est_pvsyst, 1000)
+    pvsyst = bishop88(
+        diode_voltage=vd_pvsyst, photocurrent=il_pvsyst,
+        saturation_current=io_pvsyst, resistance_series=rs_pvsyst,
+        resistance_shunt=rsh_pvsyst, nNsVth=nnsvt_pvsyst,
+        d2mutau=pvsyst_fs_495['d2mutau'],
+        NsVbi=VOLTAGE_BUILTIN*pvsyst_fs_495['cells_in_series']
+    )
+    # 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
+    voc_pvsyst = np.interp(0, pvsyst[0][::-1], pvsyst[1][::-1])
+    assert np.isclose(voc_pvsyst, expected['voc'], *tol)