diff --git a/shared/lib_irradproc.cpp b/shared/lib_irradproc.cpp index 35db1295a..53ce8d0c5 100644 --- a/shared/lib_irradproc.cpp +++ b/shared/lib_irradproc.cpp @@ -1692,7 +1692,7 @@ perez(double, double dn, double df, double alb, double inc, double tilt, double } } -void ineichen(double clearsky_results[3], double apparent_zenith, double absolute_airmass, double linke_turbidity = 1.0, double altitude = 0.0, double dni_extra = 1364.0, bool perez_enhancement = false) { +void ineichen(double clearsky_results[3], double apparent_zenith, int month, int day, double pressure = 101325.0, double linke_turbidity = 1.0, double altitude = 0.0, double dni_extra = 1364.0, bool perez_enhancement = false) { double cos_zenith = Max(cosd(apparent_zenith), 0); double tl = linke_turbidity; @@ -1701,15 +1701,25 @@ void ineichen(double clearsky_results[3], double apparent_zenith, double absolut double cg1 = 5.09e-5 * altitude + 0.868; double cg2 = 3.92e-5 * altitude + 0.0387; - double ghi = exp(-cg2 * absolute_airmass * (fh1 + fh2 * (tl - 1))); - if (perez_enhancement) ghi *= exp(0.01 * pow(absolute_airmass, 1.8)); + //double am = Min(15.25, 1.0 / (cosd(apparent_zenith) + 0.15 * (pow(93.9 - apparent_zenith, -1.253)))); // air mass + double am = Min(15.25, 1.0 / (cosd(apparent_zenith) + 0.50572 * (pow(6.07995 + (90 - apparent_zenith), -1.6364)))); // air mass kastenyoung1989 pvlib - ghi = cg1 * dni_extra * cos_zenith * tl / tl * Max(ghi, 0); + double abs_am = am * pressure / 101325.0; + + double ghi = exp(-cg2 * abs_am * (fh1 + fh2 * (tl - 1))); + if (perez_enhancement) ghi *= exp(0.01 * pow(abs_am, 1.8)); + + double Gon = 1367 * (1 + 0.033 * cos(360.0 / 365.0 * day_of_year(month, day) * M_PI / + 180)); + + if (dni_extra != 0) Gon = dni_extra; + + ghi = cg1 * Gon * cos_zenith * tl / tl * Max(ghi, 0); double b = 0.664 + 0.163 / fh1; - double bnci = b * exp(-0.09 * absolute_airmass * (tl - 1)); - bnci = dni_extra * Max(bnci, 0); + double bnci = b * exp(-0.09 * abs_am * (tl - 1)); + bnci = Gon * Max(bnci, 0); double bnci_2 = ((1 - (0.1 - 0.2 * exp(-tl)) / (0.1 + 0.882 / fh1)) / cos_zenith); bnci_2 = ghi * Min(Max(bnci_2, 0), 1e20); @@ -1958,6 +1968,12 @@ void irrad::get_irrad(double *ghi, double *dni, double *dhi) { *dhi = diffuseHorizontal; } +void irrad::get_clearsky_irrad(double* ghi_cs, double* dni_cs, double* dhi_cs) { + *ghi_cs = clearskyIrradiance[0]; + *dni_cs = clearskyIrradiance[1]; + *dhi_cs = clearskyIrradiance[2]; +} + void irrad::set_time(int y, int m, int d, int h, double min, double delt_hr) { this->year = y; this->month = m; @@ -2160,7 +2176,7 @@ int irrad::calc() { //clearsky if (enableSubhourlyClipping) { - ineichen(clearskyIrradiance, RTOD * sunAnglesRadians[1], 1.5, 1.0, elevation); + ineichen(clearskyIrradiance, RTOD * sunAnglesRadians[1], month, day, pressure * 100.0, 1.0, elevation, 0, true); } diff --git a/shared/lib_irradproc.h b/shared/lib_irradproc.h index 864aff360..dd772dc9a 100644 --- a/shared/lib_irradproc.h +++ b/shared/lib_irradproc.h @@ -778,7 +778,7 @@ void perez(double hextra, double dn, double df, double alb, double inc, double t * \param[out] clearsky_results[1] clear sky dni * \param[out] clearsky_results[2] clear sky dhi */ -void ineichen(double clearsky_results[3], double apparent_zenith, double absolute_airmass, double linke_turbidity, double altitude, double dni_extra, bool perez_enhancement); +void ineichen(double clearsky_results[3], double apparent_zenith, int month, int day, double pressure, double linke_turbidity, double altitude, double dni_extra, bool perez_enhancement); /** * Isotropic sky model for diffuse irradiance on a tilted surface, see also perez(), hdkr(). * @@ -1168,6 +1168,8 @@ class irrad /// Return the front-side irradiance components void get_irrad(double* ghi, double* dni, double* dhi); + void get_clearsky_irrad(double* ghi_cs, double* dni_cs, double* dhi_cs); + /// Return the effective hour and fraction used for the sun position calculation double get_sunpos_calc_hour(); diff --git a/shared/lib_pv_io_manager.cpp b/shared/lib_pv_io_manager.cpp index 2d7169342..9f17b5e31 100644 --- a/shared/lib_pv_io_manager.cpp +++ b/shared/lib_pv_io_manager.cpp @@ -933,6 +933,10 @@ void PVSystem_IO::AllocateOutputs(compute_module* cm) p_inverterACOutputPreLoss = cm->allocate("ac_gross", numberOfWeatherFileRecords); p_acWiringLoss = cm->allocate("ac_wiring_loss", numberOfWeatherFileRecords); p_ClippingPotential = cm->allocate("clipping_potential", numberOfWeatherFileRecords); + p_CPBin = cm->allocate("clipping_potential_bin", numberOfWeatherFileRecords); + //p_DNIIndex = cm->allocate("dni_index", numberOfWeatherFileRecords); + p_DNIIndexBin = cm->allocate("dni_index_bin", numberOfWeatherFileRecords); + p_transmissionLoss = cm->allocate("ac_transmission_loss", numberOfWeatherFileRecords); p_acPerfAdjLoss = cm->allocate("ac_perf_adj_loss", numberOfWeatherFileRecords); p_acLifetimeLoss = cm->allocate("ac_lifetime_loss", numberOfWeatherFileRecords); diff --git a/shared/lib_pv_io_manager.h b/shared/lib_pv_io_manager.h index c745f48ad..45847a9d9 100644 --- a/shared/lib_pv_io_manager.h +++ b/shared/lib_pv_io_manager.h @@ -422,6 +422,9 @@ struct PVSystem_IO ssc_number_t *p_subhourlyClippingLoss; ssc_number_t* p_subhourlyClippingLossFactor; ssc_number_t* p_ClippingPotential; + //ssc_number_t* p_DNIIndex; + ssc_number_t* p_CPBin; + ssc_number_t* p_DNIIndexBin; }; /** diff --git a/ssc/cmod_pvsamv1.cpp b/ssc/cmod_pvsamv1.cpp index c10f3961d..545f576f7 100644 --- a/ssc/cmod_pvsamv1.cpp +++ b/ssc/cmod_pvsamv1.cpp @@ -652,6 +652,7 @@ static var_info _cm_vtab_pvsamv1[] = { { SSC_OUTPUT, SSC_ARRAY, "subarray1_dc_gross", "Subarray 1 DC power gross", "kW", "", "Time Series (Subarray 1)", "*", "", "" }, { SSC_OUTPUT, SSC_ARRAY, "subarray1_voc", "Subarray 1 Open circuit DC voltage", "V", "", "Time Series (Subarray 1)", "", "", "" }, { SSC_OUTPUT, SSC_ARRAY, "subarray1_isc", "Subarray 1 String short circuit DC current", "A", "", "Time Series (Subarray 1)", "", "", "" }, + { SSC_OUTPUT, SSC_ARRAY, "subarray1_dni_index", "Subarray 1 DNI Index", "W/m2", "", "Time Series (Subarray 1)", "", "", "" }, { SSC_OUTPUT, SSC_ARRAY, "subarray2_surf_tilt", "Subarray 2 Surface tilt", "degrees", "", "Time Series (Subarray 2)", "", "", "" }, { SSC_OUTPUT, SSC_ARRAY, "subarray2_surf_azi", "Subarray 2 Surface azimuth", "degrees", "", "Time Series (Subarray 2)", "", "", "" }, @@ -805,6 +806,8 @@ static var_info _cm_vtab_pvsamv1[] = { { SSC_OUTPUT, SSC_ARRAY, "ac_gross", "Inverter AC output power", "kW", "", "Time Series (Array)", "*", "", "" }, { SSC_OUTPUT, SSC_ARRAY, "clipping_potential", "Clipping potential", "", "", "Time Series (Inverter)", "", "", "" }, { SSC_OUTPUT, SSC_ARRAY, "subhourly_clipping_loss", "Subhourly clipping correction loss", "kW", "", "Time Series (Inverter)", "", "", "" }, + { SSC_OUTPUT, SSC_ARRAY, "clipping_potential_bin", "Clipping potential bin", "", "", "Time Series (Inverter)", "", "", "" }, + { SSC_OUTPUT, SSC_ARRAY, "dni_index_bin", "DNI Index bin", "", "", "Time Series (Inverter)", "", "", "" }, // transformer model outputs { SSC_OUTPUT, SSC_ARRAY, "xfmr_nll_ts", "Transformer no load loss", "kW", "", "Time Series (Transformer)", "", "", "" }, @@ -960,6 +963,7 @@ static var_info _cm_vtab_pvsamv1[] = { { SSC_OUTPUT, SSC_NUMBER, "annual_ac_inv_eff_loss_percent", "AC inverter efficiency loss", "%", "", "Loss", "", "", "" }, { SSC_OUTPUT, SSC_NUMBER, "annual_ac_wiring_loss_percent", "AC wiring loss", "%", "", "Loss", "", "", "" }, { SSC_OUTPUT, SSC_NUMBER, "annual_subhourly_clipping_loss_percent", "Subhourly clipping correction loss percent", "%", "", "Loss", "", "", "" }, + { SSC_OUTPUT, SSC_NUMBER, "nominal_annual_clipping_output", "Test output for nominal annual clipping output", "kWh", "", "", "", "", "" }, { SSC_OUTPUT, SSC_NUMBER, "annual_transmission_loss_percent", "AC transmission loss", "%", "", "Loss", "", "", "" }, // { SSC_OUTPUT, SSC_NUMBER, "annual_ac_transformer_loss_percent", "AC step-up transformer loss", "%", "", "Loss", "", "", "" }, @@ -1532,8 +1536,11 @@ void cm_pvsamv1::exec() // beam, skydiff, and grounddiff IN THE PLANE OF ARRAY (W/m2) double ibeam, iskydiff, ignddiff; double ibeam_csky, iskydiff_csky, ignddiff_csky; + double ghi_cs, dni_cs, dhi_cs; double aoi, stilt, sazi, rot, btd; + irr.get_clearsky_irrad(&ghi_cs, &dni_cs, &dhi_cs); + // Ensure that the usePOAFromWF flag is false unless a reference cell has been used. // This will later get forced to false if any shading has been applied (in any scenario) // also this will also be forced to false if using the cec mcsp thermal model OR if using the spe module model with a diffuse util. factor < 1.0 @@ -1918,8 +1925,8 @@ void cm_pvsamv1::exec() PVSystem->p_poaBeamFrontCS[nn][idx] = (ssc_number_t)ibeam_csky; PVSystem->p_poaDiffuseFrontCS[nn][idx] = (ssc_number_t)(iskydiff_csky); PVSystem->p_poaDiffuseFrontCS[nn][idx] = (ssc_number_t)(ignddiff_csky); - if (ibeam_csky != 0) { - PVSystem->p_DNIIndex[nn][idx] = (ssc_number_t)(ibeam / ibeam_csky); + if (dni_cs != 0) { + PVSystem->p_DNIIndex[nn][idx] = (ssc_number_t)(Irradiance->p_weatherFileDNI[idx] / dni_cs); } else { PVSystem->p_DNIIndex[nn][idx] = 0; @@ -2591,6 +2598,7 @@ void cm_pvsamv1::exec() annual_subhourly_clipping_loss += sub_clipping_matrix.at(dni_row, clip_pot_col); */ } + assign("nominal_annual_clipping_output", nominal_annual_clipping_output); } for (size_t iyear = 0; iyear < nyears; iyear++) @@ -2714,6 +2722,7 @@ void cm_pvsamv1::exec() double clip_pot = (dcPower_kW_csky - paco ) / (paco); PVSystem->p_ClippingPotential[idx] = clip_pot; + //PVSystem->p_DNIIndex[idx] = dni_clearness_index; //Lookup matrix for percentage effect based on DNI index, Clipping potential @@ -2731,6 +2740,7 @@ void cm_pvsamv1::exec() } } } + PVSystem->p_DNIIndexBin[idx] = dni_row; //Clipping potential indexing if (clip_pot < sub_clipping_matrix.at(0, 1)) clip_pot_col = 1; @@ -2742,9 +2752,10 @@ void cm_pvsamv1::exec() } } } + PVSystem->p_CPBin[idx] = clip_pot_col; //acpwr_gross *= (1 - sub_clipping_matrix.at(dni_row, clip_pot_col)); - if (dcPower_kW_csky > 0.0) { + if (dcPower_kW > 0.0) { ac_subhourlyclipping_loss = sub_clipping_matrix.at(dni_row, clip_pot_col) * nominal_annual_clipping_output; //ac_subhourlyclipping_loss = sub_clipping_matrix.at(dni_row, clip_pot_col) * sharedInverter->powerAC_kW_clipping; }