Merge branch 'patch' into ortools

.gitattributes

@@ -1,11 +1,15 @@
 
-doc/* linguist-documentation=true
+doc/** linguist-documentation=true
 tcsdata/docs/* linguist-documentation=true
 
+
 *.c linguist-vendored
-lpsolve/* linguist-vendored=true
-nlopt/* linguist-vendored=true
-examples/* linguist-vendored=true
+lpsolve/** linguist-vendored=true
+nlopt/** linguist-vendored=true
+samples/** linguist-vendored=true
+test/input_docs/** linguist-vendored=true
+splinter/** linguist-vendored=true
+rapidjson/** linguist-vendored=true
 
 build_*/* linguist_generated=true
 

shared/lib_battery_dispatch_manual.cpp

@@ -116,11 +116,11 @@ void dispatch_manual_t::prepareDispatch(size_t hour_of_year, size_t )
 	m_batteryPower->canGridCharge = _gridcharge_array[iprofile - 1];
     m_batteryPower->canClipCharge = _can_clip_charge;
 
-	if (iprofile < _fuelcellcharge_array.size()) {
+	if (iprofile <= _fuelcellcharge_array.size()) {
 		m_batteryPower->canFuelCellCharge = _fuelcellcharge_array[iprofile - 1];
 	}
 
-    if (iprofile < _discharge_grid_array.size()) {
+    if (iprofile <= _discharge_grid_array.size()) {
         m_batteryPower->canDischargeToGrid = _discharge_grid_array[iprofile - 1];
     }
 

shared/lib_battery_powerflow.cpp

@@ -36,6 +36,9 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 #include "lib_power_electronics.h"
 #include "lib_shared_inverter.h"
 
+// 0.005 W when applied to power. Some inverters have a night time loss of 0.01 W
+double powerflow_tolerance = 0.000005;
+
 BatteryPower::BatteryPower(double dtHour) :
 		dtHour(dtHour),
 		powerSystem(0),
@@ -102,7 +105,7 @@ BatteryPower::BatteryPower(double dtHour) :
 		depthOfDischargeMax(1),
         currentChargeMax(0),
         currentDischargeMax(0),
-		tolerance(0.001){}
+		tolerance(powerflow_tolerance){}
 
 BatteryPower::BatteryPower(const BatteryPower& orig) {
     sharedInverter = orig.sharedInverter;

shared/lib_battery_powerflow.h

@@ -42,6 +42,9 @@ class SharedInverter;
 
 struct BatteryPower;
 
+// Allow for consistent power tolerances between the technology code and utility rate code
+extern double powerflow_tolerance; 
+
 /**
 * \class BatteryPowerFlow
 *

shared/lib_geothermal.cpp

@@ -99,7 +99,7 @@ namespace geothermal
 	/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
 	// GETEM Physics and general equations
 	/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
-	const bool IMITATE_GETEM = false;
+	const bool IMITATE_GETEM = false; 
 	const double GETEM_FT_IN_METER = (IMITATE_GETEM) ? 3.28083 : physics::FT_PER_METER; // feet per meter - largest source of discrepancy
 	//const double GETEM_PSI_PER_BAR = (IMITATE_GETEM) ? 14.50377 : physics::PSI_PER_BAR; // psi per bar
     const double GETEM_PSI_PER_BAR = 14.50377; // psi per bar

shared/lib_irradproc.cpp

@@ -2193,9 +2193,7 @@ int irrad::calc() {
     }
 
     //clearsky
-    if (enableSubhourlyClipping) {
-        ineichen(clearskyIrradiance, RTOD * sunAnglesRadians[1], month, day, pressure * 100.0, 1.0, elevation, 0, true);
-    }
+    ineichen(clearskyIrradiance, RTOD * sunAnglesRadians[1], month, day, pressure * 100.0, 1.0, elevation, 0, true);
 
 
     planeOfArrayIrradianceFront[0] = planeOfArrayIrradianceFront[1] = planeOfArrayIrradianceFront[2] = 0;
@@ -2382,6 +2380,7 @@ int irrad::calc_rear_side(double transmissionFactor, double groundClearanceHeigh
         getBackSurfaceIrradiances(pvBackShadeFraction, rowToRow, verticalHeight, clearanceGround, distanceBetweenRows,
                                   horizontalLength, rearGroundGHI, frontGroundGHI, frontReflected,
                                   rearIrradiancePerCellrow, rearAverageIrradiance);
+
         getBackSurfaceIrradiancesCS(pvBackShadeFraction, rowToRow, verticalHeight, clearanceGround, distanceBetweenRows,
             horizontalLength, rearGroundGHI, frontGroundGHI, frontReflected,
             rearIrradiancePerCellrowCS, rearAverageIrradianceCS);
@@ -3069,16 +3068,16 @@ void irrad::getBackSurfaceIrradiancesCS(double pvBackShadeFraction, double rowTo
 
     // Calculate diffuse isotropic irradiance for a horizontal surface
     perez(0, clearskyIrradiance[1], clearskyIrradiance[2], albedo, solarZenithRadians, 0, solarZenithRadians,
-        planeOfArrayIrradianceRear, diffuseIrradianceRear);
-    double isotropicSkyDiffuse = diffuseIrradianceRear[0];
+        planeOfArrayIrradianceRearCS, diffuseIrradianceRearCS);
+    double isotropicSkyDiffuse = diffuseIrradianceRearCS[0];
 
     // Calculate components for a 90 degree tilt to get horizon brightening
     double surfaceAnglesRadians90[5] = { 0, 0, 0, 0, 0 };
     incidence(0, 90.0, 180.0, 45.0, solarZenithRadians, solarAzimuthRadians, this->enableBacktrack,
         this->groundCoverageRatio, this->slopeTilt, this->slopeAzm, this->forceToStow, this->stowAngleDegrees, this->useCustomRotAngles, this->customRotAngle, surfaceAnglesRadians90);
     perez(0, clearskyIrradiance[1], clearskyIrradiance[2], albedo, surfaceAnglesRadians90[0],
-        surfaceAnglesRadians90[1], solarZenithRadians, planeOfArrayIrradianceRear, diffuseIrradianceRear);
-    double horizonDiffuse = diffuseIrradianceRear[2];
+        surfaceAnglesRadians90[1], solarZenithRadians, planeOfArrayIrradianceRearCS, diffuseIrradianceRear);
+    double horizonDiffuse = diffuseIrradianceRearCS[2];
 
     // Calculate x,y coordinates of bottom and top edges of PV row in back of desired PV row so that portions of sky and ground viewed by the
     // PV cell may be determined. Origin of x-y axis is the ground point below the lower front edge of the desired PV row. The row in back of
@@ -3092,13 +3091,13 @@ void irrad::getBackSurfaceIrradiancesCS(double pvBackShadeFraction, double rowTo
 
 // Calculate diffuse and direct component irradiances for each cell row (assuming 6 rows)
     std::vector<double> rearDirectDiffuse;              // the direct and sky diffuse irradiance incident on the rear for each cell row, before losses (shading, soiling, etc.)
-    poaRearDirectDiffuse = 0.;                          // the average direct and sky diffuse irradiance incident on the rear, before losses (shading, soiling, etc.)
+    poaRearDirectDiffuseCS = 0.;                          // the average direct and sky diffuse irradiance incident on the rear, before losses (shading, soiling, etc.)
     std::vector<double> rearRowReflections;             // the reflected irradiance from the rear row on the rear of each cell row
-    poaRearRowReflections = 0.;                         // the average reflected irradiance from the rear row on the rear
+    poaRearRowReflectionsCS = 0.;                         // the average reflected irradiance from the rear row on the rear
     std::vector<double> rearGroundReflected;            // the ground reflected irradiance onto the rear of each cell row, considering view factor
-    poaRearGroundReflected = 0.;                        // the average ground reflected irradiance onto the rear, considering view factor
+    poaRearGroundReflectedCS = 0.;                        // the average ground reflected irradiance onto the rear, considering view factor
     std::vector<double> rearSelfShaded;                 // the direct and circumsolar shaded from being incident on the rear, for each cell
-    poaRearSelfShaded = 0.;                             // the average direct and circumsolar shaded from being incident on the rear
+    poaRearSelfShadedCS = 0.;                             // the average direct and circumsolar shaded from being incident on the rear
     size_t cellRows = poaRearIrradRes;
     for (size_t i = 0; i != cellRows; i++) {
         // Calculate diffuse irradiances and reflected amounts for each cell row over its field of view of 180 degrees,
@@ -3277,9 +3276,9 @@ void irrad::getBackSurfaceIrradiancesCS(double pvBackShadeFraction, double rowTo
             solarAzimuthRadians, this->enableBacktrack,
             this->groundCoverageRatio, this->slopeTilt, this->slopeAzm, this->forceToStow, this->stowAngleDegrees, this->useCustomRotAngles, this->customRotAngle, surfaceAnglesRadians);
         perez(0, clearskyIrradiance[1], clearskyIrradiance[2], albedo, surfaceAnglesRadians[0],
-            surfaceAnglesRadians[1], solarZenithRadians, planeOfArrayIrradianceRear, diffuseIrradianceRear);
+            surfaceAnglesRadians[1], solarZenithRadians, planeOfArrayIrradianceRearCS, diffuseIrradianceRearCS);
 
-        double rear_direct_circumsolar = planeOfArrayIrradianceRear[0] + diffuseIrradianceRear[1];
+        double rear_direct_circumsolar = planeOfArrayIrradianceRearCS[0] + diffuseIrradianceRearCS[1];
         rearDirectDiffuse[i] += rear_direct_circumsolar;
 
         double cellShade = pvBackShadeFraction * cellRows - i;
@@ -3301,10 +3300,10 @@ void irrad::getBackSurfaceIrradiancesCS(double pvBackShadeFraction, double rowTo
         }
 
         rearAverageIrradiance += rearIrradiance[i] / cellRows;
-        poaRearDirectDiffuse += rearDirectDiffuse[i] / cellRows;
-        poaRearRowReflections += rearRowReflections[i] / cellRows;
-        poaRearSelfShaded += rearSelfShaded[i] / cellRows;
-        poaRearGroundReflected += rearGroundReflected[i] / cellRows;
+        poaRearDirectDiffuseCS += rearDirectDiffuse[i] / cellRows;
+        poaRearRowReflectionsCS += rearRowReflections[i] / cellRows;
+        poaRearSelfShadedCS += rearSelfShaded[i] / cellRows;
+        poaRearGroundReflectedCS += rearGroundReflected[i] / cellRows;
         double xy = 1.;
     }
 

shared/lib_irradproc.h

@@ -1024,16 +1024,22 @@ class irrad
     double poaRearRowReflections;           ///< Rear row reflected irradiance on rear (W/m2)
     double poaRearGroundReflected;          ///< Ground reflected irradiance onto the rear (W/m2)
     double poaRearSelfShaded;               ///< Irradiance shaded from being incident on the rear (W/m2)
+    double poaRearDirectDiffuseCS;            ///< Direct and sky diffuse irradiance on rear (W/m2)
+    double poaRearRowReflectionsCS;           ///< Rear row reflected irradiance on rear (W/m2)
+    double poaRearGroundReflectedCS;          ///< Ground reflected irradiance onto the rear (W/m2)
+    double poaRearSelfShadedCS;               ///< Irradiance shaded from being incident on the rear (W/m2)
 
     // Outputs
     double sunAnglesRadians[9];				///< Sun angles in radians calculated from solarpos()	
     double surfaceAnglesRadians[5];			///< Surface angles in radians calculated from incidence()
     double planeOfArrayIrradianceFront[3];	///< Front-side plane-of-array irradiance for beam, sky diffuse, ground diffuse (W/m2)
     double planeOfArrayIrradianceFrontCS[3]; ///< Front-side plane-of-array clearsky irradiance for beam, sky diffuse, ground diffuse (W/m2)
     double planeOfArrayIrradianceRear[3];	///< Rear-side plane-of-array irradiance for beam, sky diffuse, ground diffuse (W/m2)
+    double planeOfArrayIrradianceRearCS[3];
     double diffuseIrradianceFront[3];		///< Front-side diffuse irradiance for isotropic, circumsolar, and horizon (W/m2)
     double diffuseIrradianceFrontCS[3];		///< Front-side diffuse clearsky irradiance for isotropic, circumsolar, and horizon (W/m2)
     double diffuseIrradianceRear[3];		///< Rear-side diffuse irradiance for isotropic, circumsolar, and horizon (W/m2)
+    double diffuseIrradianceRearCS[3];
     int timeStepSunPosition[3];				///< [0] effective hour of day used for sun position, [1] effective minute of hour used for sun position, [2] is sun up?  (0=no, 1=midday, 2=sunup, 3=sundown)
     double planeOfArrayIrradianceRearAverage; ///< Average rear side plane-of-array irradiance (W/m2)
     double planeOfArrayIrradianceRearAverageCS; ///< Average rear side clearsky plane-of-array irradiance (W/m2)

shared/lib_pv_io_manager.cpp

@@ -505,7 +505,7 @@ Subarray_IO::Subarray_IO(compute_module* cm, const std::string& cmName, size_t s
             if (trackMode == irrad::SEASONAL_TILT)
                 throw exec_error(cmName, "Time-series tilt input may not be used with the snow model at this time: subarray " + util::to_string((int)(subarrayNumber)));
             // if tracking mode is 1-axis tracking, don't need to limit tilt angles
-            if (snowModel.setup(selfShadingInputs.nmody, (float)tiltDegrees, (trackMode != irrad::SINGLE_AXIS))) {
+            if (snowModel.setup(selfShadingInputs.nmody, (float)tiltDegrees, cm->as_double("snow_slide_coefficient"), (trackMode != irrad::SINGLE_AXIS))) {
                 if (!snowModel.good) {
                     cm->log(snowModel.msg, SSC_ERROR);
                 }
@@ -973,10 +973,15 @@ void PVSystem_IO::AllocateOutputs(compute_module* cm)
     p_dcLifetimeLoss = cm->allocate("dc_lifetime_loss", numberOfWeatherFileRecords);
     p_systemDCPower = cm->allocate("dc_net", numberOfLifetimeRecords);
     p_systemACPower = cm->allocate("gen", numberOfLifetimeRecords);
+    p_systemACPowerMax = cm->allocate("ac_csky_max", numberOfLifetimeRecords);
 
     p_systemDCPowerCS = cm->allocate("dc_net_clearsky", numberOfLifetimeRecords);
-    p_subhourlyClippingLoss = cm->allocate("subhourly_clipping_loss", numberOfLifetimeRecords);
-    p_subhourlyClippingLossFactor = cm->allocate("subhourly_clipping_loss_factor", numberOfLifetimeRecords);
+    if (cm->as_boolean("enable_subhourly_clipping")) {
+        p_subhourlyClippingLoss = cm->allocate("subhourly_clipping_loss", numberOfLifetimeRecords);
+    }
+    if (cm->as_boolean("enable_subinterval_distribution")) {
+        p_DistributionClippingLoss = cm->allocate("distribution_clipping_loss", numberOfLifetimeRecords);
+    }
 
     if (Simulation->useLifetimeOutput)
     {

shared/lib_pv_io_manager.h

@@ -419,10 +419,12 @@ struct PVSystem_IO
 	ssc_number_t *p_systemDCPower; // kWdc
     ssc_number_t* p_systemDCPowerCS; // kWdc
 	ssc_number_t *p_systemACPower; // kWac
+    ssc_number_t* p_systemACPowerMax; //kWac
 
     ssc_number_t *p_subhourlyClippingLoss;
     ssc_number_t* p_subhourlyClippingLossFactor;
     ssc_number_t* p_ClippingPotential;
+    ssc_number_t* p_DistributionClippingLoss;
     //ssc_number_t* p_DNIIndex;
     ssc_number_t* p_CPBin;
     ssc_number_t* p_DNIIndexBin;

shared/lib_pvshade.cpp

@@ -201,6 +201,7 @@ void diffuse_reduce(
 
     // sky diffuse reduction
     Fskydiff = skydiffderates.lookup(stilt);
+    Fskydiff = fmin(Fskydiff, 1.0);
     reduced_skydiff = Fskydiff * poa_sky;
 
 	double solalt = 90 - solzen;

shared/lib_shared_inverter.cpp

@@ -230,7 +230,7 @@ void SharedInverter::calculateTempDerate(double V, double tempC, double& p_dc_ra
 
 double SharedInverter::getInverterDCMaxPower(double p_dc_rated)
 {
-    double inv_dc_max_power;
+    double inv_dc_max_power = p_dc_rated * util::kilowatt_to_watt; //if the inverter type isn't one of the following, assume that max power is equal to rated power"
     if (m_inverterType == SANDIA_INVERTER || m_inverterType == DATASHEET_INVERTER || m_inverterType == COEFFICIENT_GENERATOR)
         //m_sandiaInverter->acpower(std::fabs(powerDC_Watts) / m_numInverters, DCStringVoltage, &powerAC_Watts, &P_par, &P_lr, &efficiencyAC, &powerClipLoss_kW, &powerConsumptionLoss_kW, &powerNightLoss_kW);
         inv_dc_max_power = m_sandiaInverter->Pdco;
@@ -240,9 +240,6 @@ double SharedInverter::getInverterDCMaxPower(double p_dc_rated)
     else if (m_inverterType == OND_INVERTER)
         //m_ondInverter->acpower(std::fabs(powerDC_Watts) / m_numInverters, DCStringVoltage, tempC, &powerAC_Watts, &P_par, &P_lr, &efficiencyAC, &powerClipLoss_kW, &powerConsumptionLoss_kW, &powerNightLoss_kW, &dcWiringLoss_ond_kW, &acWiringLoss_ond_kW);
         inv_dc_max_power = m_ondInverter->PMaxDC;
-    else if (m_inverterType == NONE) {
-        inv_dc_max_power = p_dc_rated * util::kilowatt_to_watt;
-    }
 
     return inv_dc_max_power;
 }

shared/lib_snowmodel.cpp

@@ -62,7 +62,7 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 pvsnowmodel::pvsnowmodel()
 {
 	mSlope = -80;
-	sSlope = (float)1.97;
+	//sSlope = (float)1.97;
 	deltaThreshold = 1.00;
 	depthThreshold = 1.00;
 	previousDepth = 0;
@@ -76,11 +76,11 @@ pvsnowmodel::pvsnowmodel()
 
 }
 
-bool pvsnowmodel::setup(int nmody_in, float baseTilt_in, bool limitTilt){
+bool pvsnowmodel::setup(int nmody_in, float baseTilt_in, float snow_slide_coeff, bool limitTilt){
 
 	nmody = nmody_in;
 	baseTilt = baseTilt_in;
-
+    sSlope = snow_slide_coeff;
 	if(limitTilt && (baseTilt>45 || baseTilt < 10)){
 		good = true;
 		msg = util::format("The snow model is designed to work for PV arrays with a tilt angle between 10 and 45 degrees, but will generate results for tilt angles outside this range. The system you are modeling includes a subarray tilt angle of %f degrees.", baseTilt);

shared/lib_snowmodel.h

@@ -41,7 +41,7 @@ class pvsnowmodel
 	pvsnowmodel();
 
 	// limitTilt requires tilt to be between 10 and 45 degrees
-	bool setup(int, float, bool limitTilt = true);
+	bool setup(int, float, float snow_slide_coeff = 1.97, bool limitTilt = true);
 
 	bool getLoss(float poa, float tilt, float wspd, float tdry, float snowDepth, int sunup, float dt, float &returnLoss);
 
@@ -63,4 +63,4 @@ class pvsnowmodel
 							//  if an error has occured
 };
 
-#endif
+#endif

shared/lib_utility_rate_equations.cpp

@@ -1352,6 +1352,7 @@ void forecast_setup::setup(rate_data* rate, std::vector<double>& P_pv_ac, std::v
     }
     for (size_t idx = 0; idx < num_recs && idx < array_size; idx++)
     {
+        size_t year_one_index = util::yearOneIndex(1.0 / _steps_per_hour, idx);
         double grid_power = P_pv_ac[idx] - P_load_ac[idx];
 
         gross_load_during_month += P_load_ac[idx] * _dt_hour;
@@ -1367,7 +1368,7 @@ void forecast_setup::setup(rate_data* rate, std::vector<double>& P_pv_ac, std::v
         }
 
         if (rate->dc_enabled) {
-            int dc_tou_period = rate->get_dc_tou_row(hour_of_year + step, curr_month - 1);
+            int dc_tou_period = rate->get_dc_tou_row(year_one_index, curr_month - 1);
             size_t month_idx = year * 12 + (curr_month - 1);
             double peak = monthly_peaks.at(month_idx, dc_tou_period) - peak_offset; // Peak for dispatch calcs in battery: peak minus battery capacity
             if (-1.0 * grid_power > peak) {

shared/lib_weatherfile.cpp

@@ -116,7 +116,7 @@ static float col_or_nan(const std::string& s)
         }
     }
     else
-        return std::numeric_limits<float>::quiet_NaN();;
+        return std::numeric_limits<float>::quiet_NaN();
 }
 
 static double conv_deg_min_sec(double degrees,
@@ -824,9 +824,11 @@ bool weatherfile::open(const std::string& file, bool header_only)
             {
                 m_hdr.lon = col_or_nan(value);
             }
-            else if (name == "tz" || name == "timezone" || name == "time zone") // require "time zone" and "local time zone" columns in NSRDB files are the same
+            else if (name == "tz" || name == "timezone" || name == "time zone" || name == "local time zone") // require "time zone" and "local time zone" columns in NSRDB files are the same
             {
-                m_hdr.tz = col_or_nan(value);
+                // some nsrdb endpoints like nsrdb-msg-v1-0-0 have both "local time zone" and "time zone" with "time zone" empty, so we need to read "local time zone" and ignore "time zone"
+                if (std::isnan(m_hdr.tz))  // only assign value if one was not assigned in an earlier pass
+                    m_hdr.tz = col_or_nan(value);
             }
             else if (name == "el" || name == "elev" || name == "elevation" || name == "site elevation" || name == "altitude")
             {