diff --git a/docs/sphinx/source/api.rst b/docs/sphinx/source/api.rst index 8ddbd40f59..4070ed56c3 100644 --- a/docs/sphinx/source/api.rst +++ b/docs/sphinx/source/api.rst @@ -201,10 +201,10 @@ Functions relevant for the single diode model. .. autosummary:: :toctree: generated/ - pvsystem.singlediode pvsystem.calcparams_desoto - pvsystem.v_from_i pvsystem.i_from_v + pvsystem.singlediode + pvsystem.v_from_i SAPM model ---------- diff --git a/docs/sphinx/source/whatsnew/v0.5.1.rst b/docs/sphinx/source/whatsnew/v0.5.1.rst index f889b676a3..bf928337dc 100644 --- a/docs/sphinx/source/whatsnew/v0.5.1.rst +++ b/docs/sphinx/source/whatsnew/v0.5.1.rst @@ -16,18 +16,35 @@ Bug fixes Enhancements ~~~~~~~~~~~~ * Improve clearsky.lookup_linke_turbidity speed. (:issue:`368`) +* Ideal devices supported in single diode model, e.g., + resistance_series = 0 and/or resistance_shunt = numpy.inf (:issue:`340`) +* `pvsystem.v_from_i` and `pvsystem.i_from_v` computations for near ideal + devices are more numerically stable. However, very, very near ideal + resistance_series and/or resistance_shunt may still cause issues with the + implicit solver (:issue:`340`) + +API Changes +~~~~~~~~~~~ +* `pvsystem.v_from_i` and `pvsystem.i_from_v` functions now accept + resistance_series = 0 and/or resistance_shunt = numpy.inf as inputs + (:issue:`340`) Documentation ~~~~~~~~~~~~~ * Doc string of modelchain.basic_chain was updated to describe args more accurately +* Doc strings of `singlediode`, `pvsystem.v_from_i`, and `pvsystem.i_from_v` + were updated to describe acceptable input arg ranges Testing ~~~~~~~ * Changed test for clearsky.haurwitz to operate on zenith angles +* Significant new test cases added for `pvsystem.v_from_i` and + `pvsystem.i_from_v` (:issue:`340`) Contributors ~~~~~~~~~~~~ * Cliff Hansen * KonstantinTr * Will Holmgren +* Mark Campanelli diff --git a/pvlib/atmosphere.py b/pvlib/atmosphere.py index a69541c93a..4878cc796d 100644 --- a/pvlib/atmosphere.py +++ b/pvlib/atmosphere.py @@ -380,7 +380,7 @@ def first_solar_spectral_correction(pw, airmass_absolute, module_type=None, The module used to calculate the spectral correction coefficients corresponds to the Mult-crystalline silicon Manufacturer 2 Model C from [3]_. Spectral Response (SR) of CIGS - and a-Si modules used to derive coefficients can be found in [4]_ + and a-Si modules used to derive coefficients can be found in [4]_ coefficients : None or array-like, default None allows for entry of user defined spectral correction diff --git a/pvlib/pvsystem.py b/pvlib/pvsystem.py index fdd7e2d9be..b19cce54f6 100644 --- a/pvlib/pvsystem.py +++ b/pvlib/pvsystem.py @@ -588,10 +588,9 @@ def __init__(self, pvsystem=None, location=None, **kwargs): Location.__init__(self, **new_kwargs) def __repr__(self): - attrs = [ - 'name', 'latitude', 'longitude', 'altitude', 'tz', 'surface_tilt', - 'surface_azimuth', 'module', 'inverter', 'albedo', 'racking_model' - ] + attrs = ['name', 'latitude', 'longitude', 'altitude', 'tz', + 'surface_tilt', 'surface_azimuth', 'module', 'inverter', + 'albedo', 'racking_model'] return ('LocalizedPVSystem: \n ' + '\n '.join( ('{}: {}'.format(attr, getattr(self, attr)) for attr in attrs))) @@ -1599,18 +1598,22 @@ def singlediode(photocurrent, saturation_current, resistance_series, photocurrent : numeric Light-generated current (photocurrent) in amperes under desired IV curve conditions. Often abbreviated ``I_L``. + 0 <= photocurrent saturation_current : numeric Diode saturation current in amperes under desired IV curve conditions. Often abbreviated ``I_0``. + 0 < saturation_current resistance_series : numeric Series resistance in ohms under desired IV curve conditions. Often abbreviated ``Rs``. + 0 <= resistance_series < numpy.inf resistance_shunt : numeric Shunt resistance in ohms under desired IV curve conditions. Often abbreviated ``Rsh``. + 0 < resistance_shunt <= numpy.inf nNsVth : numeric The product of three components. 1) The usual diode ideal factor @@ -1620,6 +1623,7 @@ def singlediode(photocurrent, saturation_current, resistance_series, ``k*temp_cell/q``, where k is Boltzmann's constant (J/K), temp_cell is the temperature of the p-n junction in Kelvin, and q is the charge of an electron (coulombs). + 0 < nNsVth ivcurve_pnts : None or int, default None Number of points in the desired IV curve. If None or 0, no @@ -1675,12 +1679,12 @@ def singlediode(photocurrent, saturation_current, resistance_series, calcparams_desoto ''' - # Find short circuit current using Lambert W - i_sc = i_from_v(resistance_shunt, resistance_series, nNsVth, 0.01, + # Compute short circuit current + i_sc = i_from_v(resistance_shunt, resistance_series, nNsVth, 0., saturation_current, photocurrent) - # Find open circuit voltage using Lambert W - v_oc = v_from_i(resistance_shunt, resistance_series, nNsVth, 0.0, + # Compute open circuit voltage + v_oc = v_from_i(resistance_shunt, resistance_series, nNsVth, 0., saturation_current, photocurrent) params = {'r_sh': resistance_shunt, @@ -1689,7 +1693,7 @@ def singlediode(photocurrent, saturation_current, resistance_series, 'i_0': saturation_current, 'i_l': photocurrent} - p_mp, v_mp = _golden_sect_DataFrame(params, 0, v_oc*1.14, _pwr_optfcn) + p_mp, v_mp = _golden_sect_DataFrame(params, 0., v_oc * 1.14, _pwr_optfcn) # Invert the Power-Current curve. Find the current where the inverted power # is minimized. This is i_mp. Start the optimization at v_oc/2 @@ -1697,11 +1701,11 @@ def singlediode(photocurrent, saturation_current, resistance_series, saturation_current, photocurrent) # Find Ix and Ixx using Lambert W - i_x = i_from_v(resistance_shunt, resistance_series, nNsVth, - 0.5*v_oc, saturation_current, photocurrent) + i_x = i_from_v(resistance_shunt, resistance_series, nNsVth, 0.5 * v_oc, + saturation_current, photocurrent) i_xx = i_from_v(resistance_shunt, resistance_series, nNsVth, - 0.5*(v_oc+v_mp), saturation_current, photocurrent) + 0.5 * (v_oc + v_mp), saturation_current, photocurrent) out = OrderedDict() out['i_sc'] = i_sc @@ -1716,9 +1720,10 @@ def singlediode(photocurrent, saturation_current, resistance_series, if ivcurve_pnts: ivcurve_v = (np.asarray(v_oc)[..., np.newaxis] * np.linspace(0, 1, ivcurve_pnts)) - ivcurve_i = i_from_v( - resistance_shunt, resistance_series, nNsVth, ivcurve_v.T, - saturation_current, photocurrent).T + + ivcurve_i = i_from_v(resistance_shunt, resistance_series, nNsVth, + ivcurve_v.T, saturation_current, photocurrent).T + out['v'] = ivcurve_v out['i'] = ivcurve_i @@ -1805,25 +1810,38 @@ def _pwr_optfcn(df, loc): Function to find power from ``i_from_v``. ''' - I = i_from_v(df['r_sh'], df['r_s'], df['nNsVth'], - df[loc], df['i_0'], df['i_l']) - return I*df[loc] + I = i_from_v(df['r_sh'], df['r_s'], df['nNsVth'], df[loc], df['i_0'], + df['i_l']) + + return I * df[loc] def v_from_i(resistance_shunt, resistance_series, nNsVth, current, saturation_current, photocurrent): ''' - Calculates voltage from current per Eq 3 Jain and Kapoor 2004 [1]. + Device voltage at the given device current for the single diode model. + + Uses the single diode model (SDM) as described in, e.g., + Jain and Kapoor 2004 [1]. + The solution is per Eq 3 of [1] except when resistance_shunt=numpy.inf, + in which case the explict solution for voltage is used. + Ideal device parameters are specified by resistance_shunt=np.inf and + resistance_series=0. + Inputs to this function can include scalars and pandas.Series, but it is + the caller's responsibility to ensure that the arguments are all float64 + and within the proper ranges. Parameters ---------- resistance_shunt : numeric Shunt resistance in ohms under desired IV curve conditions. Often abbreviated ``Rsh``. + 0 < resistance_shunt <= numpy.inf resistance_series : numeric Series resistance in ohms under desired IV curve conditions. Often abbreviated ``Rs``. + 0 <= resistance_series < numpy.inf nNsVth : numeric The product of three components. 1) The usual diode ideal factor @@ -1833,6 +1851,7 @@ def v_from_i(resistance_shunt, resistance_series, nNsVth, current, ``k*temp_cell/q``, where k is Boltzmann's constant (J/K), temp_cell is the temperature of the p-n junction in Kelvin, and q is the charge of an electron (coulombs). + 0 < nNsVth current : numeric The current in amperes under desired IV curve conditions. @@ -1840,14 +1859,16 @@ def v_from_i(resistance_shunt, resistance_series, nNsVth, current, saturation_current : numeric Diode saturation current in amperes under desired IV curve conditions. Often abbreviated ``I_0``. + 0 < saturation_current photocurrent : numeric Light-generated current (photocurrent) in amperes under desired IV curve conditions. Often abbreviated ``I_L``. + 0 <= photocurrent Returns ------- - current : np.array + current : np.ndarray or scalar References ---------- @@ -1860,51 +1881,105 @@ def v_from_i(resistance_shunt, resistance_series, nNsVth, current, except ImportError: raise ImportError('This function requires scipy') - Rsh = resistance_shunt - Rs = resistance_series - I0 = saturation_current - IL = photocurrent - I = current - - argW = I0 * Rsh / nNsVth * np.exp(Rsh * (-I + IL + I0) / nNsVth) - lambertwterm = lambertw(argW).real - - # Calculate using log(argW) in case argW is really big - logargW = (np.log(I0) + np.log(Rsh) - np.log(nNsVth) + - Rsh * (-I + IL + I0) / nNsVth) - - # Three iterations of Newton-Raphson method to solve - # w+log(w)=logargW. The initial guess is w=logargW. Where direct - # evaluation (above) results in NaN from overflow, 3 iterations - # of Newton's method gives approximately 8 digits of precision. - w = logargW - for i in range(0, 3): - w = w * (1 - np.log(w) + logargW) / (1 + w) - lambertwterm_log = w - - lambertwterm = np.where(np.isfinite(lambertwterm), lambertwterm, - lambertwterm_log) - - # Eqn. 3 in Jain and Kapoor, 2004 - V = -I*(Rs + Rsh) + IL*Rsh - nNsVth*lambertwterm + I0*Rsh - - return V + # Record if inputs were all scalar + output_is_scalar = all(map(np.isscalar, + [resistance_shunt, resistance_series, nNsVth, + current, saturation_current, photocurrent])) + + # This transforms Gsh=1/Rsh, including ideal Rsh=np.inf into Gsh=0., which + # is generally more numerically stable + conductance_shunt = 1./resistance_shunt + + # Ensure that we are working with read-only views of numpy arrays + # Turns Series into arrays so that we don't have to worry about + # multidimensional broadcasting failing + Gsh, Rs, a, I, I0, IL = \ + np.broadcast_arrays(conductance_shunt, resistance_series, nNsVth, + current, saturation_current, photocurrent) + + # Intitalize output V (I might not be float64) + V = np.full_like(I, np.nan, dtype=np.float64) + + # Determine indices where 0 < Gsh requires implicit model solution + idx_p = 0. < Gsh + + # Determine indices where 0 = Gsh allows explicit model solution + idx_z = 0. == Gsh + + # Explicit solutions where Gsh=0 + if np.any(idx_z): + V[idx_z] = a[idx_z]*np.log1p((IL[idx_z] - I[idx_z])/I0[idx_z]) - \ + I[idx_z]*Rs[idx_z] + + # Only compute using LambertW if there are cases with Gsh>0 + if np.any(idx_p): + # LambertW argument, cannot be float128, may overflow to np.inf + argW = I0[idx_p] / (Gsh[idx_p]*a[idx_p]) * \ + np.exp((-I[idx_p] + IL[idx_p] + I0[idx_p]) / + (Gsh[idx_p]*a[idx_p])) + + # lambertw typically returns complex value with zero imaginary part + # may overflow to np.inf + lambertwterm = lambertw(argW).real + + # Record indices where lambertw input overflowed output + idx_inf = np.logical_not(np.isfinite(lambertwterm)) + + # Only re-compute LambertW if it overflowed + if np.any(idx_inf): + # Calculate using log(argW) in case argW is really big + logargW = (np.log(I0[idx_p]) - np.log(Gsh[idx_p]) - + np.log(a[idx_p]) + + (-I[idx_p] + IL[idx_p] + I0[idx_p]) / + (Gsh[idx_p] * a[idx_p]))[idx_inf] + + # Three iterations of Newton-Raphson method to solve + # w+log(w)=logargW. The initial guess is w=logargW. Where direct + # evaluation (above) results in NaN from overflow, 3 iterations + # of Newton's method gives approximately 8 digits of precision. + w = logargW + for _ in range(0, 3): + w = w * (1. - np.log(w) + logargW) / (1. + w) + lambertwterm[idx_inf] = w + + # Eqn. 3 in Jain and Kapoor, 2004 + # V = -I*(Rs + Rsh) + IL*Rsh - a*lambertwterm + I0*Rsh + # Recast in terms of Gsh=1/Rsh for better numerical stability. + V[idx_p] = (IL[idx_p] + I0[idx_p] - I[idx_p])/Gsh[idx_p] - \ + I[idx_p]*Rs[idx_p] - a[idx_p]*lambertwterm + + if output_is_scalar: + return np.asscalar(V) + else: + return V def i_from_v(resistance_shunt, resistance_series, nNsVth, voltage, saturation_current, photocurrent): ''' - Calculates current from voltage per Eq 2 Jain and Kapoor 2004 [1]. + Device current at the given device voltage for the single diode model. + + Uses the single diode model (SDM) as described in, e.g., + Jain and Kapoor 2004 [1]. + The solution is per Eq 2 of [1] except when resistance_series=0, + in which case the explict solution for current is used. + Ideal device parameters are specified by resistance_shunt=np.inf and + resistance_series=0. + Inputs to this function can include scalars and pandas.Series, but it is + the caller's responsibility to ensure that the arguments are all float64 + and within the proper ranges. Parameters ---------- resistance_shunt : numeric Shunt resistance in ohms under desired IV curve conditions. Often abbreviated ``Rsh``. + 0 < resistance_shunt <= numpy.inf resistance_series : numeric Series resistance in ohms under desired IV curve conditions. Often abbreviated ``Rs``. + 0 <= resistance_series < numpy.inf nNsVth : numeric The product of three components. 1) The usual diode ideal factor @@ -1914,6 +1989,7 @@ def i_from_v(resistance_shunt, resistance_series, nNsVth, voltage, ``k*temp_cell/q``, where k is Boltzmann's constant (J/K), temp_cell is the temperature of the p-n junction in Kelvin, and q is the charge of an electron (coulombs). + 0 < nNsVth voltage : numeric The voltage in Volts under desired IV curve conditions. @@ -1921,14 +1997,16 @@ def i_from_v(resistance_shunt, resistance_series, nNsVth, voltage, saturation_current : numeric Diode saturation current in amperes under desired IV curve conditions. Often abbreviated ``I_0``. + 0 < saturation_current photocurrent : numeric Light-generated current (photocurrent) in amperes under desired IV curve conditions. Often abbreviated ``I_L``. + 0 <= photocurrent Returns ------- - current : np.array + current : np.ndarray or scalar References ---------- @@ -1941,23 +2019,58 @@ def i_from_v(resistance_shunt, resistance_series, nNsVth, voltage, except ImportError: raise ImportError('This function requires scipy') - # asarray turns Series into arrays so that we don't have to worry - # about multidimensional broadcasting failing - Rsh = np.asarray(resistance_shunt) - Rs = np.asarray(resistance_series) - I0 = np.asarray(saturation_current) - IL = np.asarray(photocurrent) - V = np.asarray(voltage) - - argW = (Rs*I0*Rsh * - np.exp(Rsh*(Rs*(IL+I0)+V) / (nNsVth*(Rs+Rsh))) / - (nNsVth*(Rs + Rsh))) - lambertwterm = lambertw(argW).real - - # Eqn. 4 in Jain and Kapoor, 2004 - I = -V/(Rs + Rsh) - (nNsVth/Rs)*lambertwterm + Rsh*(IL + I0)/(Rs + Rsh) - - return I + # Record if inputs were all scalar + output_is_scalar = all(map(np.isscalar, + [resistance_shunt, resistance_series, nNsVth, + voltage, saturation_current, photocurrent])) + + # This transforms Gsh=1/Rsh, including ideal Rsh=np.inf into Gsh=0., which + # is generally more numerically stable + conductance_shunt = 1./resistance_shunt + + # Ensure that we are working with read-only views of numpy arrays + # Turns Series into arrays so that we don't have to worry about + # multidimensional broadcasting failing + Gsh, Rs, a, V, I0, IL = \ + np.broadcast_arrays(conductance_shunt, resistance_series, nNsVth, + voltage, saturation_current, photocurrent) + + # Intitalize output I (V might not be float64) + I = np.full_like(V, np.nan, dtype=np.float64) + + # Determine indices where 0 < Rs requires implicit model solution + idx_p = 0. < Rs + + # Determine indices where 0 = Rs allows explicit model solution + idx_z = 0. == Rs + + # Explicit solutions where Rs=0 + if np.any(idx_z): + I[idx_z] = IL[idx_z] - I0[idx_z]*np.expm1(V[idx_z]/a[idx_z]) - \ + Gsh[idx_z]*V[idx_z] + + # Only compute using LambertW if there are cases with Rs>0 + # Does NOT handle possibility of overflow, github issue 298 + if np.any(idx_p): + # LambertW argument, cannot be float128, may overflow to np.inf + argW = Rs[idx_p]*I0[idx_p]/(a[idx_p]*(Rs[idx_p]*Gsh[idx_p] + 1.)) * \ + np.exp((Rs[idx_p]*(IL[idx_p] + I0[idx_p]) + V[idx_p]) / + (a[idx_p]*(Rs[idx_p]*Gsh[idx_p] + 1.))) + + # lambertw typically returns complex value with zero imaginary part + # may overflow to np.inf + lambertwterm = lambertw(argW).real + + # Eqn. 2 in Jain and Kapoor, 2004 + # I = -V/(Rs + Rsh) - (a/Rs)*lambertwterm + Rsh*(IL + I0)/(Rs + Rsh) + # Recast in terms of Gsh=1/Rsh for better numerical stability. + I[idx_p] = (IL[idx_p] + I0[idx_p] - V[idx_p]*Gsh[idx_p]) / \ + (Rs[idx_p]*Gsh[idx_p] + 1.) - (a[idx_p]/Rs[idx_p])*lambertwterm + + if output_is_scalar: + return np.asscalar(I) + else: + return I def snlinverter(v_dc, p_dc, inverter): diff --git a/pvlib/test/test_irradiance.py b/pvlib/test/test_irradiance.py index ba4a1ba74c..4a6fe224db 100755 --- a/pvlib/test/test_irradiance.py +++ b/pvlib/test/test_irradiance.py @@ -425,14 +425,15 @@ def test_dni(): @pytest.mark.parametrize( - 'surface_tilt,surface_azimuth,solar_zenith,solar_azimuth,aoi_expected,aoi_proj_expected', [ - (0, 0, 0, 0, 0, 1), - (30, 180, 30, 180, 0, 1), - (30, 180, 150, 0, 180, -1), - (90, 0, 30, 60, 75.5224878, 0.25), - (90, 0, 30, 170, 119.4987042, -0.4924038)]) + 'surface_tilt,surface_azimuth,solar_zenith,' + + 'solar_azimuth,aoi_expected,aoi_proj_expected', + [(0, 0, 0, 0, 0, 1), + (30, 180, 30, 180, 0, 1), + (30, 180, 150, 0, 180, -1), + (90, 0, 30, 60, 75.5224878, 0.25), + (90, 0, 30, 170, 119.4987042, -0.4924038)]) def test_aoi_and_aoi_projection(surface_tilt, surface_azimuth, solar_zenith, - solar_azimuth, aoi_expected, + solar_azimuth, aoi_expected, aoi_proj_expected): aoi = irradiance.aoi(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth) diff --git a/pvlib/test/test_modelchain.py b/pvlib/test/test_modelchain.py index aa5726f2cc..9c13e4b2f8 100644 --- a/pvlib/test/test_modelchain.py +++ b/pvlib/test/test_modelchain.py @@ -444,6 +444,7 @@ def test_ModelChain___repr__(system, location, strategy, strategy_str): assert mc.__repr__() == expected + @requires_scipy def test_weather_irradiance_input(system, location): """Test will raise a warning and should be removed in future versions.""" diff --git a/pvlib/test/test_pvsystem.py b/pvlib/test/test_pvsystem.py index 9283536846..c7e7345ca4 100644 --- a/pvlib/test/test_pvsystem.py +++ b/pvlib/test/test_pvsystem.py @@ -353,37 +353,245 @@ def test_PVSystem_calcparams_desoto(cec_module_params): assert_allclose(nNsVth, 0.473) -@requires_scipy -def test_v_from_i(): - output = pvsystem.v_from_i(20, .1, .5, 3, 6e-7, 7) - assert_allclose(7.5049875193450521, output, atol=1e-5) - - -@requires_scipy -def test_v_from_i_big(): - output = pvsystem.v_from_i(500, 10, 4.06, 0, 6e-10, 1.2) - assert_allclose(86.320000493521079, output, atol=1e-5) +@pytest.fixture(params=[ + { # Can handle all python scalar inputs + 'Rsh': 20., + 'Rs': 0.1, + 'nNsVth': 0.5, + 'I': 3., + 'I0': 6.e-7, + 'IL': 7., + 'V_expected': 7.5049875193450521 + }, + { # Can handle all rank-0 array inputs + 'Rsh': np.array(20.), + 'Rs': np.array(0.1), + 'nNsVth': np.array(0.5), + 'I': np.array(3.), + 'I0': np.array(6.e-7), + 'IL': np.array(7.), + 'V_expected': np.array(7.5049875193450521) + }, + { # Can handle all rank-1 singleton array inputs + 'Rsh': np.array([20.]), + 'Rs': np.array([0.1]), + 'nNsVth': np.array([0.5]), + 'I': np.array([3.]), + 'I0': np.array([6.e-7]), + 'IL': np.array([7.]), + 'V_expected': np.array([7.5049875193450521]) + }, + { # Can handle all rank-1 non-singleton array inputs with infinite shunt + # resistance, Rsh=inf gives V=Voc=nNsVth*(np.log(IL + I0) - np.log(I0) + # at I=0 + 'Rsh': np.array([np.inf, 20.]), + 'Rs': np.array([0.1, 0.1]), + 'nNsVth': np.array([0.5, 0.5]), + 'I': np.array([0., 3.]), + 'I0': np.array([6.e-7, 6.e-7]), + 'IL': np.array([7., 7.]), + 'V_expected': np.array([0.5*(np.log(7. + 6.e-7) - np.log(6.e-7)), + 7.5049875193450521]) + }, + { # Can handle mixed inputs with a rank-2 array with infinite shunt + # resistance, Rsh=inf gives V=Voc=nNsVth*(np.log(IL + I0) - np.log(I0) + # at I=0 + 'Rsh': np.array([[np.inf, np.inf], [np.inf, np.inf]]), + 'Rs': np.array([0.1]), + 'nNsVth': np.array(0.5), + 'I': 0., + 'I0': np.array([6.e-7]), + 'IL': np.array([7.]), + 'V_expected': 0.5*(np.log(7. + 6.e-7) - np.log(6.e-7))*np.ones((2, 2)) + }, + { # Can handle ideal series and shunt, Rsh=inf and Rs=0 give + # V = nNsVth*(np.log(IL - I + I0) - np.log(I0)) + 'Rsh': np.inf, + 'Rs': 0., + 'nNsVth': 0.5, + 'I': np.array([7., 7./2., 0.]), + 'I0': 6.e-7, + 'IL': 7., + 'V_expected': np.array([0., 0.5*(np.log(7. - 7./2. + 6.e-7) - + np.log(6.e-7)), 0.5*(np.log(7. + 6.e-7) - + np.log(6.e-7))]) + }, + { # Can handle only ideal series resistance, no closed form solution + 'Rsh': 20., + 'Rs': 0., + 'nNsVth': 0.5, + 'I': 3., + 'I0': 6.e-7, + 'IL': 7., + 'V_expected': 7.804987519345062 + }, + { # Can handle all python scalar inputs with big LambertW arg + 'Rsh': 500., + 'Rs': 10., + 'nNsVth': 4.06, + 'I': 0., + 'I0': 6.e-10, + 'IL': 1.2, + 'V_expected': 86.320000493521079 + }, + { # Can handle all python scalar inputs with bigger LambertW arg + # 1000 W/m^2 on a Canadian Solar 220M with 20 C ambient temp + # github issue 225 (this appears to be from PR 226 not issue 225) + 'Rsh': 190., + 'Rs': 1.065, + 'nNsVth': 2.89, + 'I': 0., + 'I0': 7.05196029e-08, + 'IL': 10.491262, + 'V_expected': 54.303958833791455 + }, + { # Can handle all python scalar inputs with bigger LambertW arg + # 1000 W/m^2 on a Canadian Solar 220M with 20 C ambient temp + # github issue 225 + 'Rsh': 381.68, + 'Rs': 1.065, + 'nNsVth': 2.681527737715915, + 'I': 0., + 'I0': 1.8739027472625636e-09, + 'IL': 5.1366949999999996, + 'V_expected': 58.19323124611128 + }, + { # Verify mixed solution type indexing logic + 'Rsh': np.array([np.inf, 190., 381.68]), + 'Rs': 1.065, + 'nNsVth': np.array([2.89, 2.89, 2.681527737715915]), + 'I': 0., + 'I0': np.array([7.05196029e-08, 7.05196029e-08, 1.8739027472625636e-09]), + 'IL': np.array([10.491262, 10.491262, 5.1366949999999996]), + 'V_expected': np.array([2.89*np.log1p(10.491262/7.05196029e-08), + 54.303958833791455, 58.19323124611128]) + }]) +def fixture_v_from_i(request): + return request.param @requires_scipy -def test_v_from_i_bigger(): - # 1000 W/m^2 on a Canadian Solar 220M with 20 C ambient temp - # github issue 225 - output = pvsystem.v_from_i(190, 1.065, 2.89, 0, 7.05196029e-08, 10.491262) - assert_allclose(54.303958833791455, output, atol=1e-5) +def test_v_from_i(fixture_v_from_i): + # Solution set loaded from fixture + Rsh = fixture_v_from_i['Rsh'] + Rs = fixture_v_from_i['Rs'] + nNsVth = fixture_v_from_i['nNsVth'] + I = fixture_v_from_i['I'] + I0 = fixture_v_from_i['I0'] + IL = fixture_v_from_i['IL'] + V_expected = fixture_v_from_i['V_expected'] + + # Convergence criteria + atol = 1.e-11 + + V = pvsystem.v_from_i(Rsh, Rs, nNsVth, I, I0, IL) + assert(isinstance(V, type(V_expected))) + if isinstance(V, type(np.ndarray)): + assert(isinstance(V.dtype, type(V_expected.dtype))) + assert(V.shape == V_expected.shape) + assert_allclose(V, V_expected, atol=atol) + + +@pytest.fixture(params=[ + { # Can handle all python scalar inputs + 'Rsh': 20., + 'Rs': 0.1, + 'nNsVth': 0.5, + 'V': 40., + 'I0': 6.e-7, + 'IL': 7., + 'I_expected': -299.746389916 + }, + { # Can handle all rank-0 array inputs + 'Rsh': np.array(20.), + 'Rs': np.array(0.1), + 'nNsVth': np.array(0.5), + 'V': np.array(40.), + 'I0': np.array(6.e-7), + 'IL': np.array(7.), + 'I_expected': np.array(-299.746389916) + }, + { # Can handle all rank-1 singleton array inputs + 'Rsh': np.array([20.]), + 'Rs': np.array([0.1]), + 'nNsVth': np.array([0.5]), + 'V': np.array([40.]), + 'I0': np.array([6.e-7]), + 'IL': np.array([7.]), + 'I_expected': np.array([-299.746389916]) + }, + { # Can handle all rank-1 non-singleton array inputs with a zero + # series resistance, Rs=0 gives I=IL=Isc at V=0 + 'Rsh': np.array([20., 20.]), + 'Rs': np.array([0., 0.1]), + 'nNsVth': np.array([0.5, 0.5]), + 'V': np.array([0., 40.]), + 'I0': np.array([6.e-7, 6.e-7]), + 'IL': np.array([7., 7.]), + 'I_expected': np.array([7., -299.746389916]) + }, + { # Can handle mixed inputs with a rank-2 array with zero series + # resistance, Rs=0 gives I=IL=Isc at V=0 + 'Rsh': np.array([20.]), + 'Rs': np.array([[0., 0.], [0., 0.]]), + 'nNsVth': np.array(0.5), + 'V': 0., + 'I0': np.array([6.e-7]), + 'IL': np.array([7.]), + 'I_expected': np.array([[7., 7.], [7., 7.]]) + }, + { # Can handle ideal series and shunt, Rsh=inf and Rs=0 give + # V_oc = nNsVth*(np.log(IL + I0) - np.log(I0)) + 'Rsh': np.inf, + 'Rs': 0., + 'nNsVth': 0.5, + 'V': np.array([0., 0.5*(np.log(7. + 6.e-7) - np.log(6.e-7))/2., + 0.5*(np.log(7. + 6.e-7) - np.log(6.e-7))]), + 'I0': 6.e-7, + 'IL': 7., + 'I_expected': np.array([7., 7. - 6.e-7*np.expm1((np.log(7. + 6.e-7) - + np.log(6.e-7))/2.), 0.]) + }, + { # Can handle only ideal shunt resistance, no closed form solution + 'Rsh': np.inf, + 'Rs': 0.1, + 'nNsVth': 0.5, + 'V': 40., + 'I0': 6.e-7, + 'IL': 7., + 'I_expected': -299.7383436645412 + }]) +def fixture_i_from_v(request): + return request.param @requires_scipy -def test_i_from_v(): - output = pvsystem.i_from_v(20, .1, .5, 40, 6e-7, 7) - assert_allclose(-299.746389916, output, atol=1e-5) +def test_i_from_v(fixture_i_from_v): + # Solution set loaded from fixture + Rsh = fixture_i_from_v['Rsh'] + Rs = fixture_i_from_v['Rs'] + nNsVth = fixture_i_from_v['nNsVth'] + V = fixture_i_from_v['V'] + I0 = fixture_i_from_v['I0'] + IL = fixture_i_from_v['IL'] + I_expected = fixture_i_from_v['I_expected'] + + # Convergence criteria + atol = 1.e-11 + + I = pvsystem.i_from_v(Rsh, Rs, nNsVth, V, I0, IL) + assert(isinstance(I, type(I_expected))) + if isinstance(I, type(np.ndarray)): + assert(isinstance(I.dtype, type(I_expected.dtype))) + assert(I.shape == I_expected.shape) + assert_allclose(I, I_expected, atol=atol) @requires_scipy def test_PVSystem_i_from_v(): system = pvsystem.PVSystem() output = system.i_from_v(20, .1, .5, 40, 6e-7, 7) - assert_allclose(-299.746389916, output, atol=1e-5) + assert_allclose(output, -299.746389916, atol=1e-5) @requires_scipy @@ -440,7 +648,7 @@ def test_singlediode_floats(sam_data): if k in ['i', 'v']: assert v is None else: - assert_allclose(expected[k], v, atol=3) + assert_allclose(v, expected[k], atol=3) @requires_scipy @@ -453,11 +661,11 @@ def test_singlediode_floats_ivcurve(): 'i_x': 6.7556075876880621, 'i_sc': 6.9646747613963198, 'v_mp': 6.221535886625464, - 'i': np.array([6.965172e+00, 6.755882e+00, 2.575717e-14]), - 'v': np.array([0. , 4.05315, 8.1063])} + 'i': np.array([6.965172e+00, 6.755882e+00, 2.575717e-14]), + 'v': np.array([0., 4.05315, 8.1063])} assert isinstance(out, dict) for k, v in out.items(): - assert_allclose(expected[k], v, atol=3) + assert_allclose(v, expected[k], atol=3) @requires_scipy @@ -465,33 +673,31 @@ def test_singlediode_series_ivcurve(cec_module_params): times = pd.DatetimeIndex(start='2015-06-01', periods=3, freq='6H') poa_data = pd.Series([0, 400, 800], index=times) IL, I0, Rs, Rsh, nNsVth = pvsystem.calcparams_desoto( - poa_data, - temp_cell=25, - alpha_isc=cec_module_params['alpha_sc'], - module_parameters=cec_module_params, - EgRef=1.121, - dEgdT=-0.0002677) + poa_data, temp_cell=25, + alpha_isc=cec_module_params['alpha_sc'], + module_parameters=cec_module_params, + EgRef=1.121, dEgdT=-0.0002677) out = pvsystem.singlediode(IL, I0, Rs, Rsh, nNsVth, ivcurve_pnts=3) - expected = OrderedDict([('i_sc', array([ nan, 3.01054475, 6.00675648])), - ('v_oc', array([ nan, 9.96886962, 10.29530483])), - ('i_mp', array([ nan, 2.65191983, 5.28594672])), - ('v_mp', array([ nan, 8.33392491, 8.4159707 ])), - ('p_mp', array([ nan, 22.10090078, 44.48637274])), - ('i_x', array([ nan, 2.88414114, 5.74622046])), - ('i_xx', array([ nan, 2.04340914, 3.90007956])), - ('v', - array([[ nan, nan, nan], - [ 0. , 4.98443481, 9.96886962], - [ 0. , 5.14765242, 10.29530483]])), - ('i', - array([[ nan, nan, nan], - [ 3.01079860e+00, 2.88414114e+00, 3.10862447e-14], - [ 6.00726296e+00, 5.74622046e+00, 0.00000000e+00]]))]) + expected = OrderedDict([('i_sc', array([0., 3.01054475, 6.00675648])), + ('v_oc', array([0., 9.96886962, 10.29530483])), + ('i_mp', array([0., 2.65191983, 5.28594672])), + ('v_mp', array([0., 8.33392491, 8.4159707])), + ('p_mp', array([0., 22.10090078, 44.48637274])), + ('i_x', array([0., 2.88414114, 5.74622046])), + ('i_xx', array([0., 2.04340914, 3.90007956])), + ('v', array([[0., 0., 0.], + [0., 4.98443481, 9.96886962], + [0., 5.14765242, 10.29530483]])), + ('i', array([[0., 0., 0.], + [3.01079860e+00, 2.88414114e+00, + 3.10862447e-14], + [6.00726296e+00, 5.74622046e+00, + 0.00000000e+00]]))]) for k, v in out.items(): - assert_allclose(expected[k], v, atol=1e-2) + assert_allclose(v, expected[k], atol=1e-2) def test_scale_voltage_current_power(sam_data): @@ -521,17 +727,17 @@ def test_PVSystem_scale_voltage_current_power(): def test_sapm_celltemp(): default = pvsystem.sapm_celltemp(900, 5, 20) - assert_allclose(43.509, default['temp_cell'], 3) - assert_allclose(40.809, default['temp_module'], 3) + assert_allclose(default['temp_cell'], 43.509, 3) + assert_allclose(default['temp_module'], 40.809, 3) assert_frame_equal(default, pvsystem.sapm_celltemp(900, 5, 20, [-3.47, -.0594, 3])) def test_sapm_celltemp_dict_like(): default = pvsystem.sapm_celltemp(900, 5, 20) - assert_allclose(43.509, default['temp_cell'], 3) - assert_allclose(40.809, default['temp_module'], 3) - model = {'a':-3.47, 'b':-.0594, 'deltaT':3} + assert_allclose(default['temp_cell'], 43.509, 3) + assert_allclose(default['temp_module'], 40.809, 3) + model = {'a': -3.47, 'b': -.0594, 'deltaT': 3} assert_frame_equal(default, pvsystem.sapm_celltemp(900, 5, 20, model)) model = pd.Series(model) assert_frame_equal(default, pvsystem.sapm_celltemp(900, 5, 20, model)) @@ -759,7 +965,7 @@ def test_LocalizedPVSystem___repr__(): def test_pvwatts_dc_scalars(): expected = 88.65 out = pvsystem.pvwatts_dc(900, 30, 100, -0.003) - assert_allclose(expected, out) + assert_allclose(out, expected) @needs_numpy_1_10 @@ -767,11 +973,11 @@ def test_pvwatts_dc_arrays(): irrad_trans = np.array([np.nan, 900, 900]) temp_cell = np.array([30, np.nan, 30]) irrad_trans, temp_cell = np.meshgrid(irrad_trans, temp_cell) - expected = np.array([[ nan, 88.65, 88.65], - [ nan, nan, nan], - [ nan, 88.65, 88.65]]) + expected = np.array([[nan, 88.65, 88.65], + [nan, nan, nan], + [nan, 88.65, 88.65]]) out = pvsystem.pvwatts_dc(irrad_trans, temp_cell, 100, -0.003) - assert_allclose(expected, out, equal_nan=True) + assert_allclose(out, expected, equal_nan=True) def test_pvwatts_dc_series(): @@ -785,18 +991,18 @@ def test_pvwatts_dc_series(): def test_pvwatts_ac_scalars(): expected = 85.58556604752516 out = pvsystem.pvwatts_ac(90, 100, 0.95) - assert_allclose(expected, out) + assert_allclose(out, expected) @needs_numpy_1_10 def test_pvwatts_ac_arrays(): pdc = np.array([[np.nan], [50], [100]]) pdc0 = 100 - expected = np.array([[ nan], - [ 47.60843624], - [ 95. ]]) + expected = np.array([[nan], + [47.60843624], + [95.]]) out = pvsystem.pvwatts_ac(pdc, pdc0, 0.95) - assert_allclose(expected, out, equal_nan=True) + assert_allclose(out, expected, equal_nan=True) def test_pvwatts_ac_series(): @@ -810,7 +1016,7 @@ def test_pvwatts_ac_series(): def test_pvwatts_losses_default(): expected = 14.075660688264469 out = pvsystem.pvwatts_losses() - assert_allclose(expected, out) + assert_allclose(out, expected) @needs_numpy_1_10 @@ -818,7 +1024,7 @@ def test_pvwatts_losses_arrays(): expected = np.array([nan, 14.934904]) age = np.array([nan, 1]) out = pvsystem.pvwatts_losses(age=age) - assert_allclose(expected, out) + assert_allclose(out, expected) def test_pvwatts_losses_series():