Fix array bounds error for interzone windows and fix convexity of mirrored surfaces

doc/input-output-reference/src/overview/group-thermal-zone-description-geometry.tex

@@ -1983,7 +1983,7 @@ \subsection{Surface Vertices}\label{surface-vertices}
 use the same vertex input. The numeric parameters indicated below are taken from the \hyperref[buildingsurfacedetailed]{BuildingSurface:Detailed} definition; the others may not be exactly the same but are identical in configuration. They are also ``extensible'' -- so, to define more vertices for these surfaces, simply add the required number of vertices (X, Y, and Z coordinates for each vertex) to the input file. Note that \hyperref[fenestrationsurfacedetailed]{FenestrationSurface:Detailed} is not extensible and is limited to 4 (max) vertices. If the Number of Surface Vertex groups is left blank or entered as \textbf{autocalculate}, EnergyPlus looks at the number of groups entered and figures out how many coordinate groups are entered.
 
 \begin{callout}
-\warning{Note that the resolution on the surface vertex input is 1 millimeter (.001 meter). Therefore, using vertices that are very close together (\textless{} 1 mm) may result in invalid dot product and fatal errors during shading calculations.}
+\warning{Note that the resolution for surface vertex input is 1 centimeter (0.01 meter). Vertices that are \textless{} 1 cm apart will be combined.}
 \end{callout}
 
 \begin{figure}[hbtp] % fig 29

src/EnergyPlus/DataGlobalConstants.hh

@@ -571,7 +571,12 @@ namespace Constant {
     Real64 constexpr OneFifth = 1.0 / 5.0;   // 1/5 in highest precision
     Real64 constexpr OneSixth = 1.0 / 6.0;   // 1/6 in highest precision
     Real64 constexpr FourFifths = 4.0 / 5.0; // 4/5 in highest precision
+    Real64 constexpr OneThousandth = 1.0e-3; // Used as a tolerance in various places
+    Real64 constexpr OneMillionth = 1.0e-6;  // Used as a tolerance in various places
 
+    Real64 constexpr OneCentimeter = 0.01;     // Geometric tolerance in meters
+    Real64 constexpr TwoCentimeters = 0.02;    // Geometric tolerance in meters
+    Real64 constexpr SmallDistance = 1.0e-4;   // Geometric tolerance in meters
     Real64 constexpr MaxEXPArg = 709.78;       // maximum exponent in EXP() function
     Real64 constexpr Pi = 3.14159265358979324; // Pi 3.1415926535897932384626435
     Real64 constexpr PiOvr2 = Pi / 2.0;        // Pi/2

src/EnergyPlus/DataSurfaces.cc

@@ -74,7 +74,6 @@ namespace EnergyPlus::DataSurfaces {
 //                      Dec 2006, DJS (PSU) added logical ecoroof variable
 //                      Dec 2008, TH added new properties to SurfaceWindowCalc for thermochromic windows
 //                      Jul 2011, M.J. Witte and C.O. Pedersen, add new fields to OSC for last T, max and min
-//       RE-ENGINEERED  na
 
 // Using/Aliasing
 using namespace DataVectorTypes;
@@ -88,15 +87,6 @@ using namespace WindowManager;
 
 Array1D_string const cExtBoundCondition({-6, 0}, {"KivaFoundation", "FCGround", "OSCM", "OSC", "OSC", "Ground", "ExternalEnvironment"});
 
-// Parameters to indicate surface classes
-// Surface Class (FLOOR, WALL, ROOF (incl's CEILING), WINDOW, DOOR, GLASSDOOR,
-// SHADING (includes OVERHANG, WING), DETACHED, INTMASS),
-// TDD:DOME, TDD:DIFFUSER (for tubular daylighting device)
-// (Note: GLASSDOOR and TDD:DIFFUSER get overwritten as WINDOW
-// in SurfaceGeometry.cc, SurfaceWindow%OriginalClass holds the true value)
-// why aren't these sequential (LKL - 13 Aug 2007)
-
-// Constructor
 Surface2D::Surface2D(ShapeCat const shapeCat, int const axis, Vertices const &v, Vector2D const &vl, Vector2D const &vu)
     : axis(axis), vertices(v), vl(vl), vu(vu)
 {
@@ -180,7 +170,7 @@ Surface2D::Surface2D(ShapeCat const shapeCat, int const axis, Vertices const &v,
                 xt = xte;
             }
 #endif
-            assert((shapeCat == ShapeCat::Nonconvex) || (crossEdges.size() == 2));
+            assert((shapeCat == ShapeCat::Nonconvex) || (crossEdges.size() == 2u));
             for (auto const &edge : crossEdges) {
                 size_type const iEdge(std::get<2>(edge));
                 slab.edges.push_back(iEdge); // Add edge to slab
@@ -470,7 +460,7 @@ Surface2D SurfaceData::computed_surface2d() const
 
 Real64 SurfaceData::get_average_height(EnergyPlusData &state) const
 {
-    if (std::abs(SinTilt) < 1.e-4) {
+    if (std::abs(SinTilt) < Constant::SmallDistance) {
         return 0.0;
     }
     using Vertex2D = ObjexxFCL::Vector2<Real64>;

src/EnergyPlus/HeatBalanceSurfaceManager.cc

@@ -409,7 +409,7 @@ void InitSurfaceHeatBalance(EnergyPlusData &state)
         if (state.dataHeatBalSurfMgr->InitSurfaceHeatBalancefirstTime) {
             DisplayString(state, "Computing Interior Diffuse Solar Exchange through Interzone Windows");
         }
-        ComputeDifSolExcZonesWIZWindows(state, state.dataGlobal->NumOfZones);
+        ComputeDifSolExcZonesWIZWindows(state);
     }
 
     Dayltg::initDaylighting(state, state.dataHeatBalSurfMgr->InitSurfaceHeatBalancefirstTime);
@@ -4397,21 +4397,22 @@ void ComputeIntSWAbsorpFactors(EnergyPlusData &state)
     } // End of enclosure loop
 }
 
-void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state, int const NumberOfEnclosures) // Number of solar enclosures
+void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state)
 {
 
     // SUBROUTINE INFORMATION:
     //       MODIFIED       Jun 2007 - Lawrie - Speed enhancements.
     //       RE-ENGINEERED  Winkelmann, Lawrie
 
     // PURPOSE OF THIS SUBROUTINE:
-    // This subroutine computes the diffuse solar exchange factors between zones with
+    // This subroutine computes the diffuse solar exchange factors between enclosures with
     // interzone windows.
 
+    int const numEnclosures = state.dataViewFactor->NumOfSolarEnclosures;
     if (!allocated(state.dataHeatBalSurf->ZoneFractDifShortZtoZ)) {
-        state.dataHeatBalSurf->ZoneFractDifShortZtoZ.allocate(NumberOfEnclosures, NumberOfEnclosures);
-        state.dataHeatBalSurf->EnclSolRecDifShortFromZ.allocate(NumberOfEnclosures);
-        state.dataHeatBalSurfMgr->DiffuseArray.allocate(NumberOfEnclosures, NumberOfEnclosures);
+        state.dataHeatBalSurf->ZoneFractDifShortZtoZ.allocate(numEnclosures, numEnclosures);
+        state.dataHeatBalSurf->EnclSolRecDifShortFromZ.allocate(numEnclosures);
+        state.dataHeatBalSurfMgr->DiffuseArray.allocate(numEnclosures, numEnclosures);
     }
 
     state.dataHeatBalSurf->EnclSolRecDifShortFromZ = false;
@@ -4438,10 +4439,10 @@ void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state, int const NumberOfEn
     //          Compute fractions for multiple passes.
 
     Array2D<Real64>::size_type l(0u), m(0u), d(0u);
-    for (int NZ = 1; NZ <= NumberOfEnclosures; ++NZ, d += NumberOfEnclosures + 1) {
+    for (int NZ = 1; NZ <= numEnclosures; ++NZ, d += numEnclosures + 1) {
         m = NZ - 1;
         Real64 D_d(0.0); // Local accumulator
-        for (int MZ = 1; MZ <= NumberOfEnclosures; ++MZ, ++l, m += NumberOfEnclosures) {
+        for (int MZ = 1; MZ <= numEnclosures; ++MZ, ++l, m += numEnclosures) {
             if (MZ == NZ) continue;
             state.dataHeatBalSurfMgr->DiffuseArray[l] =
                 state.dataHeatBalSurf->ZoneFractDifShortZtoZ[l] /
@@ -4454,13 +4455,12 @@ void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state, int const NumberOfEn
 
     state.dataHeatBalSurf->ZoneFractDifShortZtoZ = state.dataHeatBalSurfMgr->DiffuseArray;
     // added for CR 7999 & 7869
-    assert(state.dataHeatBalSurf->ZoneFractDifShortZtoZ.isize1() == NumberOfEnclosures);
-    assert(state.dataHeatBalSurf->ZoneFractDifShortZtoZ.isize2() == NumberOfEnclosures);
-    l = 0u;
-    for (int NZ = 1; NZ <= NumberOfEnclosures; ++NZ) {
-        for (int MZ = 1; MZ <= NumberOfEnclosures; ++MZ, ++l) {
+    assert(state.dataHeatBalSurf->ZoneFractDifShortZtoZ.isize1() == numEnclosures);
+    assert(state.dataHeatBalSurf->ZoneFractDifShortZtoZ.isize2() == numEnclosures);
+    for (int NZ = 1; NZ <= numEnclosures; ++NZ) {
+        for (int MZ = 1; MZ <= numEnclosures; ++MZ) {
             if (MZ == NZ) continue;
-            if (state.dataHeatBalSurf->ZoneFractDifShortZtoZ[l] > 0.0) { // [ l ] == ( MZ, NZ )
+            if (state.dataHeatBalSurf->ZoneFractDifShortZtoZ(MZ, NZ) > 0.0) {
                 state.dataHeatBalSurf->EnclSolRecDifShortFromZ(NZ) = true;
                 break;
             }
@@ -4469,23 +4469,23 @@ void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state, int const NumberOfEn
 
     //           Compute fractions for multiple zones.
 
-    for (int IZ = 1; IZ <= NumberOfEnclosures; ++IZ) {
+    for (int IZ = 1; IZ <= numEnclosures; ++IZ) {
         if (!state.dataHeatBalSurf->EnclSolRecDifShortFromZ(IZ)) continue;
 
-        for (int JZ = 1; JZ <= NumberOfEnclosures; ++JZ) {
+        for (int JZ = 1; JZ <= numEnclosures; ++JZ) {
             if (!state.dataHeatBalSurf->EnclSolRecDifShortFromZ(JZ)) continue;
             if (IZ == JZ) continue;
             if (state.dataHeatBalSurfMgr->DiffuseArray(IZ, JZ) == 0.0) continue;
 
-            for (int KZ = 1; KZ <= NumberOfEnclosures; ++KZ) {
+            for (int KZ = 1; KZ <= numEnclosures; ++KZ) {
                 if (!state.dataHeatBalSurf->EnclSolRecDifShortFromZ(KZ)) continue;
                 if (IZ == KZ) continue;
                 if (JZ == KZ) continue;
                 if (state.dataHeatBalSurfMgr->DiffuseArray(JZ, KZ) == 0.0) continue;
                 state.dataHeatBalSurf->ZoneFractDifShortZtoZ(IZ, KZ) +=
                     state.dataHeatBalSurfMgr->DiffuseArray(JZ, KZ) * state.dataHeatBalSurfMgr->DiffuseArray(IZ, JZ);
 
-                for (int LZ = 1; LZ <= NumberOfEnclosures; ++LZ) {
+                for (int LZ = 1; LZ <= numEnclosures; ++LZ) {
                     if (!state.dataHeatBalSurf->EnclSolRecDifShortFromZ(LZ)) continue;
                     if (IZ == LZ) continue;
                     if (JZ == LZ) continue;
@@ -4495,7 +4495,7 @@ void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state, int const NumberOfEn
                                                                             state.dataHeatBalSurfMgr->DiffuseArray(JZ, KZ) *
                                                                             state.dataHeatBalSurfMgr->DiffuseArray(IZ, JZ);
 
-                    for (int MZ = 1; MZ <= NumberOfEnclosures; ++MZ) {
+                    for (int MZ = 1; MZ <= numEnclosures; ++MZ) {
                         if (!state.dataHeatBalSurf->EnclSolRecDifShortFromZ(MZ)) continue;
                         if (IZ == MZ) continue;
                         if (JZ == MZ) continue;

src/EnergyPlus/HeatBalanceSurfaceManager.hh

@@ -114,7 +114,7 @@ namespace HeatBalanceSurfaceManager {
 
     void ComputeIntSWAbsorpFactors(EnergyPlusData &state);
 
-    void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state, int NumberOfEnclosures); // Number of solar enclosures
+    void ComputeDifSolExcZonesWIZWindows(EnergyPlusData &state);
 
     void InitEMSControlledSurfaceProperties(EnergyPlusData &state);
 

src/EnergyPlus/SolarReflectionManager.cc

@@ -311,7 +311,7 @@ namespace SolarReflectionManager {
                     ObsVec = state.dataSurface->Surface(ObsSurfNum).Vertex(loop);
                     DotProd = dot(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec, ObsVec - RecVec);
                     // CR8251      IF(DotProd > 0.01d0) THEN  ! This obstructing-surface vertex is not behind receiving surface
-                    if (DotProd > 1.0e-6) { // This obstructing-surface vertex is not behind receiving surface
+                    if (DotProd > Constant::OneMillionth) { // This obstructing-surface vertex is not behind receiving surface
                         ObsBehindRec = false;
                         break;
                     }

src/EnergyPlus/SolarShading.cc

@@ -5746,7 +5746,7 @@ void SHADOW(EnergyPlusData &state,
                      state.dataSurface->Surface(NGRS).lcsz.y * state.dataSolarShading->SUNCOS(2) +
                      state.dataSurface->Surface(NGRS).lcsz.z * state.dataSolarShading->SUNCOS(3);
 
-                if (std::abs(ZS) > 1.e-4) {
+                if (std::abs(ZS) > Constant::SmallDistance) {
                     state.dataSolarShading->XShadowProjection = XS / ZS;
                     state.dataSolarShading->YShadowProjection = YS / ZS;
                     if (std::abs(state.dataSolarShading->XShadowProjection) < 1.e-8) state.dataSolarShading->XShadowProjection = 0.0;
@@ -11231,7 +11231,7 @@ void CalcBeamSolarOnWinRevealSurface(EnergyPlusData &state)
                     if (A2ill < 0.0) A2ill = 0.0;
 
                     // Quantities related to outside reveal
-                    if (A1ill > 1.0e-6) {
+                    if (A1ill > Constant::OneMillionth) {
 
                         state.dataSurface->SurfWinBmSolAbsdOutsReveal(SurfNum) +=
                             A1ill * state.dataSurface->SurfWinOutsideRevealSolAbs(SurfNum) * CosBeta * tmp_SunlitFracWithoutReveal;
@@ -11259,7 +11259,7 @@ void CalcBeamSolarOnWinRevealSurface(EnergyPlusData &state)
 
                     if (NOT_SHADED(ShadeFlag) || ShadeFlag == WinShadingType::SwitchableGlazing) {
 
-                        if (A2ill > 1.0e-6) {
+                        if (A2ill > Constant::OneMillionth) {
 
                             DiffReflGlass = state.dataConstruction->Construct(ConstrNum).ReflectSolDiffBack;
                             if (ShadeFlag == WinShadingType::SwitchableGlazing) {
@@ -12751,7 +12751,7 @@ void CalcComplexWindowOverlap(EnergyPlusData &state,
         ZSp = -SdotZ;
 
         // Projection of shadows for current basis direciton
-        if (std::abs(ZSp) > 1.e-4) {
+        if (std::abs(ZSp) > Constant::SmallDistance) {
             XShadowProjection = XSp / ZSp;
             YShadowProjection = YSp / ZSp;
             if (std::abs(XShadowProjection) < 1.e-8) XShadowProjection = 0.0;