#!/usr/bin/env python
# coding: utf-8

# # Explicit methods for solving IV curves of a PV system
# 
# This is the last notebook in a series exploring how to model PV systems using the single diode model (SDM). The first two notebooks presented implicit approaches that return either the corresponding current for a specified voltage or the current and voltage of the max power point. This notebook uses an explicit method which parameterizes the SDM using the diode voltage, and returns a current and voltage corresponding to the value of the parameter. The diode voltage is defined as the voltage across the diode in the SDM and is related to the current and voltage by the following relation, $V_{diode} = V + I R_s$, in which $R_s$ is the series resistance.
# 
# ## PVMismatch
# [PVMismatch](https://sunpower.github.io/PVMismatch/) is a Python package that solves IV curves using this explicit approach which was discussed in by J. Bishop in _Solar_ (1980). There is another notebook in this series that is a repeat of this notebook, but using PVMismatch.
# 
# ## Conclusions
# Overall, it appears that the explicit approach is about 5x faster than the implicit approaches.

# In[1]:


get_ipython().run_line_magic('matplotlib', 'inline')

import numpy as np
from matplotlib import pyplot as plt


# In[2]:


# CONSTANTS
IL = 7.0  # [A] photogenerated current
I0 = 1.0e-6  # [A] dark current
RSH = 20.0  # [ohms] shunt resistance
RS = 0.001  # [ohms] series resistance
GAMMA = 1.23  # diode ideality
VT = 0.026  # [V] thermal voltage at 300[K]
VBYPASS = -0.5  # bypass diode trigger voltage [V]
NSERIES_CELLS = 72  # number of series cells per module
NMODS = 8  # number of modules per string
NSTRINGS = 3  # number of strings in system
NBYPASS = 3  # number of bypass diodes per module, assumes cells divided evenly


# In[3]:


def sdm(il, i0, rsh, rs, gamma, vt, vdiode):
    """
    Residual for single diode model (SDM) and derivative.

    Args:
        il (numeric): photogenerated current [A]
        i0 (numeric): dark current [A]
        rsh (numeric): shunt resistance [Ohms]
        rs (numeric): series resistance [Ohms]
        gamma (numeric): diode ideality factor
        vt (numeric): thermal voltage [V]
        vdiode (numeric): diode voltage [V]

    Returns:
        vcell (numeric): cell voltages [V]
        icell (numeric): cell currents [A]
    """
    icell = il - i0*(np.exp(vdiode / gamma / vt) - 1.0) - vdiode/rsh
    vcell = vdiode - icell*rs  # voltage drop across diode in SDM [V]
    
    return icell, vcell


# In[4]:


# add some noise to the irradiance and temperature
noise = np.random.randn(NSERIES_CELLS, NMODS, NSTRINGS) / 1000


# In[5]:


# estimate upper limit of system voltage
VOC_EST = (GAMMA*VT * np.log(IL / I0 + 1.0))*NMODS*NSERIES_CELLS
VOC_EST


# In[6]:


# make a distribution that is more closely spaced around Vmp to Voc
# and less spaced near Isc, also in reverse order, so current increases
NPTS = 1000
VDIODE = (VOC_EST-np.logspace(0, np.log10(VOC_EST), NPTS))/(NMODS*NSERIES_CELLS)
plt.plot(VDIODE)
plt.title('diode voltage is log spaced')
plt.xlabel('index')
plt.ylabel('voltage [V]')
plt.grid();


# In[7]:


# make the cells with some noise
iv = []
for pvstr in range(NSTRINGS):
    for pvmod in range(NMODS):
        for pvcell in range(NSERIES_CELLS):
            i, v = sdm(IL+noise[pvcell, pvmod, pvstr], I0, RSH, RS, GAMMA, VT, VDIODE)
            iv.append([i, v])
pvcells = np.array(iv).reshape(NSTRINGS, NMODS, NSERIES_CELLS, 2, NPTS)  # reshape


# In[8]:


# make a plot of the the cell-0, mod-0, str-0
plt.plot(pvcells[0, 0, 0][1, :], pvcells[0, 0, 0][0, :])
plt.title('iv curve cell-0, mod-0, str-0')
plt.ylabel('current [A]')
plt.xlabel('voltage [V]')
plt.grid();


# In[9]:


# combine substrings and check bypass diode activation
iv = []
ncells_per_sub = NSERIES_CELLS/NBYPASS
for pvstr in range(NSTRINGS):
    for pvmod in range(NMODS):
        for pvsub in range(NBYPASS):
            idx0 = int(pvsub*ncells_per_sub)  # 1st index of substring
            idx1 = int(idx0+ncells_per_sub)  # last indest of substring
            isub = pvcells[pvstr, pvmod, idx0][0, :]
            vsub = np.zeros_like(isub)
            for ivcell in pvcells[pvstr, pvmod, idx0:idx1]:
                vsub += np.interp(isub, ivcell[0, :], ivcell[1, :])
            vsub[vsub < VBYPASS] = VBYPASS
            iv.append([isub, vsub])
pvsubs = np.array(iv).reshape(NSTRINGS, NMODS, NBYPASS, 2, NPTS)


# In[10]:


# make a plot of the the substring-0, mod-0, str-0
plt.plot(pvsubs[0, 0, 0][1, :], pvsubs[0, 0, 0][0, :])
plt.title('iv curve substring-0, mod-0, str-0')
plt.ylabel('current [A]')
plt.xlabel('voltage [V]')
plt.grid();


# In[11]:


# combine substrings into modules
iv = []
for pvstr in range(NSTRINGS):
    for pvmod in range(NMODS):
        imod = pvsubs[pvstr, pvmod, 0][0, :]
        vmod = np.zeros_like(imod)
        for ivsub in pvsubs[pvstr, pvmod, :]:
            vmod += np.interp(imod, ivsub[0, :], ivsub[1, :])
        iv.append([imod, vmod])
pvmods = np.array(iv).reshape(NSTRINGS, NMODS, 2, NPTS)


# In[12]:


# make a plot of the the mod-0, str-0
plt.plot(pvmods[0, 0][1, :], pvmods[0, 0][0, :])
plt.title('iv curve mod-0, str-0')
plt.ylabel('current [A]')
plt.xlabel('voltage [V]')
plt.grid();


# In[13]:


# combine modules into strings
iv = []
for pvstr in range(NSTRINGS):
    istr = pvmods[pvstr, 0][0, :]
    vstr = np.zeros_like(istr)
    for ivmod in pvmods[pvstr, :]:
        vstr += np.interp(istr, ivmod[0, :], ivmod[1, :])
    iv.append([istr, vstr])
pvstrs = np.array(iv).reshape(NSTRINGS, 2, NPTS)


# In[14]:


# make a plot of the the str-0
plt.plot(pvstrs[0][1, :], pvstrs[0][0, :])
plt.title('iv curve str-0')
plt.ylabel('current [A]')
plt.xlabel('voltage [V]')
plt.grid();


# In[15]:


# combine string into system
vsys = pvstrs[0][1, ::-1]
isys = np.zeros_like(vsys)
for ivstr in pvstrs:
    isys += np.interp(vsys, ivstr[1, ::-1], ivstr[0, ::-1])
pvsys = np.array([isys, vsys, isys*vsys])


# In[16]:


# make a plot of the the system
f, ax = plt.subplots(1, 2, figsize=(12, 4))
ax[0].plot(pvsys[1, :], pvsys[0, :])
ax[0].set_title('iv curve system')
ax[0].set_ylabel('current [A]')
ax[0].set_xlabel('voltage [V]')
ax[0].grid()
# power
ax[1].plot(pvsys[1, :], pvsys[2, :])
ax[1].set_title('iv curve system')
ax[1].set_ylabel('power [W]')
ax[1].set_xlabel('voltage [V]')
ax[1].grid();
plt.tight_layout()


# In[17]:


# approximate max power, depends on NPTS and spacing
Pmp = np.max(pvsys[2, :])
Pmp


# In[18]:


# approximate location of MPP, depends on NPTS and spacing
mpp = np.argmax(pvsys[2, :])
Imp, Vmp, mpp = pvsys[:, mpp]
Imp, Vmp, mpp