pvlib-python
pvlib-python copied to clipboard
IAM that supports AR coating like Fresnel
Problem
Currently pvlib supports the DeSoto physical model (similar to normal glass), ASHRAE, Martin & Ruiz, and SAPM polynomial, but it doesn't have a pure Fresnel model that allows additional interfaces like an AR coating.
- DeSoto physical model is most similar to the Fresnel for normal glass but only has one interface, so is limited to IAM curves below it only, while an AR coating would have a greater ρ
- Martin & Ruiz could be used to approximate an AR coated glass if the correct
a_r
were known. The default ofa_r=0.16
is slightly above the normal glass Fresnel IAM, but ana_r=0.14
seems to match an AR coating with index of refraction of 1.2 most closely.
Proposal
a new method in pvl.iam.fresnel_ar(aoi, n_ar=1.2, n_air=1.0, n_glass=1.56)
that implements the Fresnel equation
Alternative
Suggest readers to use Martin & Ruiz with a_r=0.14
instead of default.
additional content
PVsyst has switched to Fresnel equations. We can duplicate their methods ignoring additional reflections and the encapsulant layer:
import numpy as np
import matplotlib.pyplot as plt
plt.ion()
# constants
n_glass = 1.56
n_air = 1.0
theta_inc = np.linspace(0, 88, 100)
def snell(theta_1, n1, n2):
"""Snell's equation"""
sintheta_2 = n1/n2 * np.sin(np.radians(theta_1))
return sintheta_2, np.degrees(np.arcsin(sintheta_2))
def refl_s(theta_1, theta_2, n1, n2):
"""Fresnel's equation"""
n1_costheta_1 = n1*np.cos(np.radians(theta_1))
n2_costheta_2 = n2*np.cos(np.radians(theta_2))
return np.abs((n1_costheta_1 - n2_costheta_2)/(n1_costheta_1 + n2_costheta_2))**2
def refl_p(theta_1, theta_2, n1, n2):
"""Fresnel's equation"""
n1_costheta_2 = n1*np.cos(np.radians(theta_2))
n2_costheta_1 = n2*np.cos(np.radians(theta_1))
return np.abs((n1_costheta_2 - n2_costheta_1)/(n1_costheta_2 + n2_costheta_1))**2
def refl_eff(rs, rp):
"""effective reflectivity"""
return (rs+rp)/2
def trans(refl):
"""transmissivity"""
return 1-refl
def refl0(n1, n2):
"""reflectivity at normal incidence"""
return np.abs((n1-n2)/(n1+n2))**2
def fresnel(theta_inc, n1=n_air, n2=n_glass):
"""calculate IAM using Fresnel's Law"""
_, theta_tr = snell(theta_inc, n1, n2)
rs = refl_s(theta_inc, theta_tr, n1, n2)
rp = refl_p(theta_inc, theta_tr, n1, n2)
reff = refl_eff(rs, rp)
r0 = refl0(n1, n2)
return trans(reff)/trans(r0)
def ashrae(theta_inc, b0=0.05):
"""ASHRAE equation"""
return 1 - b0*(1/np.cos(np.radians(theta_inc)) - 1)
def fresnel_ar(theta_inc, n_ar, n1=n_air, n2=n_glass):
"""calculate IAM using Fresnel's law with AR"""
# use fresnel() for n2=n_ar
_, theta_ar = snell(theta_inc, n1, n_ar)
rs_ar1 = refl_s(theta_inc, theta_ar, n1, n_ar)
rp_ar1 = refl_p(theta_inc, theta_ar, n1, n_ar)
r0_ar1 = refl0(n1, n_ar)
# repeat with fresnel() with n1=n_ar
_, theta_tr = snell(theta_ar, n_ar, n2)
rs = refl_s(theta_ar, theta_tr, n_ar, n2)
rp = refl_p(theta_ar, theta_tr, n_ar, n2)
# note that combined reflectivity is product of transmissivity!
# so... rho12 = 1 - (1-rho1)(1-rho2)
reff = refl_eff(1-(1-rs_ar1)*(1-rs), 1-(1-rp_ar1)*(1-rp))
r0 = 1-(1-refl0(n_ar, n2))*(1-r0_ar1)
return trans(reff)/trans(r0)
# plot Fresnel for normal glass and ASHRAE
plt.plot(theta_inc, fresnel(theta_inc))
plt.plot(theta_inc, ashrae(theta_inc))
# calculate IAM for AR with n=1.1 and plot
iam_ar11 = fresnel_ar(theta_inc, n_ar=1.1)
plt.plot(theta_inc, iam_ar11)
# repeat for AR with n=1.2
iam_ar12 = fresnel_ar(theta_inc, n_ar=1.2)
plt.plot(theta_inc, iam_ar12)
# make plot pretty
plt.legend(['Fresnel, normal glass', 'ASHRAE, $b_0=0.05$', 'Fresnel $n_{AR}=1.1$', 'Fresnel $n_{AR}=1.2$'])
plt.title("IAM correction, Fresnel vs. ASHRAE, using basic eqn's")
plt.ylabel('IAM')
plt.xlabel(r'incidence angle $\theta_{inc} [\degree]$')
plt.grid()
plt.ylim([0.55,1.05])
+1. This reference might be relevant.
I seem to recall from somewhere that PVsyst actually interpolates from a fixed set of pre-calculated values when simulating.
PVsyst allows a user specified custom IAM v AOI lookup table in the module PAN file, but that presupposes there exist qualified IAM measurements either from a lab or the manufacturer. Otherwise they use Fresnel as of v6.67. See https://www.pvsyst.com/help/iam_loss.htm
Yes, what I meant is that they use the Fresnel equations to populate the table for interpolation. At least this is my recollection.
In #1616 I submitted a proposal to extend the ìam.physical()
function with an extra optional argument n_ar
to support AR coating.
Main differences with the code proposed above by @mikofski :
- instead of adding new functions
fresnel
andfresnel_ar
, the existing functionphysical
(which is equivalent withfresnel
) is extended with support for AR coating through the optional argumentn_ar
. - internal transmissions are no longer ignored but taken into account.
Main differences with the previous code of iam.physical()
(besides support for AR coating):
- the Fresnel equations for the reflectance for parallel and perpendicular polarized light has been adapted from the tan formulation (as used by Desoto) to the much more common cos formulation (as used everywhere else). As far as I can tell they are mathematically fully equivalent (in combination with Snell's law), but the cos formulation doesn't suffer from 0/0 division for aoi equal to zero. That's why no special handling is needed anymore for angles near zero.
- the result for 90° is no longer exactly 0, but (on my laptop) 3.47262923e-16, which I think is fully fine. But for that I had to include
a_tol
in the existingtest_physical
unit test. - (probably not relevant) though aoi is normally expected to be between 0° and 180°, there were apparently also tests for negative aoi. The new code is now also consistent between positive and negative angles for angles outside the range -360°, 360° (though I didn't include tests for it as I doubt whether we should really support that).
Main differences with PVSyst:
- I don't know the reason but there is still a (very minor) difference between the values in PVSyst compared to the values returned by this function.
- the result for 90° is no longer exactly 0, but (on my laptop) 3.47262923e-16, which I think is fully fine. But for that I had to include
a_tol
in the existingtest_physical
unit test.
For fun you could try using scipy.special.cosdg()
- I don't know the reason but there is still a (very minor) difference between the values in PVSyst compared to the values returned by this function.
I don't think they do the internal reflections.
- the result for 90° is no longer exactly 0, but (on my laptop) 3.47262923e-16, which I think is fully fine. But for that I had to include
a_tol
in the existingtest_physical
unit test.For fun you could try using
scipy.special.cosdg()
I didn't know about scipy.special.cosdg()
. Using that function one gets indeed exactly 0 for aoi of 90°.
I still don't think it matters as:
- 3.47262923e-16 is really very close to zero
- for aoi > 90° we already get exactly 0
- timing of sunset, sunrise has its uncertainties and depends on atmospheric conditions
- the sun moves by 1° every 4 minutes
- the diameter of the sun is 0.53°.
Nevertheless, maybe pvlib should in general use these scipy.special
functions instead of the pvlib.tools
functions.
- I don't know the reason but there is still a (very minor) difference between the values in PVSyst compared to the values returned by this function.
I don't think they do the internal reflections.
Results are actually closer to PVSyst when taking internal reflections into account compared to when not taking them into account. Compare e.g. the results for n_ar = 1.29
(and L=0
):
aoi | pvsyst | with reflections | without reflections |
---|---|---|---|
0 | 1 | 1 | 1 |
30 | 0.999 | 0.999 | 0.999 |
50 | 0.987 | 0.987 | 0.987 |
60 | 0.962 | 0.963 | 0.962 |
70 | 0.892 | 0.892 | 0.890 |
75 | 0.816 | 0.814 | 0.812 |
80 | 0.681 | 0.679 | 0.675 |
85 | 0.440 | 0.437 | 0.434 |
90 | 0 | 0 | 0 |
So, I'd guess also PVSyst takes internal reflections into account.
Are all the relevant parameter values documented in PVsyst? I can't find them on their web site.
I’m -1 on scipy.special.cosdg
see previous thread.
great to see closer agreement with PVsyst, not that it’s the reference, but still it’s de facto. Their online help says they consider reflections, but in my pseudo-code I ignored it. What is L parameter?
great to see closer agreement with PVsyst, not that it’s the reference, but still it’s de facto. Their online help says they consider reflections, but in my pseudo-code I ignored it. What is L parameter?
L is the thickness of the glass. It's an optional parameter of iam.physical
, which is by default 2 mm (though most PV panels are actually 3.2 mm thick). It has a minor impact on the IAM values. You considered it zero as well in your calculation above. I don't know what PVSyst considers but taking it into account reduces the IAM values and thus causes an increase of the deviations with regard to PVSyst.
Are all the relevant parameter values documented in PVsyst? I can't find them on their web site.
The refractive index for the glass and AR coating can be selected in the tool and its default values (n_glass=1.526, n_ar=1.29) for the Incidence Angle mode 'Fresnel, Glass with AF' can be seen in some threads on the PVsyst forum, e.g. here and here. I couldn't find any values for the other two parameters (glazing thickness L
and glazing extinction coefficient K
), but that's possibly because PVsyst ignores the absorption in the IAM calculation. That's why I assumed L=0
for the comparison table (and the glazing extinction coefficient doesn't matter any more when the glazing thickness is zero). Moreover, assuming a non-zero glazing thickness would only lower the values, and we need higher values to get closer to the PVsyst values.
Well, from my perspective if it is clear you implement these formulas correctly for a two-layer window then it doesn't need to match PVsyst or claim that it is doing the same as PVsyst.
PVWatts v5 and SAM both have the option of a second snell & fresnel iteration for AR coatings, although I don't think either of them consider repeated internal reflections. I point these out as potentially useful inclusions in the references and possible influences on the name of the new function.
L is the thickness of the glass. It's an optional parameter of
iam.physical
, which is by default 2 mm (though most PV panels are actually 3.2 mm thick).
This is one default that has always bothered me out of principle. For the incorrect value there is an academic paper that can be referenced, thus giving someone else the responsibility, but for the correct value we haven't got one!
L is the thickness of the glass. It's an optional parameter of
iam.physical
, which is by default 2 mm (though most PV panels are actually 3.2 mm thick).This is one default that has always bothered me out of principle. For the incorrect value there is an academic paper that can be referenced, thus giving someone else the responsibility, but for the correct value we haven't got one!
I'd be OK changing it, assuming explanatory comments are added.
The only reason to stick with L=0.002 m is to stay consistent with the reference. The basis for that value: it was what W. De Soto used in his Master's thesis (same title as the paper), with no further discussion provided. So we shouldn't regard 2mm as the result of any concerted research.
L is the thickness of the glass. It's an optional parameter of
iam.physical
, which is by default 2 mm (though most PV panels are actually 3.2 mm thick).This is one default that has always bothered me out of principle. For the incorrect value there is an academic paper that can be referenced, thus giving someone else the responsibility, but for the correct value we haven't got one!
FYI, the yearly 'International Technology Roadmap for Photovoltaics (ITRPV)' (donwloadable from https://www.vdma.org/international-technology-roadmap-photovoltaic) has a figure about the "world market share of front glass thickness in modules". Below you find the relevant figure in the latest report from November 2022.
Good reference to have. I suspect that the increasing fraction of thinner front glass is related to an increasing fraction of glass-glass modules.