inverter.py

"""
This module contains functions for inverter modeling and for fitting inverter
models to data.

Inverter models calculate AC power output from DC input. Model parameters
should be passed as a single dict.

Functions for estimating parameters for inverter models should follow the
naming pattern 'fit_<model name>', e.g., fit_sandia.

"""

import numpy as np
import pandas as pd
from numpy.polynomial.polynomial import polyfit  # different than np.polyfit


def _sandia_eff(v_dc, p_dc, inverter):
    r'''
    Calculate the inverter AC power without clipping
    '''
    Paco = inverter['Paco']
    Pdco = inverter['Pdco']
    Vdco = inverter['Vdco']
    C0 = inverter['C0']
    C1 = inverter['C1']
    C2 = inverter['C2']
    C3 = inverter['C3']
    Pso = inverter['Pso']

    A = Pdco * (1 + C1 * (v_dc - Vdco))
    B = Pso * (1 + C2 * (v_dc - Vdco))
    C = C0 * (1 + C3 * (v_dc - Vdco))

    return (Paco / (A - B) - C * (A - B)) * (p_dc - B) + C * (p_dc - B)**2


def _sandia_limits(power_ac, p_dc, Paco, Pnt, Pso):
    r'''
    Applies minimum and maximum power limits to `power_ac`
    '''
    power_ac = np.minimum(Paco, power_ac)
    min_ac_power = -1.0 * abs(Pnt)
    below_limit = p_dc < Pso
    try:
        power_ac[below_limit] = min_ac_power
    except TypeError:  # power_ac is a float
        if below_limit:
            power_ac = min_ac_power
    return power_ac


def sandia(v_dc, p_dc, inverter):
    r'''
    Convert DC power and voltage to AC power using Sandia's
    Grid-Connected PV Inverter model.

    Parameters
    ----------
    v_dc : numeric
        DC voltage input to the inverter. [V]

    p_dc : numeric
        DC power input to the inverter. [W]

    inverter : dict-like
        Defines parameters for the inverter model in [1]_.

    Returns
    -------
    power_ac : numeric
        AC power output. [W]

    Notes
    -----

    Determines the AC power output of an inverter given the DC voltage and DC
    power. Output AC power is bounded above by the parameter ``Paco``, to
    represent inverter "clipping".  When `power_ac` would be less than
    parameter ``Pso`` (startup power required), then `power_ac` is set to
    ``-Pnt``, representing self-consumption. `power_ac` is not adjusted for
    maximum power point tracking (MPPT) voltage windows or maximum current
    limits of the inverter.

    Required model parameters are:

    ======   ============================================================
    Column   Description
    ======   ============================================================
    Paco     AC power rating of the inverter. [W]
    Pdco     DC power input that results in Paco output at reference
             voltage Vdco. [W]
    Vdco     DC voltage at which the AC power rating is achieved
             with Pdco power input. [V]
    Pso      DC power required to start the inversion process, or
             self-consumption by inverter, strongly influences inverter
             efficiency at low power levels. [W]
    C0       Parameter defining the curvature (parabolic) of the
             relationship between AC power and DC power at the reference
             operating condition. [1/W]
    C1       Empirical coefficient allowing ``Pdco`` to vary linearly
             with DC voltage input. [1/V]
    C2       Empirical coefficient allowing ``Pso`` to vary linearly with
             DC voltage input. [1/V]
    C3       Empirical coefficient allowing ``C0`` to vary linearly with
             DC voltage input. [1/V]
    Pnt      AC power consumed by the inverter at night (night tare). [W]
    ======   ============================================================

    A copy of the parameter database from the System Advisor Model (SAM) [2]_
    is provided with pvlib and may be read  using
    :py:func:`pvlib.pvsystem.retrieve_sam`.

    References
    ----------
    .. [1] D. King, S. Gonzalez, G. Galbraith, W. Boyson, "Performance Model
       for Grid-Connected Photovoltaic Inverters", SAND2007-5036, Sandia
       National Laboratories.

    .. [2] System Advisor Model web page. https://sam.nrel.gov.

    See Also
    --------
    pvlib.pvsystem.retrieve_sam
    '''

    Paco = inverter['Paco']
    Pnt = inverter['Pnt']
    Pso = inverter['Pso']

    power_ac = _sandia_eff(v_dc, p_dc, inverter)
    power_ac = _sandia_limits(power_ac, p_dc, Paco, Pnt, Pso)

    if isinstance(p_dc, pd.Series):
        power_ac = pd.Series(power_ac, index=p_dc.index)

    return power_ac


def sandia_multi(v_dc, p_dc, inverter):
    r'''
    Convert DC power and voltage to AC power for an inverter with multiple
    MPPT inputs.

    Uses Sandia's Grid-Connected PV Inverter model [1]_. Extension of [1]_
    to inverters with multiple, unbalanced inputs as described in [2]_.

    Parameters
    ----------
    v_dc : tuple, list or array of numeric
        DC voltage on each MPPT input of the inverter. If type is array, must
        be 2d with axis 0 being the MPPT inputs. [V]

    p_dc : tuple, list or array of numeric
        DC power on each MPPT input of the inverter. If type is array, must
        be 2d with axis 0 being the MPPT inputs. [W]

    inverter : dict-like
        Defines parameters for the inverter model in [1]_.

    Returns
    -------
    power_ac : numeric
        AC power output for the inverter. [W]

    Raises
    ------
    ValueError
        If v_dc and p_dc have different lengths.

    Notes
    -----
    See :py:func:`pvlib.inverter.sandia` for definition of the parameters in
    `inverter`.

    References
    ----------
    .. [1] D. King, S. Gonzalez, G. Galbraith, W. Boyson, "Performance Model
       for Grid-Connected Photovoltaic Inverters", SAND2007-5036, Sandia
       National Laboratories.
    .. [2] C. Hansen, J. Johnson, R. Darbali-Zamora, N. Gurule. "Modeling
       Efficiency Of Inverters With Multiple Inputs", 49th IEEE Photovoltaic
       Specialist Conference, Philadelphia, PA, USA. June 2022.

    See Also
    --------
    pvlib.inverter.sandia
    '''

    if len(p_dc) != len(v_dc):
        raise ValueError('p_dc and v_dc have different lengths')
    power_dc = sum(p_dc)
    power_ac = 0. * power_dc

    for vdc, pdc in zip(v_dc, p_dc):
        power_ac += pdc / power_dc * _sandia_eff(vdc, power_dc, inverter)

    return _sandia_limits(power_ac, power_dc, inverter['Paco'],
                          inverter['Pnt'], inverter['Pso'])


def adr(v_dc, p_dc, inverter, vtol=0.10):
    r'''
    Converts DC power and voltage to AC power using Anton Driesse's
    grid-connected inverter efficiency model.

    Parameters
    ----------
    v_dc : numeric
        DC voltage input to the inverter, should be >= 0. [V]

    p_dc : numeric
        DC power input to the inverter, should be >= 0. [W]

    inverter : dict-like
        Defines parameters for the inverter model in [1]_.  See Notes for
        required model parameters. A parameter database is provided with pvlib
        and may be read using :py:func:`pvlib.pvsystem.retrieve_sam`.

    vtol : numeric, default 0.1
        Fraction of DC voltage that determines how far the efficiency model is
        extrapolated beyond the inverter's normal input voltage operating
        range. 0.0 <= vtol <= 1.0. [unitless]

    Returns
    -------
    power_ac : numeric
        AC power output. [W]

    Notes
    -----
    Determines the AC power output of an inverter given the DC voltage and DC
    power. Output AC power is bounded above by the parameter ``Pacmax``, to
    represent inverter "clipping". AC power is bounded below by ``-Pnt``
    (negative when power is consumed rather than produced) which represents
    self-consumption. `power_ac` is not adjusted for maximum power point
    tracking (MPPT) voltage windows or maximum current limits of the inverter.

    Required model parameters are:

    ================ ==========================================================
    Column           Description
    ================ ==========================================================
    Pnom             Nominal DC power, typically the DC power needed to produce
                     maximum AC power output. [W]
    Vnom             Nominal DC input voltage. Typically the level at which the
                     highest efficiency is achieved. [V]
    Vmax             Maximum DC input voltage. [V]
    Vmin             Minimum DC input voltage. [V]
    Vdcmax           Maximum voltage supplied from DC array. [V]
    MPPTHi           Maximum DC voltage for MPPT range. [V]
    MPPTLow          Minimum DC voltage for MPPT range. [V]
    Pacmax           Maximum AC output power, used to clip the output power
                     if needed. [W]
    ADRCoefficients  A list of 9 coefficients that capture the influence
                     of input voltage and power on inverter losses, and thereby
                     efficiency. Corresponds to terms from [1]_ (in order):
                     :math: `b_{0,0}, b_{1,0}, b_{2,0}, b_{0,1}, b_{1,1},
                     b_{2,1}, b_{0,2}, b_{1,2},  b_{2,2}`. See [1]_ for the
                     use of each coefficient and its associated unit.
    Pnt              AC power consumed by inverter at night (night tare) to
                     maintain circuitry required to sense the PV array
                     voltage. [W]
    ================ ==========================================================

    AC power output is set to NaN where the input DC voltage exceeds a limit
    M = max(Vmax, Vdcmax, MPPTHi) x (1 + vtol), and where the input DC voltage
    is less than a limit m = max(Vmin, MPPTLow) x (1 - vtol)

    References
    ----------
    .. [1] A. Driesse, "Beyond the Curves: Modeling the Electrical Efficiency
       of Photovoltaic Inverters", 33rd IEEE Photovoltaic Specialist
       Conference (PVSC), June 2008

    See Also
    --------
    pvlib.inverter.sandia
    pvlib.pvsystem.retrieve_sam
    '''

    p_nom = inverter['Pnom']
    v_nom = inverter['Vnom']
    pac_max = inverter['Pacmax']
    p_nt = inverter['Pnt']
    ce_list = inverter['ADRCoefficients']
    v_max = inverter['Vmax']
    v_min = inverter['Vmin']
    vdc_max = inverter['Vdcmax']
    mppt_hi = inverter['MPPTHi']
    mppt_low = inverter['MPPTLow']

    v_lim_upper = float(np.nanmax([v_max, vdc_max, mppt_hi]) * (1 + vtol))
    v_lim_lower = float(np.nanmax([v_min, mppt_low]) * (1 - vtol))

    pdc = p_dc / p_nom
    vdc = v_dc / v_nom
    # zero voltage will lead to division by zero, but since power is
    # set to night time value later, these errors can be safely ignored
    with np.errstate(invalid='ignore', divide='ignore'):
        poly = np.array([pdc**0,  # replace with np.ones_like?
                         pdc,
                         pdc**2,
                         vdc - 1,
                         pdc * (vdc - 1),
                         pdc**2 * (vdc - 1),
                         1. / vdc - 1,  # divide by 0
                         pdc * (1. / vdc - 1),  # invalid 0./0. --> nan
                         pdc**2 * (1. / vdc - 1)])  # divide by 0
    p_loss = np.dot(np.array(ce_list), poly)
    power_ac = p_nom * (pdc - p_loss)
    p_nt = -1 * np.absolute(p_nt)

    # set output to nan where input is outside of limits
    # errstate silences case where input is nan
    with np.errstate(invalid='ignore'):
        invalid = (v_lim_upper < v_dc) | (v_dc < v_lim_lower)
    power_ac = np.where(invalid, np.nan, power_ac)

    # set night values
    power_ac = np.where(vdc == 0, p_nt, power_ac)
    power_ac = np.maximum(power_ac, p_nt)

    # set max ac output
    power_ac = np.minimum(power_ac, pac_max)

    if isinstance(p_dc, pd.Series):
        power_ac = pd.Series(power_ac, index=pdc.index)

    return power_ac


def pvwatts(pdc, pdc0, eta_inv_nom=0.96, eta_inv_ref=0.9637):
    r"""
    NREL's PVWatts inverter model.

    The PVWatts inverter model [1]_ calculates inverter efficiency :math:`\eta`
    as a function of input DC power :math:`P_{dc}`

    .. math::

        \eta = \frac{\eta_{nom}}{\eta_{ref}} (-0.0162\zeta - \frac{0.0059}
        {\zeta} + 0.9858)

    where :math:`\zeta=P_{dc}/P_{dc0}` and :math:`P_{dc0}=P_{ac0}/\eta_{nom}`.

    Output AC power is then given by

    .. math::

        P_{ac} = \min(\eta P_{dc}, P_{ac0})

    Parameters
    ----------
    pdc : numeric
        DC power. Same unit as ``pdc0``.
    pdc0 : numeric
        DC input limit of the inverter.  Same unit as ``pdc``.
    eta_inv_nom : numeric, default 0.96
        Nominal inverter efficiency. [unitless]
    eta_inv_ref : numeric, default 0.9637
        Reference inverter efficiency. PVWatts defines it to be 0.9637
        and is included here for flexibility. [unitless]

    Returns
    -------
    power_ac : numeric
        AC power.  Same unit as ``pdc0``.

    Notes
    -----
    When sourcing ``pdc`` from pvlib functions
    (e.g. :py:func:`pvlib.pvsystem.pvwatts_dc`) their DC power output is in W,
    and ``pdc0`` should have the same unit (W).

    Note that ``pdc0`` is also used as a symbol in
    :py:func:`pvlib.pvsystem.pvwatts_dc`. ``pdc0`` in this function refers to
    the DC power input limit of the inverter. ``pdc0`` in
    :py:func:`pvlib.pvsystem.pvwatts_dc` refers to the DC power of the modules
    at reference conditions.

    See Also
    --------
    pvlib.inverter.pvwatts_multi

    References
    ----------
    .. [1] A. P. Dobos, "PVWatts Version 5 Manual,"
       http://pvwatts.nrel.gov/downloads/pvwattsv5.pdf (2014).
    """

    pac0 = eta_inv_nom * pdc0
    zeta = pdc / pdc0

    # arrays to help avoid divide by 0 for scalar and array
    eta = np.zeros_like(pdc, dtype=float)
    pdc_neq_0 = ~np.equal(pdc, 0)

    # eta < 0 if zeta < 0.006. power_ac is forced to be >= 0 below. GH 541
    # In some published versions of [1] the parentheses are missing
    eta = eta_inv_nom / eta_inv_ref * (
        -0.0162 * zeta - np.divide(0.0059, zeta, out=eta, where=pdc_neq_0)
        + 0.9858)  # noQA: W503

    power_ac = eta * pdc
    power_ac = np.minimum(pac0, power_ac)
    power_ac = np.maximum(0, power_ac)     # GH 541

    return power_ac


def pvwatts_multi(pdc, pdc0, eta_inv_nom=0.96, eta_inv_ref=0.9637):
    r"""
    Extend NREL's PVWatts inverter model for multiple MPP inputs.

    DC input power is summed over MPP inputs to obtain the DC power
    input to the PVWatts inverter model. See :py:func:`pvlib.inverter.pvwatts`
    for details.

    Parameters
    ----------
    pdc : tuple, list or array of numeric
        DC power on each MPPT input of the inverter. If type is array, must
        be 2d with axis 0 being the MPPT inputs. Same unit as ``pdc0``.
    pdc0 : numeric
        Total DC power limit of the inverter.  Same unit as ``pdc``.
    eta_inv_nom : numeric, default 0.96
        Nominal inverter efficiency. [unitless]
    eta_inv_ref : numeric, default 0.9637
        Reference inverter efficiency. PVWatts defines it to be 0.9637
        and is included here for flexibility. [unitless]

    Returns
    -------
    power_ac : numeric
        AC power.  Same unit as ``pdc0``.

    See Also
    --------
    pvlib.inverter.pvwatts
    """
    return pvwatts(sum(pdc), pdc0, eta_inv_nom, eta_inv_ref)


def fit_sandia(ac_power, dc_power, dc_voltage, dc_voltage_level, p_ac_0, p_nt):
    r'''
    Determine parameters for the Sandia inverter model.

    Parameters
    ----------
    ac_power : array_like
        AC power output at each data point [W].
    dc_power : array_like
        DC power input at each data point [W].
    dc_voltage : array_like
        DC input voltage at each data point [V].
    dc_voltage_level : array_like
        DC input voltage level at each data point. Values must be 'Vmin',
        'Vnom' or 'Vmax'.
    p_ac_0 : float
        Rated AC power of the inverter [W].
    p_nt : float
        Night tare, i.e., power consumed while inverter is not delivering
        AC power. [W]

    Returns
    -------
    dict
        A set of parameters for the Sandia inverter model [1]_. See
        :py:func:`pvlib.inverter.sandia` for a description of keys and values.

    See Also
    --------
    pvlib.inverter.sandia

    Notes
    -----
    The fitting procedure to estimate parameters is described at [2]_.
    A data point is a pair of values (dc_power, ac_power). Typically, inverter
    performance is measured or described at three DC input voltage levels,
    denoted 'Vmin', 'Vnom' and 'Vmax' and at each level, inverter efficiency
    is determined at various output power levels. For example,
    the CEC inverter test protocol [3]_ specifies measurement of input DC
    power that delivers AC output power of 0.1, 0.2, 0.3, 0.5, 0.75 and 1.0 of
    the inverter's AC power rating.

    References
    ----------
    .. [1] D. King, S. Gonzalez, G. Galbraith, W. Boyson, "Performance Model
       for Grid-Connected Photovoltaic Inverters", SAND2007-5036, Sandia
       National Laboratories.
    .. [2] Sandia Inverter Model page, PV Performance Modeling Collaborative
       https://pvpmc.sandia.gov/modeling-steps/dc-to-ac-conversion/sandia-inverter-model/
    .. [3] W. Bower, et al., "Performance Test Protocol for Evaluating
       Inverters Used in Grid-Connected Photovoltaic Systems", available at
       https://www.energy.ca.gov/sites/default/files/2020-06/2004-11-22_Sandia_Test_Protocol_ada.pdf
    '''  # noqa: E501

    voltage_levels = ['Vmin', 'Vnom', 'Vmax']

    # average dc input voltage at each voltage level
    v_d = np.array(
        [dc_voltage[dc_voltage_level == 'Vmin'].mean(),
         dc_voltage[dc_voltage_level == 'Vnom'].mean(),
         dc_voltage[dc_voltage_level == 'Vmax'].mean()])
    v_nom = v_d[1]  # model parameter
    # independent variable for regressions, x_d
    x_d = v_d - v_nom

    # empty dataframe to contain intermediate variables
    coeffs = pd.DataFrame(index=voltage_levels,
                          columns=['a', 'b', 'c', 'p_dc', 'p_s0'], data=np.nan)

    def solve_quad(a, b, c):
        return (-b + (b**2 - 4 * a * c)**.5) / (2 * a)

    # [2] STEP 3E, fit a line to (DC voltage, model_coefficient)
    def extract_c(x_d, add):
        beta0, beta1 = polyfit(x_d, add, 1)
        c = beta1 / beta0
        return beta0, beta1, c

    for d in voltage_levels:
        x = dc_power[dc_voltage_level == d]
        y = ac_power[dc_voltage_level == d]
        # [2] STEP 3B
        # fit a quadratic to (DC power, AC power)
        c, b, a = polyfit(x, y, 2)

        # [2] STEP 3D, solve for p_dc and p_s0
        p_dc = solve_quad(a, b, (c - p_ac_0))
        p_s0 = solve_quad(a, b, c)

        # Add values to dataframe at index d
        coeffs['a'][d] = a
        coeffs['p_dc'][d] = p_dc
        coeffs['p_s0'][d] = p_s0

    b_dc0, b_dc1, c1 = extract_c(x_d, coeffs['p_dc'])
    b_s0, b_s1, c2 = extract_c(x_d, coeffs['p_s0'])
    b_c0, b_c1, c3 = extract_c(x_d, coeffs['a'])

    p_dc0 = b_dc0
    p_s0 = b_s0
    c0 = b_c0

    # prepare dict and return
    return {'Paco': p_ac_0, 'Pdco': p_dc0, 'Vdco': v_nom, 'Pso': p_s0,
            'C0': c0, 'C1': c1, 'C2': c2, 'C3': c3, 'Pnt': p_nt}