Line data Source code
1 : // EnergyPlus, Copyright (c) 1996-2025, The Board of Trustees of the University of Illinois,
2 : // The Regents of the University of California, through Lawrence Berkeley National Laboratory
3 : // (subject to receipt of any required approvals from the U.S. Dept. of Energy), Oak Ridge
4 : // National Laboratory, managed by UT-Battelle, Alliance for Sustainable Energy, LLC, and other
5 : // contributors. All rights reserved.
6 : //
7 : // NOTICE: This Software was developed under funding from the U.S. Department of Energy and the
8 : // U.S. Government consequently retains certain rights. As such, the U.S. Government has been
9 : // granted for itself and others acting on its behalf a paid-up, nonexclusive, irrevocable,
10 : // worldwide license in the Software to reproduce, distribute copies to the public, prepare
11 : // derivative works, and perform publicly and display publicly, and to permit others to do so.
12 : //
13 : // Redistribution and use in source and binary forms, with or without modification, are permitted
14 : // provided that the following conditions are met:
15 : //
16 : // (1) Redistributions of source code must retain the above copyright notice, this list of
17 : // conditions and the following disclaimer.
18 : //
19 : // (2) Redistributions in binary form must reproduce the above copyright notice, this list of
20 : // conditions and the following disclaimer in the documentation and/or other materials
21 : // provided with the distribution.
22 : //
23 : // (3) Neither the name of the University of California, Lawrence Berkeley National Laboratory,
24 : // the University of Illinois, U.S. Dept. of Energy nor the names of its contributors may be
25 : // used to endorse or promote products derived from this software without specific prior
26 : // written permission.
27 : //
28 : // (4) Use of EnergyPlus(TM) Name. If Licensee (i) distributes the software in stand-alone form
29 : // without changes from the version obtained under this License, or (ii) Licensee makes a
30 : // reference solely to the software portion of its product, Licensee must refer to the
31 : // software as "EnergyPlus version X" software, where "X" is the version number Licensee
32 : // obtained under this License and may not use a different name for the software. Except as
33 : // specifically required in this Section (4), Licensee shall not use in a company name, a
34 : // product name, in advertising, publicity, or other promotional activities any name, trade
35 : // name, trademark, logo, or other designation of "EnergyPlus", "E+", "e+" or confusingly
36 : // similar designation, without the U.S. Department of Energy's prior written consent.
37 : //
38 : // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
39 : // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY
40 : // AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
41 : // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
42 : // CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
43 : // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
44 : // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
45 : // OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
46 : // POSSIBILITY OF SUCH DAMAGE.
47 :
48 : // C++ Headers
49 : #include <cassert>
50 :
51 : // EnergyPlus Headers
52 : #include <EnergyPlus/Data/EnergyPlusData.hh>
53 : #include <EnergyPlus/DataGlobals.hh>
54 : #include <EnergyPlus/TARCOGArgs.hh>
55 : #include <EnergyPlus/TARCOGCommon.hh>
56 : #include <EnergyPlus/TARCOGGasses90.hh>
57 : #include <EnergyPlus/TARCOGGassesParams.hh>
58 : #include <EnergyPlus/TARCOGParams.hh>
59 : #include <EnergyPlus/TarcogShading.hh>
60 : #include <EnergyPlus/ThermalISO15099Calc.hh>
61 :
62 : namespace EnergyPlus::ThermalISO15099Calc {
63 : //***********************************************************************
64 : // ThermalISO15099Calc: a TARCOG module
65 : // module For Calculation of Thermal Performance Indices For Center
66 : // of Glass According to ISO 15099
67 : // History of Revisions:
68 : // Revision: 6.0.36 (June/22/2010)
69 : // - Initial setup, extracted and refactored from TARCOG.for
70 : //***********************************************************************
71 :
72 : // MODULE INFORMATION:
73 : // AUTHOR D. Charlie Curcija
74 : // DATE WRITTEN July/2000
75 : // MODIFIED na
76 : // RE-ENGINEERED March/27/2012, Simon Vidanovic
77 :
78 : // Revision: 7.0.13 (March/27/2012), Simon Vidanovic
79 : // - feature: New set of equations is set instead of hhat coefficients and new approach to solution which improves
80 : // speed and stability. Note that this solution does not include laminates
81 :
82 : // PURPOSE OF THIS MODULE:
83 : // Module For Calculation of Thermal Performance Indices For Center
84 : // of Glass According to ISO 15099
85 :
86 : // Using/Aliasing
87 : using namespace TARCOGGassesParams;
88 : using namespace TARCOGParams;
89 : using namespace TARCOGArgs;
90 : using namespace TARCOGCommon;
91 : using namespace TARCOGOutput;
92 : using namespace TARCOGGasses90;
93 : using namespace TarcogShading;
94 :
95 : // Calculation outcome
96 : enum class CalculationOutcome
97 : {
98 : Invalid = -1,
99 : OK,
100 : Num
101 : };
102 :
103 0 : void film(Real64 const tex, Real64 const tw, Real64 const ws, int const iwd, Real64 &hcout, int const ibc)
104 : {
105 : //***********************************************************************
106 : // purpose - to find outdoor film coeff
107 : //***********************************************************************
108 : // Inputs -
109 : // tex - outdoor air temp [k]
110 : // tw - outside surface temp
111 : // ws - wind speed [m/s]
112 : // iwd - wind direction [0 - windward; 1 - leeward]
113 : // Outputs
114 : // hcout - convective film coeff [w m-2 k-1]
115 :
116 : // Locals
117 0 : Real64 constexpr conv(5.6783);
118 :
119 : Real64 vc;
120 : Real64 acoef;
121 : Real64 bexp;
122 :
123 : // calculation of convection component of exterior film coefficient using the :
124 0 : switch (ibc) {
125 0 : case 0: { // ISO 15099
126 0 : hcout = 4.0 + 4.0 * ws;
127 0 : } break;
128 0 : case -1: { // old ASHRAE SPC142 correlation
129 0 : if (iwd == 0) { // windward
130 0 : if (ws > 2.0) {
131 0 : vc = 0.25 * ws;
132 : } else {
133 0 : vc = 0.5;
134 : }
135 : } else { // leeward
136 0 : vc = 0.3 + 0.05 * ws;
137 : }
138 0 : hcout = 3.28 * std::pow(vc, 0.605);
139 0 : hcout *= conv; // convert to metric
140 0 : } break;
141 0 : case -2: { // Yazdanian-Klems correlation:
142 0 : if (iwd == 0) { // windward
143 0 : acoef = 2.38;
144 0 : bexp = 0.89;
145 : } else { // leeward
146 0 : acoef = 2.86;
147 0 : bexp = 0.617;
148 : }
149 0 : hcout = std::sqrt(pow_2(0.84 * std::pow(tw - tex, 0.33)) + pow_2(acoef * std::pow(ws, bexp)));
150 0 : } break;
151 0 : case -3: { // Kimura correlation (Section 8.4.2.3 in ISO 15099-2001):
152 0 : if (iwd == 0) { // windward
153 0 : if (ws > 2.0) {
154 0 : vc = 0.25 * ws;
155 : } else {
156 0 : vc = 0.5 * ws;
157 : }
158 : } else { // leeward
159 0 : vc = 0.3 + 0.05 * ws;
160 : }
161 0 : hcout = 4.7 + 7.6 * vc;
162 0 : } break;
163 0 : default:
164 0 : break;
165 : }
166 0 : }
167 :
168 0 : void Calc_ISO15099(EnergyPlusData &state,
169 : Files &files,
170 : int const nlayer,
171 : int const iwd,
172 : Real64 &tout,
173 : Real64 &tind,
174 : Real64 &trmin,
175 : Real64 const wso,
176 : Real64 const wsi,
177 : Real64 const dir,
178 : Real64 const outir,
179 : int const isky,
180 : Real64 const tsky,
181 : Real64 &esky,
182 : Real64 const fclr,
183 : Real64 const VacuumPressure,
184 : Real64 const VacuumMaxGapThickness,
185 : Array1D<Real64> &gap,
186 : Array1D<Real64> &thick,
187 : Array1D<Real64> &scon,
188 : const Array1D<Real64> &tir,
189 : const Array1D<Real64> &emis,
190 : Real64 const totsol,
191 : Real64 const tilt,
192 : const Array1D<Real64> &asol,
193 : Real64 const height,
194 : Real64 const heightt,
195 : Real64 const width,
196 : const Array1D<Real64> &presure,
197 : Array2A_int const iprop,
198 : Array2A<Real64> const frct,
199 : Array2A<Real64> const xgcon,
200 : Array2A<Real64> const xgvis,
201 : Array2A<Real64> const xgcp,
202 : const Array1D<Real64> &xwght,
203 : const Array1D<Real64> &gama,
204 : const Array1D_int &nmix,
205 : const Array1D_int &SupportPillar,
206 : const Array1D<Real64> &PillarSpacing,
207 : const Array1D<Real64> &PillarRadius,
208 : Array1D<Real64> &theta,
209 : Array1D<Real64> &q,
210 : Array1D<Real64> &qv,
211 : Real64 &ufactor,
212 : Real64 &sc,
213 : Real64 &hflux,
214 : Real64 &hcin,
215 : Real64 &hcout,
216 : Real64 &hrin,
217 : Real64 &hrout,
218 : Real64 &hin,
219 : Real64 &hout,
220 : Array1D<Real64> &hcgas,
221 : Array1D<Real64> &hrgas,
222 : Real64 &shgc,
223 : int &nperr,
224 : std::string &ErrorMessage,
225 : Real64 &shgct,
226 : Real64 &tamb,
227 : Real64 &troom,
228 : const Array1D_int &ibc,
229 : const Array1D<Real64> &Atop,
230 : const Array1D<Real64> &Abot,
231 : const Array1D<Real64> &Al,
232 : const Array1D<Real64> &Ar,
233 : const Array1D<Real64> &Ah,
234 : const Array1D<Real64> &SlatThick,
235 : const Array1D<Real64> &SlatWidth,
236 : const Array1D<Real64> &SlatAngle,
237 : const Array1D<Real64> &SlatCond,
238 : const Array1D<Real64> &SlatSpacing,
239 : const Array1D<Real64> &SlatCurve,
240 : const Array1D<Real64> &vvent,
241 : const Array1D<Real64> &tvent,
242 : const Array1D<TARCOGLayerType> &LayerType,
243 : const Array1D_int &nslice,
244 : const Array1D<Real64> &LaminateA,
245 : const Array1D<Real64> &LaminateB,
246 : const Array1D<Real64> &sumsol,
247 : Array1D<Real64> &Ra,
248 : Array1D<Real64> &Nu,
249 : TARCOGThermalModel const ThermalMod,
250 : int const Debug_mode,
251 : Real64 &ShadeEmisRatioOut,
252 : Real64 &ShadeEmisRatioIn,
253 : Real64 &ShadeHcRatioOut,
254 : Real64 &ShadeHcRatioIn,
255 : Real64 &HcUnshadedOut,
256 : Real64 &HcUnshadedIn,
257 : Array1D<Real64> &Keff,
258 : Array1D<Real64> &ShadeGapKeffConv,
259 : Real64 const SDScalar,
260 : int const SHGCCalc,
261 : int &NumOfIterations,
262 : Real64 const edgeGlCorrFac)
263 : {
264 :
265 : // Argument array dimensioning
266 0 : EP_SIZE_CHECK(gap, MaxGap);
267 0 : EP_SIZE_CHECK(thick, maxlay);
268 0 : EP_SIZE_CHECK(scon, maxlay);
269 0 : EP_SIZE_CHECK(tir, maxlay2);
270 0 : EP_SIZE_CHECK(emis, maxlay2);
271 0 : EP_SIZE_CHECK(asol, maxlay);
272 0 : EP_SIZE_CHECK(presure, maxlay1);
273 0 : iprop.dim(maxgas, maxlay1);
274 0 : frct.dim(maxgas, maxlay1);
275 0 : xgcon.dim(3, maxgas);
276 0 : xgvis.dim(3, maxgas);
277 0 : xgcp.dim(3, maxgas);
278 0 : EP_SIZE_CHECK(xwght, maxgas);
279 0 : EP_SIZE_CHECK(gama, maxgas);
280 0 : EP_SIZE_CHECK(nmix, maxlay1);
281 0 : EP_SIZE_CHECK(SupportPillar, maxlay);
282 0 : EP_SIZE_CHECK(PillarSpacing, maxlay);
283 0 : EP_SIZE_CHECK(PillarRadius, maxlay);
284 0 : EP_SIZE_CHECK(theta, maxlay2);
285 0 : EP_SIZE_CHECK(q, maxlay3);
286 0 : EP_SIZE_CHECK(qv, maxlay1);
287 0 : EP_SIZE_CHECK(hcgas, maxlay1);
288 0 : EP_SIZE_CHECK(hrgas, maxlay1);
289 0 : EP_SIZE_CHECK(ibc, 2);
290 0 : EP_SIZE_CHECK(Atop, maxlay);
291 0 : EP_SIZE_CHECK(Abot, maxlay);
292 0 : EP_SIZE_CHECK(Al, maxlay);
293 0 : EP_SIZE_CHECK(Ar, maxlay);
294 0 : EP_SIZE_CHECK(Ah, maxlay);
295 0 : EP_SIZE_CHECK(SlatThick, maxlay);
296 0 : EP_SIZE_CHECK(SlatWidth, maxlay);
297 0 : EP_SIZE_CHECK(SlatAngle, maxlay);
298 0 : EP_SIZE_CHECK(SlatCond, maxlay);
299 0 : EP_SIZE_CHECK(SlatSpacing, maxlay);
300 0 : EP_SIZE_CHECK(SlatCurve, maxlay);
301 0 : EP_SIZE_CHECK(vvent, maxlay1);
302 0 : EP_SIZE_CHECK(tvent, maxlay1);
303 0 : EP_SIZE_CHECK(LayerType, maxlay);
304 0 : EP_SIZE_CHECK(nslice, maxlay);
305 0 : EP_SIZE_CHECK(LaminateA, maxlay);
306 0 : EP_SIZE_CHECK(LaminateB, maxlay);
307 0 : EP_SIZE_CHECK(sumsol, maxlay);
308 0 : EP_SIZE_CHECK(Ra, maxlay);
309 0 : EP_SIZE_CHECK(Nu, maxlay);
310 0 : EP_SIZE_CHECK(Keff, maxlay);
311 0 : EP_SIZE_CHECK(ShadeGapKeffConv, MaxGap);
312 :
313 : // REAL(r64) :: grho(maxgas,3)
314 : Real64 shgct_NOSD;
315 : Real64 trmout;
316 :
317 : Real64 Gout;
318 : Real64 Gin;
319 : Real64 AchievedErrorTolerance;
320 : Real64 AchievedErrorToleranceSolar;
321 : int NumOfIter;
322 : int NumOfIterSolar;
323 :
324 : Real64 tgg;
325 : Real64 qc1;
326 : Real64 qc2;
327 : Real64 qcgg;
328 :
329 : Real64 ShadeHcModifiedOut;
330 : Real64 ShadeHcModifiedIn;
331 :
332 : // cbi...Variables for "unshaded" run:
333 :
334 : bool NeedUnshadedRun;
335 : int nlayer_NOSD;
336 : Real64 AchievedErrorTolerance_NOSD;
337 : int NumOfIter_NOSD;
338 : Real64 hin_NOSD;
339 : Real64 flux_NOSD;
340 : Real64 hcin_NOSD;
341 : Real64 hrin_NOSD;
342 : Real64 hcout_NOSD;
343 : Real64 hrout_NOSD;
344 : Real64 tamb_NOSD;
345 : Real64 troom_NOSD;
346 : Real64 ufactor_NOSD;
347 : Real64 sc_NOSD;
348 : Real64 hflux_NOSD;
349 : Real64 shgc_NOSD;
350 : Real64 hout_NOSD;
351 : // REAL(r64) :: rs_NOSD(maxlay3)!,sol(maxlay)
352 : Real64 ShadeEmisRatioOut_NOSD;
353 : Real64 ShadeEmisRatioIn_NOSD;
354 : Real64 ShadeHcRatioOut_NOSD;
355 : Real64 ShadeHcRatioIn_NOSD;
356 : Real64 ShadeHcModifiedOut_NOSD;
357 : Real64 ShadeHcModifiedIn_NOSD;
358 :
359 : int FirstSpecularLayer;
360 : int LastSpecularLayer;
361 :
362 : // cbi...Other variables:
363 : Real64 flux;
364 : Real64 hint;
365 : Real64 houtt;
366 : Real64 ebsky;
367 : Real64 ebroom;
368 : int i;
369 : int j;
370 : int OriginalIndex;
371 : int UnshadedDebug;
372 :
373 : // Autodesk:Uninit Initialize variables used uninitialized
374 0 : shgc_NOSD = 0.0; // Autodesk:Uninit Force default initialization
375 0 : sc_NOSD = 0.0; // Autodesk:Uninit Force default initialization
376 0 : hflux_NOSD = 0.0; // Autodesk:Uninit Force default initialization
377 0 : ShadeHcRatioIn_NOSD = 0.0; // Autodesk:Uninit Force default initialization
378 0 : ShadeHcRatioOut_NOSD = 0.0; // Autodesk:Uninit Force default initialization
379 :
380 0 : AchievedErrorTolerance = 0.0;
381 0 : AchievedErrorToleranceSolar = 0.0;
382 0 : AchievedErrorTolerance_NOSD = 0.0;
383 :
384 0 : PrepVariablesISO15099(nlayer,
385 : tout,
386 : tind,
387 : trmin,
388 : isky,
389 : outir,
390 : tsky,
391 : esky,
392 : fclr,
393 : gap,
394 : thick,
395 : scon,
396 : tir,
397 : emis,
398 : tilt,
399 : hin,
400 : hout,
401 : ibc,
402 : SlatThick,
403 : SlatWidth,
404 : SlatAngle,
405 : SlatCond,
406 : LayerType,
407 : ThermalMod,
408 : SDScalar,
409 : ShadeEmisRatioOut,
410 : ShadeEmisRatioIn,
411 : ShadeHcRatioOut,
412 : ShadeHcRatioIn,
413 : Keff,
414 : ShadeGapKeffConv,
415 : sc,
416 : shgc,
417 : ufactor,
418 : flux,
419 0 : state.dataThermalISO15099Calc->LaminateAU,
420 0 : state.dataThermalISO15099Calc->sumsolU,
421 0 : state.dataThermalISO15099Calc->sol0,
422 : hint,
423 : houtt,
424 : trmout,
425 : ebsky,
426 : ebroom,
427 : Gout,
428 : Gin,
429 0 : state.dataThermalISO15099Calc->rir,
430 0 : state.dataThermalISO15099Calc->vfreevent,
431 : nperr,
432 : ErrorMessage);
433 :
434 0 : for (int i = 1; i <= nlayer; ++i) {
435 0 : state.dataThermalISO15099Calc->EffectiveOpenness(i) = Ah(i) / (width * height);
436 : }
437 :
438 0 : updateEffectiveMultipliers(nlayer,
439 : width,
440 : height,
441 : Atop,
442 : Abot,
443 : Al,
444 : Ar,
445 : Ah,
446 0 : state.dataThermalISO15099Calc->Atop_eff,
447 0 : state.dataThermalISO15099Calc->Abot_eff,
448 0 : state.dataThermalISO15099Calc->Al_eff,
449 0 : state.dataThermalISO15099Calc->Ar_eff,
450 0 : state.dataThermalISO15099Calc->Ah_eff,
451 : LayerType,
452 : SlatAngle);
453 :
454 : // No option to take hardcoded variables. All gas coefficients are now passed from outside.
455 : // if (GoAhead(nperr)) call propcon90(ISO15099,mgas,xgcon,xgvis,xgcp,xgrho,xwght,nperr)
456 :
457 : // exit on error
458 0 : if (!(GoAhead(nperr))) {
459 0 : return;
460 : }
461 :
462 : // bi...Write intermediate results to output file:
463 0 : if (files.WriteDebugOutput) {
464 0 : WriteModifiedArguments(files.DebugOutputFile,
465 0 : files.DBGD,
466 : esky,
467 : trmout,
468 : trmin,
469 : ebsky,
470 : ebroom,
471 : Gout,
472 : Gin,
473 : nlayer,
474 : LayerType,
475 : nmix,
476 : frct,
477 : thick,
478 : scon,
479 : gap,
480 : xgcon,
481 : xgvis,
482 : xgcp,
483 : xwght);
484 : }
485 :
486 : // cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
487 : // This is "solar radiation" pass
488 : // cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
489 :
490 : // This is main calculation in case UFactor calculations will not be performed
491 0 : if ((dir > 0.0) || (SHGCCalc == 0)) {
492 : // call therm1d to calculate heat flux with solar radiation
493 :
494 0 : therm1d(state,
495 : files,
496 : nlayer,
497 : iwd,
498 : tout,
499 : tind,
500 : wso,
501 : wsi,
502 : VacuumPressure,
503 : VacuumMaxGapThickness,
504 : dir,
505 : ebsky,
506 : Gout,
507 : trmout,
508 : trmin,
509 : ebroom,
510 : Gin,
511 : tir,
512 0 : state.dataThermalISO15099Calc->rir,
513 : emis,
514 : gap,
515 : thick,
516 : scon,
517 : tilt,
518 : asol,
519 : height,
520 : heightt,
521 : width,
522 : iprop,
523 : frct,
524 : presure,
525 : nmix,
526 : xwght,
527 : xgcon,
528 : xgvis,
529 : xgcp,
530 : gama,
531 : SupportPillar,
532 : PillarSpacing,
533 : PillarRadius,
534 : theta,
535 : q,
536 : qv,
537 : flux,
538 : hcin,
539 : hrin,
540 : hcout,
541 : hrout,
542 : hin,
543 : hout,
544 : hcgas,
545 : hrgas,
546 : ufactor,
547 : nperr,
548 : ErrorMessage,
549 : tamb,
550 : troom,
551 : ibc,
552 0 : state.dataThermalISO15099Calc->Atop_eff,
553 0 : state.dataThermalISO15099Calc->Abot_eff,
554 0 : state.dataThermalISO15099Calc->Al_eff,
555 0 : state.dataThermalISO15099Calc->Ar_eff,
556 0 : state.dataThermalISO15099Calc->Ah_eff,
557 0 : state.dataThermalISO15099Calc->EffectiveOpenness,
558 : vvent,
559 : tvent,
560 : LayerType,
561 : Ra,
562 : Nu,
563 0 : state.dataThermalISO15099Calc->vfreevent,
564 0 : state.dataThermalISO15099Calc->qcgas,
565 0 : state.dataThermalISO15099Calc->qrgas,
566 0 : state.dataThermalISO15099Calc->Ebf,
567 0 : state.dataThermalISO15099Calc->Ebb,
568 0 : state.dataThermalISO15099Calc->Rf,
569 0 : state.dataThermalISO15099Calc->Rb,
570 : ShadeEmisRatioOut,
571 : ShadeEmisRatioIn,
572 : ShadeHcModifiedOut,
573 : ShadeHcModifiedIn,
574 : ThermalMod,
575 : Debug_mode,
576 : AchievedErrorToleranceSolar,
577 : NumOfIterSolar,
578 : edgeGlCorrFac);
579 :
580 0 : NumOfIterations = NumOfIterSolar;
581 : // exit on error:
582 :
583 0 : if (nlayer > 1) {
584 0 : for (i = 1; i <= nlayer - 1; ++i) {
585 0 : Keff(i) = gap(i) * q(2 * i + 1) / (theta(2 * i + 1) - theta(2 * i));
586 0 : if (IsShadingLayer(LayerType(i))) {
587 0 : Keff(i) = gap(i) * q(2 * i + 1) / (theta(2 * i + 1) - theta(2 * i));
588 : }
589 0 : if (IsShadingLayer(LayerType(i + 1))) {
590 0 : Keff(i) = gap(i) * q(2 * i + 1) / (theta(2 * i + 1) - theta(2 * i));
591 : }
592 : }
593 : }
594 :
595 0 : if (!(GoAhead(nperr))) {
596 0 : return;
597 : }
598 :
599 : // No need to store results in case of non-ufactor run
600 0 : if ((SHGCCalc > 0) && (dir > 0.0)) {
601 0 : solarISO15099(
602 0 : totsol, state.dataThermalISO15099Calc->rtot, state.dataThermalISO15099Calc->rs, nlayer, asol, state.dataThermalISO15099Calc->sft);
603 0 : shgct = state.dataThermalISO15099Calc->sft;
604 0 : shgct_NOSD = 0.0;
605 0 : state.dataThermalISO15099Calc->hcins = hcin;
606 0 : state.dataThermalISO15099Calc->hrins = hrin;
607 0 : state.dataThermalISO15099Calc->hins = hin;
608 0 : state.dataThermalISO15099Calc->hcouts = hcout;
609 0 : state.dataThermalISO15099Calc->hrouts = hrout;
610 0 : state.dataThermalISO15099Calc->houts = hout;
611 0 : state.dataThermalISO15099Calc->ufactors = ufactor;
612 0 : state.dataThermalISO15099Calc->fluxs = flux;
613 0 : for (i = 1; i <= nlayer; ++i) {
614 0 : state.dataThermalISO15099Calc->thetas(2 * i - 1) = theta(2 * i - 1);
615 0 : state.dataThermalISO15099Calc->thetas(2 * i) = theta(2 * i);
616 0 : state.dataThermalISO15099Calc->Ebbs(i) = state.dataThermalISO15099Calc->Ebb(i);
617 0 : state.dataThermalISO15099Calc->Ebfs(i) = state.dataThermalISO15099Calc->Ebf(i);
618 0 : state.dataThermalISO15099Calc->Rbs(i) = state.dataThermalISO15099Calc->Rb(i);
619 0 : state.dataThermalISO15099Calc->Rfs(i) = state.dataThermalISO15099Calc->Rf(i);
620 0 : state.dataThermalISO15099Calc->qs(2 * i - 1) = q(2 * i - 1);
621 0 : state.dataThermalISO15099Calc->qs(2 * i) = q(2 * i);
622 : // qprims(2*i - 1) = qprim(2*i - 1)
623 : // qprims(2*i) = qprim(2*i)
624 0 : state.dataThermalISO15099Calc->qvs(2 * i - 1) = qv(2 * i - 1);
625 0 : state.dataThermalISO15099Calc->qvs(2 * i) = qv(2 * i);
626 0 : state.dataThermalISO15099Calc->hcgass(i) = hcgas(i);
627 0 : state.dataThermalISO15099Calc->hrgass(i) = hrgas(i);
628 0 : state.dataThermalISO15099Calc->qrgaps(i) = state.dataThermalISO15099Calc->qrgas(i);
629 0 : state.dataThermalISO15099Calc->qcgaps(i) = state.dataThermalISO15099Calc->qcgas(i);
630 : }
631 : // CHECK THIS!
632 0 : state.dataThermalISO15099Calc->qs(2 * nlayer + 1) = q(2 * nlayer + 1);
633 : } // if (UFactorCalc.gt.0) then
634 : }
635 :
636 : // No solar radiation pass is not needed to be calculated
637 : // if ((SHGCCalc.gt.0).or.(dir.eq.0)) then
638 0 : if (SHGCCalc > 0) {
639 :
640 : // cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
641 : // This is "no solar radiation" pass
642 : // cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
643 :
644 0 : hin = hint;
645 0 : hout = houtt;
646 :
647 : // call therm1d to calculate heat flux without solar radiation
648 0 : therm1d(state,
649 : files,
650 : nlayer,
651 : iwd,
652 : tout,
653 : tind,
654 : wso,
655 : wsi,
656 : VacuumPressure,
657 : VacuumMaxGapThickness,
658 : 0.0,
659 : ebsky,
660 : Gout,
661 : trmout,
662 : trmin,
663 : ebroom,
664 : Gin,
665 : tir,
666 0 : state.dataThermalISO15099Calc->rir,
667 : emis,
668 : gap,
669 : thick,
670 : scon,
671 : tilt,
672 0 : state.dataThermalISO15099Calc->sol0,
673 : height,
674 : heightt,
675 : width,
676 : iprop,
677 : frct,
678 : presure,
679 : nmix,
680 : xwght,
681 : xgcon,
682 : xgvis,
683 : xgcp,
684 : gama,
685 : SupportPillar,
686 : PillarSpacing,
687 : PillarRadius,
688 : theta,
689 : q,
690 : qv,
691 : flux,
692 : hcin,
693 : hrin,
694 : hcout,
695 : hrout,
696 : hin,
697 : hout,
698 : hcgas,
699 : hrgas,
700 : ufactor,
701 : nperr,
702 : ErrorMessage,
703 : tamb,
704 : troom,
705 : ibc,
706 0 : state.dataThermalISO15099Calc->Atop_eff,
707 0 : state.dataThermalISO15099Calc->Abot_eff,
708 0 : state.dataThermalISO15099Calc->Al_eff,
709 0 : state.dataThermalISO15099Calc->Ar_eff,
710 0 : state.dataThermalISO15099Calc->Ah_eff,
711 0 : state.dataThermalISO15099Calc->EffectiveOpenness,
712 : vvent,
713 : tvent,
714 : LayerType,
715 : Ra,
716 : Nu,
717 0 : state.dataThermalISO15099Calc->vfreevent,
718 0 : state.dataThermalISO15099Calc->qcgas,
719 0 : state.dataThermalISO15099Calc->qrgas,
720 0 : state.dataThermalISO15099Calc->Ebf,
721 0 : state.dataThermalISO15099Calc->Ebb,
722 0 : state.dataThermalISO15099Calc->Rf,
723 0 : state.dataThermalISO15099Calc->Rb,
724 : ShadeEmisRatioOut,
725 : ShadeEmisRatioIn,
726 : ShadeHcModifiedOut,
727 : ShadeHcModifiedIn,
728 : ThermalMod,
729 : Debug_mode,
730 : AchievedErrorTolerance,
731 : NumOfIter,
732 : edgeGlCorrFac);
733 :
734 0 : NumOfIterations = NumOfIter;
735 :
736 : // exit on error:
737 0 : if (!(GoAhead(nperr))) {
738 0 : return;
739 : }
740 :
741 : // bi...Keep hcout, hcin in case this is an unshaded system:
742 0 : HcUnshadedOut = hcout;
743 0 : HcUnshadedIn = hcin;
744 :
745 : // bi...do an Unshaded run if necessary (Uvalue/Winter conditions):
746 : // bi...Prepare variables for UNSHADED (NO SD) run:
747 :
748 0 : NeedUnshadedRun = false;
749 0 : FirstSpecularLayer = 1;
750 0 : LastSpecularLayer = nlayer;
751 0 : nlayer_NOSD = nlayer;
752 0 : if (IsShadingLayer(LayerType(1))) {
753 0 : --nlayer_NOSD;
754 0 : FirstSpecularLayer = 2;
755 0 : NeedUnshadedRun = true;
756 : }
757 :
758 : // if (LayerType(nlayer).eq.VENETBLIND) then
759 0 : if (IsShadingLayer(LayerType(nlayer))) {
760 0 : --nlayer_NOSD;
761 0 : LastSpecularLayer = nlayer - 1;
762 0 : NeedUnshadedRun = true;
763 : }
764 :
765 : // no unshaded run for now
766 0 : NeedUnshadedRun = false;
767 : // bi...Set outdoor & indoor gas properties:
768 0 : if (NeedUnshadedRun) {
769 0 : state.dataThermalISO15099Calc->nmix_NOSD(1) = nmix(1);
770 0 : state.dataThermalISO15099Calc->presure_NOSD(1) = presure(1);
771 0 : state.dataThermalISO15099Calc->nmix_NOSD(nlayer_NOSD + 1) = nmix(nlayer + 1);
772 0 : state.dataThermalISO15099Calc->presure_NOSD(nlayer_NOSD + 1) = presure(nlayer + 1);
773 0 : for (j = 1; j <= nmix(1); ++j) {
774 0 : state.dataThermalISO15099Calc->iprop_NOSD(j, 1) = iprop(j, 1);
775 0 : state.dataThermalISO15099Calc->frct_NOSD(j, 1) = frct(j, 1);
776 : }
777 0 : for (j = 1; j <= nmix(nlayer_NOSD + 1); ++j) {
778 0 : state.dataThermalISO15099Calc->iprop_NOSD(j, nlayer_NOSD + 1) = iprop(j, nlayer + 1);
779 0 : state.dataThermalISO15099Calc->frct_NOSD(j, nlayer_NOSD + 1) = frct(j, nlayer + 1);
780 : }
781 0 : for (i = 1; i <= nlayer_NOSD; ++i) {
782 0 : OriginalIndex = FirstSpecularLayer + i - 1;
783 0 : state.dataThermalISO15099Calc->Atop_NOSD(i) = state.dataThermalISO15099Calc->Atop_eff(OriginalIndex);
784 0 : state.dataThermalISO15099Calc->Abot_NOSD(i) = state.dataThermalISO15099Calc->Abot_eff(OriginalIndex);
785 0 : state.dataThermalISO15099Calc->Al_NOSD(i) = state.dataThermalISO15099Calc->Al_eff(OriginalIndex);
786 0 : state.dataThermalISO15099Calc->Ar_NOSD(i) = state.dataThermalISO15099Calc->Ar_eff(OriginalIndex);
787 0 : state.dataThermalISO15099Calc->Ah_NOSD(i) = state.dataThermalISO15099Calc->Ah_eff(OriginalIndex);
788 :
789 0 : state.dataThermalISO15099Calc->SlatThick_NOSD(i) = SlatThick(OriginalIndex);
790 0 : state.dataThermalISO15099Calc->SlatWidth_NOSD(i) = SlatWidth(OriginalIndex);
791 0 : state.dataThermalISO15099Calc->SlatAngle_NOSD(i) = SlatAngle(OriginalIndex);
792 0 : state.dataThermalISO15099Calc->SlatCond_NOSD(i) = SlatCond(OriginalIndex);
793 0 : state.dataThermalISO15099Calc->SlatSpacing_NOSD(i) = SlatSpacing(OriginalIndex);
794 0 : state.dataThermalISO15099Calc->SlatCurve_NOSD(i) = SlatCurve(OriginalIndex);
795 :
796 : // cbi... TO do when Forced Ventilation is implemented: take care of appropriate arguments!!!
797 : // vvent_NOSD
798 : // tvent_NOSD
799 :
800 0 : state.dataThermalISO15099Calc->LayerType_NOSD(i) = LayerType(OriginalIndex);
801 :
802 0 : state.dataThermalISO15099Calc->thick_NOSD(i) = thick(OriginalIndex);
803 0 : state.dataThermalISO15099Calc->scon_NOSD(i) = scon(OriginalIndex);
804 0 : state.dataThermalISO15099Calc->tir_NOSD(2 * i - 1) = tir(2 * OriginalIndex - 1);
805 0 : state.dataThermalISO15099Calc->emis_NOSD(2 * i - 1) = emis(2 * OriginalIndex - 1);
806 0 : state.dataThermalISO15099Calc->emis_NOSD(2 * i) = emis(2 * OriginalIndex);
807 0 : state.dataThermalISO15099Calc->rir_NOSD(2 * i - 1) = state.dataThermalISO15099Calc->rir(2 * OriginalIndex - 1);
808 0 : state.dataThermalISO15099Calc->rir_NOSD(2 * i) = state.dataThermalISO15099Calc->rir(2 * OriginalIndex);
809 :
810 0 : state.dataThermalISO15099Calc->gap_NOSD(i) = gap(OriginalIndex);
811 :
812 0 : if (i < nlayer_NOSD) {
813 0 : state.dataThermalISO15099Calc->nmix_NOSD(i + 1) = nmix(OriginalIndex + 1);
814 0 : state.dataThermalISO15099Calc->presure_NOSD(i + 1) = presure(OriginalIndex + 1);
815 0 : for (j = 1; j <= state.dataThermalISO15099Calc->nmix_NOSD(i + 1); ++j) {
816 0 : state.dataThermalISO15099Calc->iprop_NOSD(j, i + 1) = iprop(j, OriginalIndex + 1);
817 0 : state.dataThermalISO15099Calc->frct_NOSD(j, i + 1) = frct(j, OriginalIndex + 1);
818 : }
819 : }
820 :
821 0 : state.dataThermalISO15099Calc->LaminateA_NOSD(i) = LaminateA(OriginalIndex);
822 0 : state.dataThermalISO15099Calc->LaminateB_NOSD(i) = LaminateB(OriginalIndex);
823 0 : state.dataThermalISO15099Calc->sumsol_NOSD(i) = sumsol(OriginalIndex);
824 :
825 0 : state.dataThermalISO15099Calc->nslice_NOSD(i) = nslice(OriginalIndex);
826 : }
827 :
828 : // This is UNSHADED pass - no solar radiation:
829 0 : hin_NOSD = hint;
830 0 : hout_NOSD = houtt;
831 :
832 : // Simon: Removed unshaded debug output for now
833 0 : UnshadedDebug = 0;
834 0 : if (files.WriteDebugOutput && (UnshadedDebug == 1)) {
835 0 : print(files.DebugOutputFile, "\n");
836 0 : print(files.DebugOutputFile, "UNSHADED RUN:\n");
837 0 : print(files.DebugOutputFile, "\n");
838 :
839 0 : WriteInputArguments(state,
840 0 : files.DebugOutputFile,
841 0 : files.DBGD,
842 : tout,
843 : tind,
844 : trmin,
845 : wso,
846 : iwd,
847 : wsi,
848 : dir,
849 : outir,
850 : isky,
851 : tsky,
852 : esky,
853 : fclr,
854 : VacuumPressure,
855 : VacuumMaxGapThickness,
856 : ibc,
857 : hout_NOSD,
858 : hin_NOSD,
859 : TARCOGGassesParams::Stdrd::ISO15099,
860 : ThermalMod,
861 : SDScalar,
862 : height,
863 : heightt,
864 : width,
865 : tilt,
866 : totsol,
867 : nlayer_NOSD,
868 0 : state.dataThermalISO15099Calc->LayerType_NOSD,
869 0 : state.dataThermalISO15099Calc->thick_NOSD,
870 0 : state.dataThermalISO15099Calc->scon_NOSD,
871 : asol,
872 0 : state.dataThermalISO15099Calc->tir_NOSD,
873 0 : state.dataThermalISO15099Calc->emis_NOSD,
874 0 : state.dataThermalISO15099Calc->Atop_NOSD,
875 0 : state.dataThermalISO15099Calc->Abot_NOSD,
876 0 : state.dataThermalISO15099Calc->Al_NOSD,
877 0 : state.dataThermalISO15099Calc->Ar_NOSD,
878 0 : state.dataThermalISO15099Calc->Ah_NOSD,
879 0 : state.dataThermalISO15099Calc->SlatThick_NOSD,
880 0 : state.dataThermalISO15099Calc->SlatWidth_NOSD,
881 0 : state.dataThermalISO15099Calc->SlatAngle_NOSD,
882 0 : state.dataThermalISO15099Calc->SlatCond_NOSD,
883 0 : state.dataThermalISO15099Calc->SlatSpacing_NOSD,
884 0 : state.dataThermalISO15099Calc->SlatCurve_NOSD,
885 0 : state.dataThermalISO15099Calc->nslice_NOSD,
886 0 : state.dataThermalISO15099Calc->LaminateA_NOSD,
887 0 : state.dataThermalISO15099Calc->LaminateB_NOSD,
888 0 : state.dataThermalISO15099Calc->sumsol_NOSD,
889 0 : state.dataThermalISO15099Calc->gap_NOSD,
890 0 : state.dataThermalISO15099Calc->vvent_NOSD,
891 0 : state.dataThermalISO15099Calc->tvent_NOSD,
892 0 : state.dataThermalISO15099Calc->presure_NOSD,
893 0 : state.dataThermalISO15099Calc->nmix_NOSD,
894 0 : state.dataThermalISO15099Calc->iprop_NOSD,
895 0 : state.dataThermalISO15099Calc->frct_NOSD,
896 : xgcon,
897 : xgvis,
898 : xgcp,
899 : xwght);
900 :
901 : } // end if UnshadedDebug = 1
902 :
903 : // cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
904 : // This is "Unshaded, No solar radiation" pass
905 : // cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
906 : // call therm1d to calculate heat flux with solar radiation
907 0 : therm1d(state,
908 : files,
909 : nlayer_NOSD,
910 : iwd,
911 : tout,
912 : tind,
913 : wso,
914 : wsi,
915 : VacuumPressure,
916 : VacuumMaxGapThickness,
917 : 0.0,
918 : ebsky,
919 : Gout,
920 : trmout,
921 : trmin,
922 : ebroom,
923 : Gin,
924 0 : state.dataThermalISO15099Calc->tir_NOSD,
925 0 : state.dataThermalISO15099Calc->rir_NOSD,
926 0 : state.dataThermalISO15099Calc->emis_NOSD,
927 0 : state.dataThermalISO15099Calc->gap_NOSD,
928 0 : state.dataThermalISO15099Calc->thick_NOSD,
929 0 : state.dataThermalISO15099Calc->scon_NOSD,
930 : tilt,
931 0 : state.dataThermalISO15099Calc->sol0,
932 : height,
933 : heightt,
934 : width,
935 0 : state.dataThermalISO15099Calc->iprop_NOSD,
936 0 : state.dataThermalISO15099Calc->frct_NOSD,
937 0 : state.dataThermalISO15099Calc->presure_NOSD,
938 0 : state.dataThermalISO15099Calc->nmix_NOSD,
939 : xwght,
940 : xgcon,
941 : xgvis,
942 : xgcp,
943 : gama,
944 : SupportPillar,
945 : PillarSpacing,
946 : PillarRadius,
947 0 : state.dataThermalISO15099Calc->theta_NOSD,
948 0 : state.dataThermalISO15099Calc->q_NOSD,
949 0 : state.dataThermalISO15099Calc->qv_NOSD,
950 : flux_NOSD,
951 : hcin_NOSD,
952 : hrin_NOSD,
953 : hcout_NOSD,
954 : hrout_NOSD,
955 : hin_NOSD,
956 : hout_NOSD,
957 0 : state.dataThermalISO15099Calc->hcgas_NOSD,
958 0 : state.dataThermalISO15099Calc->hrgas_NOSD,
959 : ufactor_NOSD,
960 : nperr,
961 : ErrorMessage,
962 : tamb_NOSD,
963 : troom_NOSD,
964 : ibc,
965 0 : state.dataThermalISO15099Calc->Atop_NOSD,
966 0 : state.dataThermalISO15099Calc->Abot_NOSD,
967 0 : state.dataThermalISO15099Calc->Al_NOSD,
968 0 : state.dataThermalISO15099Calc->Ar_NOSD,
969 0 : state.dataThermalISO15099Calc->Ah_NOSD,
970 0 : state.dataThermalISO15099Calc->EffectiveOpenness_NOSD,
971 0 : state.dataThermalISO15099Calc->vvent_NOSD,
972 0 : state.dataThermalISO15099Calc->tvent_NOSD,
973 0 : state.dataThermalISO15099Calc->LayerType_NOSD,
974 0 : state.dataThermalISO15099Calc->Ra_NOSD,
975 0 : state.dataThermalISO15099Calc->Nu_NOSD,
976 0 : state.dataThermalISO15099Calc->vfreevent_NOSD,
977 0 : state.dataThermalISO15099Calc->qcgas_NOSD,
978 0 : state.dataThermalISO15099Calc->qrgas_NOSD,
979 0 : state.dataThermalISO15099Calc->Ebf_NOSD,
980 0 : state.dataThermalISO15099Calc->Ebb_NOSD,
981 0 : state.dataThermalISO15099Calc->Rf_NOSD,
982 0 : state.dataThermalISO15099Calc->Rb_NOSD,
983 : ShadeEmisRatioOut_NOSD,
984 : ShadeEmisRatioIn_NOSD,
985 : ShadeHcModifiedOut_NOSD,
986 : ShadeHcModifiedIn_NOSD,
987 : ThermalMod,
988 : Debug_mode,
989 : AchievedErrorTolerance_NOSD,
990 : NumOfIter_NOSD,
991 : edgeGlCorrFac);
992 :
993 0 : NumOfIterations = NumOfIter_NOSD;
994 : // exit on error
995 0 : if (!(GoAhead(nperr))) {
996 0 : return;
997 : }
998 :
999 : // bi... Keep these values:
1000 0 : HcUnshadedOut = hcout_NOSD;
1001 0 : HcUnshadedIn = hcin_NOSD;
1002 :
1003 0 : ShadeHcRatioOut = ShadeHcModifiedOut / HcUnshadedOut;
1004 0 : ShadeHcRatioIn = ShadeHcModifiedIn / HcUnshadedIn;
1005 :
1006 : // bi...unshaded results:
1007 0 : if (files.WriteDebugOutput && (UnshadedDebug == 1)) {
1008 0 : WriteOutputArguments(files.DebugOutputFile,
1009 0 : files.DBGD,
1010 : nlayer_NOSD,
1011 : tamb,
1012 0 : state.dataThermalISO15099Calc->q_NOSD,
1013 0 : state.dataThermalISO15099Calc->qv_NOSD,
1014 0 : state.dataThermalISO15099Calc->qcgas_NOSD,
1015 0 : state.dataThermalISO15099Calc->qrgas_NOSD,
1016 0 : state.dataThermalISO15099Calc->theta_NOSD,
1017 0 : state.dataThermalISO15099Calc->vfreevent_NOSD,
1018 0 : state.dataThermalISO15099Calc->vvent_NOSD,
1019 0 : state.dataThermalISO15099Calc->Keff_NOSD,
1020 0 : state.dataThermalISO15099Calc->ShadeGapKeffConv_NOSD,
1021 : troom_NOSD,
1022 : ufactor_NOSD,
1023 : shgc_NOSD,
1024 : sc_NOSD,
1025 : hflux_NOSD,
1026 : shgct_NOSD,
1027 : hcin_NOSD,
1028 : hrin_NOSD,
1029 : hcout_NOSD,
1030 : hrout_NOSD,
1031 0 : state.dataThermalISO15099Calc->Ra_NOSD,
1032 0 : state.dataThermalISO15099Calc->Nu_NOSD,
1033 0 : state.dataThermalISO15099Calc->LayerType_NOSD,
1034 0 : state.dataThermalISO15099Calc->Ebf_NOSD,
1035 0 : state.dataThermalISO15099Calc->Ebb_NOSD,
1036 0 : state.dataThermalISO15099Calc->Rf_NOSD,
1037 0 : state.dataThermalISO15099Calc->Rb_NOSD,
1038 : ebsky,
1039 : Gout,
1040 : ebroom,
1041 : Gin,
1042 : ShadeEmisRatioIn_NOSD,
1043 : ShadeEmisRatioOut_NOSD,
1044 : ShadeHcRatioIn_NOSD,
1045 : ShadeHcRatioOut_NOSD,
1046 : hcin_NOSD,
1047 : hcout_NOSD,
1048 0 : state.dataThermalISO15099Calc->hcgas_NOSD,
1049 0 : state.dataThermalISO15099Calc->hrgas_NOSD,
1050 : AchievedErrorTolerance_NOSD,
1051 : NumOfIter_NOSD); // Autodesk:Uninit shgc_NOSD, sc_NOSD, hflux_NOSD,
1052 : // ShadeHcRatioIn_NOSD, ShadeHcRatioOut_NOSD were
1053 : // uninitialized
1054 : } // end if UnshadedDebug = 1
1055 : } // end if NeedUnshadedRun...
1056 :
1057 : // bi Set T6-related quantities keff, keffc: (using non-solar pass results)
1058 0 : if (nlayer > 1) {
1059 0 : for (i = 1; i <= nlayer - 1; ++i) {
1060 0 : Keff(i) = gap(i) * q(2 * i + 1) / (theta(2 * i + 1) - theta(2 * i));
1061 0 : if (IsShadingLayer(LayerType(i))) {
1062 0 : Keff(i) = gap(i) * q(2 * i + 1) / (theta(2 * i + 1) - theta(2 * i));
1063 : }
1064 0 : if (IsShadingLayer(LayerType(i + 1))) {
1065 0 : Keff(i) = gap(i) * q(2 * i + 1) / (theta(2 * i + 1) - theta(2 * i));
1066 : }
1067 0 : if (IsShadingLayer(LayerType(i))) {
1068 : // Keff(i) = gap(i) * qprim(2*i+1) / (theta(2*i+1) - theta(2*i))
1069 0 : if ((i > 1) && (i < nlayer)) {
1070 0 : tgg = gap(i - 1) + gap(i) + thick(i);
1071 0 : qc1 = state.dataThermalISO15099Calc->qcgas(i - 1);
1072 0 : qc2 = state.dataThermalISO15099Calc->qcgas(i);
1073 0 : qcgg = (qc1 + qc2) / 2.0;
1074 0 : ShadeGapKeffConv(i) = tgg * qcgg / (theta(2 * i + 1) - theta(2 * i - 2));
1075 : }
1076 : }
1077 : }
1078 : }
1079 :
1080 : } // if (UFactorCalc.ne.0) then
1081 :
1082 : // bi... For debugging purposes:
1083 0 : state.dataThermalISO15099Calc->qeff = ufactor * std::abs(tout - tind);
1084 0 : state.dataThermalISO15099Calc->flux_nonsolar = flux;
1085 :
1086 0 : if ((SHGCCalc > 0) && (dir > 0.0)) {
1087 0 : shgc = totsol - (state.dataThermalISO15099Calc->fluxs - flux) / dir;
1088 0 : sc = shgc / 0.87;
1089 0 : hcin = state.dataThermalISO15099Calc->hcins;
1090 0 : hrin = state.dataThermalISO15099Calc->hrins;
1091 0 : hin = state.dataThermalISO15099Calc->hins;
1092 0 : hcout = state.dataThermalISO15099Calc->hcouts;
1093 0 : hrout = state.dataThermalISO15099Calc->hrouts;
1094 0 : hout = state.dataThermalISO15099Calc->houts;
1095 0 : flux = state.dataThermalISO15099Calc->fluxs; // <--- ???
1096 0 : for (i = 1; i <= nlayer; ++i) {
1097 0 : theta(2 * i - 1) = state.dataThermalISO15099Calc->thetas(2 * i - 1);
1098 0 : theta(2 * i) = state.dataThermalISO15099Calc->thetas(2 * i);
1099 0 : state.dataThermalISO15099Calc->Ebb(i) = state.dataThermalISO15099Calc->Ebbs(i);
1100 0 : state.dataThermalISO15099Calc->Ebf(i) = state.dataThermalISO15099Calc->Ebfs(i);
1101 0 : state.dataThermalISO15099Calc->Rb(i) = state.dataThermalISO15099Calc->Rbs(i);
1102 0 : state.dataThermalISO15099Calc->Rf(i) = state.dataThermalISO15099Calc->Rfs(i);
1103 0 : q(2 * i - 1) = state.dataThermalISO15099Calc->qs(2 * i - 1);
1104 0 : q(2 * i) = state.dataThermalISO15099Calc->qs(2 * i);
1105 : // qprim(2*i - 1) = qprims(2*i - 1)
1106 : // qprim(2*i) = qprims(2*i)
1107 0 : qv(2 * i - 1) = state.dataThermalISO15099Calc->qvs(2 * i - 1);
1108 0 : qv(2 * i) = state.dataThermalISO15099Calc->qvs(2 * i);
1109 0 : hcgas(i) = state.dataThermalISO15099Calc->hcgass(i);
1110 0 : hrgas(i) = state.dataThermalISO15099Calc->hrgass(i);
1111 0 : state.dataThermalISO15099Calc->qcgas(i) = state.dataThermalISO15099Calc->qcgaps(i);
1112 0 : state.dataThermalISO15099Calc->qrgas(i) = state.dataThermalISO15099Calc->qrgaps(i);
1113 0 : AchievedErrorTolerance = AchievedErrorToleranceSolar;
1114 0 : NumOfIter = NumOfIterSolar;
1115 : }
1116 :
1117 : // bi CHECK THIS!
1118 0 : q(2 * nlayer + 1) = state.dataThermalISO15099Calc->qs(2 * nlayer + 1);
1119 : }
1120 :
1121 0 : hflux = flux; // save flux value for output table
1122 :
1123 : // bi... Write results to debug output file:
1124 0 : if (files.WriteDebugOutput) {
1125 0 : WriteOutputArguments(files.DebugOutputFile,
1126 0 : files.DBGD,
1127 : nlayer,
1128 : tamb,
1129 : q,
1130 : qv,
1131 0 : state.dataThermalISO15099Calc->qcgas,
1132 0 : state.dataThermalISO15099Calc->qrgas,
1133 : theta,
1134 0 : state.dataThermalISO15099Calc->vfreevent,
1135 : vvent,
1136 : Keff,
1137 : ShadeGapKeffConv,
1138 : troom,
1139 : ufactor,
1140 : shgc,
1141 : sc,
1142 : hflux,
1143 : shgct,
1144 : hcin,
1145 : hrin,
1146 : hcout,
1147 : hrout,
1148 : Ra,
1149 : Nu,
1150 : LayerType,
1151 0 : state.dataThermalISO15099Calc->Ebf,
1152 0 : state.dataThermalISO15099Calc->Ebb,
1153 0 : state.dataThermalISO15099Calc->Rf,
1154 0 : state.dataThermalISO15099Calc->Rb,
1155 : ebsky,
1156 : Gout,
1157 : ebroom,
1158 : Gin,
1159 : ShadeEmisRatioIn,
1160 : ShadeEmisRatioOut,
1161 : ShadeHcRatioIn,
1162 : ShadeHcRatioOut,
1163 : HcUnshadedIn,
1164 : HcUnshadedOut,
1165 : hcgas,
1166 : hrgas,
1167 : AchievedErrorTolerance,
1168 : NumOfIter);
1169 : } // if WriteDebugOutput.eq.true - writing output file
1170 : }
1171 :
1172 0 : void therm1d(EnergyPlusData &state,
1173 : Files &files,
1174 : int const nlayer,
1175 : int const iwd,
1176 : Real64 &tout,
1177 : Real64 &tind,
1178 : Real64 const wso,
1179 : Real64 const wsi,
1180 : Real64 const VacuumPressure,
1181 : Real64 const VacuumMaxGapThickness,
1182 : Real64 const dir,
1183 : Real64 &ebsky,
1184 : Real64 const Gout,
1185 : Real64 const trmout,
1186 : Real64 const trmin,
1187 : Real64 &ebroom,
1188 : Real64 const Gin,
1189 : const Array1D<Real64> &tir,
1190 : const Array1D<Real64> &rir,
1191 : const Array1D<Real64> &emis,
1192 : const Array1D<Real64> &gap,
1193 : const Array1D<Real64> &thick,
1194 : const Array1D<Real64> &scon,
1195 : Real64 const tilt,
1196 : const Array1D<Real64> &asol,
1197 : Real64 const height,
1198 : Real64 const heightt,
1199 : Real64 const width,
1200 : Array2_int const &iprop,
1201 : Array2<Real64> const &frct,
1202 : const Array1D<Real64> &presure,
1203 : const Array1D_int &nmix,
1204 : const Array1D<Real64> &wght,
1205 : Array2<Real64> const &gcon,
1206 : Array2<Real64> const &gvis,
1207 : Array2<Real64> const &gcp,
1208 : const Array1D<Real64> &gama,
1209 : const Array1D_int &SupportPillar,
1210 : const Array1D<Real64> &PillarSpacing,
1211 : const Array1D<Real64> &PillarRadius,
1212 : Array1D<Real64> &theta,
1213 : Array1D<Real64> &q,
1214 : Array1D<Real64> &qv,
1215 : Real64 &flux,
1216 : Real64 &hcin,
1217 : Real64 &hrin,
1218 : Real64 &hcout,
1219 : Real64 &hrout,
1220 : Real64 const hin,
1221 : Real64 const hout,
1222 : Array1D<Real64> &hcgas,
1223 : Array1D<Real64> &hrgas,
1224 : Real64 &ufactor,
1225 : int &nperr,
1226 : std::string &ErrorMessage,
1227 : Real64 &tamb,
1228 : Real64 &troom,
1229 : const Array1D_int &ibc,
1230 : const Array1D<Real64> &Atop,
1231 : const Array1D<Real64> &Abot,
1232 : const Array1D<Real64> &Al,
1233 : const Array1D<Real64> &Ar,
1234 : const Array1D<Real64> &Ah,
1235 : const Array1D<Real64> &EffectiveOpenness,
1236 : const Array1D<Real64> &vvent,
1237 : const Array1D<Real64> &tvent,
1238 : const Array1D<TARCOGLayerType> &LayerType,
1239 : Array1D<Real64> &Ra,
1240 : Array1D<Real64> &Nu,
1241 : Array1D<Real64> &vfreevent,
1242 : Array1D<Real64> &qcgas,
1243 : Array1D<Real64> &qrgas,
1244 : Array1D<Real64> &Ebf,
1245 : Array1D<Real64> &Ebb,
1246 : Array1D<Real64> &Rf,
1247 : Array1D<Real64> &Rb,
1248 : Real64 &ShadeEmisRatioOut,
1249 : Real64 &ShadeEmisRatioIn,
1250 : Real64 &ShadeHcModifiedOut,
1251 : Real64 &ShadeHcModifiedIn,
1252 : TARCOGThermalModel ThermalMod,
1253 : [[maybe_unused]] int const Debug_mode,
1254 : Real64 &AchievedErrorTolerance,
1255 : int &TotalIndex,
1256 : Real64 const edgeGlCorrFac)
1257 : {
1258 : //********************************************************************************
1259 : // Main subroutine for calculation of 1-D heat transfer in the center of glazing.
1260 : //********************************************************************************
1261 : // Inputs
1262 : // nlayer number of solid layers
1263 : // iwd wind direction
1264 : // tout outside temp in k
1265 : // tind inside temp in k
1266 : // wso wind speed in m/s
1267 : // wsi inside forced air speed m/s
1268 : // Ebsky ir flux from outside
1269 : // Gout back facing radiosity from outside
1270 : // Trmout Mean outdoor radiant temperature
1271 : // Trmin Mean indoor radiant temperature
1272 : // Ebroom ir flux from room
1273 : // Gin front facing radiosity from room
1274 : // tir ir transmittance of each layer
1275 : // rir ir reflectance of each surface
1276 : // emis ir emittances of each surface
1277 : // gap array of gap widths in meters
1278 : // thick thickness of glazing layers (m)
1279 : // scon Vector of conductivities of 'glazing' layers
1280 : // tilt Window tilt (deg). vert: tilt=90, hor out up: tilt=0, hor out down: tilt=180
1281 : // sol absorbed solar energy for each layer in w/m2
1282 : // height glazing cavity height
1283 : // heightt
1284 : // iprop
1285 : // frct
1286 : // presure
1287 : // nmix vector of number of gasses in a mixture for each gap
1288 : // hin convective indoor film coefficient (if non-zero hin input)
1289 : // hout convective outdoor film coeff. (if non-zero hout input)
1290 : // outputs
1291 : // theta temp distribution in k
1292 : // flux net heat flux between room and window
1293 : // rtot overall thermal resistance
1294 : // rs ?
1295 : // hcin convective indoor film coeff.
1296 : // hrin radiative part of indoor film coeff.
1297 : // hcout convective outdoor film coeff.
1298 : // hrout radiative part of outdoor film coeff.
1299 : // hin convective indoor film coefficient
1300 : // hout convective outdoor film coeff.
1301 : // ufactor overall u-factor
1302 : // qcgap vector of convective/conductive parts of flux in gaps
1303 : // qrgap vector of radiative parts of flux in gaps
1304 : // nperr
1305 : // *Inactives**
1306 : // wa - window azimuth (degrees, clockwise from south)
1307 : // hgas matrix of gap film coefficients
1308 : // Locals
1309 : // Ebb Vector
1310 : // Ebf Vector
1311 : // Rb Vector
1312 : // Rf Vector
1313 : // a Array
1314 : // b Array
1315 : // hhat Vector
1316 : // err iteration tolerance
1317 : // dtmax max temp difference after iteration
1318 : // index iteration step
1319 :
1320 : // Using
1321 : // Locals
1322 : // 0 - don't create debug output files
1323 : // 1 - append results to existing debug output file
1324 : // 2 - store results in new debug output file
1325 : // 3 - save in-between results (in all iterations) to existing debug file
1326 :
1327 0 : Array2D<Real64> a(4 * nlayer, 4 * nlayer);
1328 0 : Array1D<Real64> b(4 * nlayer);
1329 : // REAL(r64) :: hhatv(maxlay3),hcv(maxlay3), Ebgap(maxlay3), Tgap(maxlay1)
1330 :
1331 : // REAL(r64) :: alpha
1332 : int maxiter;
1333 :
1334 : Real64 qr_gap_out;
1335 : Real64 qr_gap_in;
1336 :
1337 0 : Array1D<Real64> told(2 * nlayer);
1338 :
1339 : // Simon: parameters used in case of JCFN iteration method
1340 0 : Array1D<Real64> FRes({1, 4 * nlayer}); // store function results from current iteration
1341 0 : Array1D<Real64> FResOld({1, 4 * nlayer}); // store function results from previous iteration
1342 0 : Array1D<Real64> FResDiff({1, 4 * nlayer}); // save difference in results between iterations
1343 0 : Array1D<Real64> Radiation({1, 2 * nlayer}); // radiation on layer surfaces. used as temporary storage during iterations
1344 :
1345 0 : Array1D<Real64> x({1, 4 * nlayer}); // temporary vector for storing results (theta and Radiation). used for easier handling
1346 0 : Array1D<Real64> dX({1, 4 * nlayer}, 0.0); // difference in results
1347 0 : Array2D<Real64> Jacobian({1, 4 * nlayer}, {1, 4 * nlayer}); // diagonal vector for jacobian computation-free newton method
1348 0 : Array1D<Real64> DRes({1, 4 * nlayer}); // used in jacobian forward-difference approximation
1349 :
1350 : // This is used to store matrix before equation solver. It is important because solver destroys
1351 : // content of matrices
1352 0 : Array2D<Real64> LeftHandSide({1, 4 * nlayer}, {1, 4 * nlayer});
1353 0 : Array1D<Real64> RightHandSide({1, 4 * nlayer});
1354 :
1355 : // Simon: Keep best achieved convergence
1356 : Real64 prevDifference;
1357 : Real64 Relaxation;
1358 0 : Array1D<Real64> RadiationSave({1, 2 * nlayer});
1359 0 : Array1D<Real64> thetaSave({1, 2 * nlayer});
1360 : int currentTry;
1361 :
1362 : int CSMFlag;
1363 : int i;
1364 : int j;
1365 : int k;
1366 : Real64 curDifference;
1367 : int index;
1368 : int curTempCorrection;
1369 :
1370 : Real64 qc_gap_in;
1371 : Real64 hc_modified_in;
1372 :
1373 : CalculationOutcome CalcOutcome;
1374 :
1375 : bool iterationsFinished; // To mark whether or not iterations are finished
1376 : bool saveIterationResults;
1377 : bool updateGapTemperature;
1378 : // logical :: TurnOnNewton
1379 :
1380 0 : int SDLayerIndex = -1;
1381 :
1382 0 : Array1D<Real64> sconScaled(maxlay);
1383 :
1384 : // Simon: This is set to zero until it is resolved what to do with modifier
1385 0 : ShadeHcModifiedOut = 0.0;
1386 0 : CSMFlag = 0;
1387 0 : CalcOutcome = CalculationOutcome::Invalid;
1388 0 : curTempCorrection = 0;
1389 0 : AchievedErrorTolerance = 0.0;
1390 0 : curDifference = 0.0;
1391 0 : currentTry = 0;
1392 0 : index = 0;
1393 0 : TotalIndex = 0;
1394 0 : iterationsFinished = false;
1395 0 : qv = 0.0;
1396 0 : Ebb = 0.0;
1397 0 : Ebf = 0.0;
1398 0 : Rb = 0.0;
1399 0 : Rf = 0.0;
1400 0 : a = 0.0;
1401 0 : b = 0.0;
1402 :
1403 0 : FRes = 0.0;
1404 0 : FResOld = 0.0;
1405 0 : FResDiff = 0.0;
1406 0 : Radiation = 0.0;
1407 0 : Relaxation = RelaxationStart;
1408 :
1409 0 : maxiter = NumOfIterations;
1410 :
1411 0 : saveIterationResults = false;
1412 :
1413 0 : for (i = 1; i <= nlayer; ++i) {
1414 0 : k = 2 * i;
1415 0 : Radiation(k) = Ebb(i);
1416 0 : Radiation(k - 1) = Ebf(i);
1417 0 : told(k - 1) = 0.0;
1418 0 : told(k) = 0.0;
1419 : }
1420 :
1421 : // bi...Set LayerTypeSpec array - need to treat venetians AND woven shades as glass:
1422 0 : if (ThermalMod == TARCOGThermalModel::CSM) {
1423 0 : for (i = 1; i <= nlayer; ++i) {
1424 0 : if (IsShadingLayer(LayerType(i))) {
1425 : // LayerTypeSpec( i ) = 0; //Unused
1426 0 : SDLayerIndex = i;
1427 : } else {
1428 : // LayerTypeSpec( i ) = LayerType( i ); //Unused
1429 : }
1430 : }
1431 : }
1432 :
1433 : // first store results before iterations begin
1434 0 : if (saveIterationResults) {
1435 0 : storeIterationResults(state,
1436 : files,
1437 : nlayer,
1438 : index,
1439 : theta,
1440 : trmout,
1441 : tamb,
1442 : trmin,
1443 : troom,
1444 : ebsky,
1445 : ebroom,
1446 : hcin,
1447 : hcout,
1448 : hrin,
1449 : hrout,
1450 : hin,
1451 : hout,
1452 : Ebb,
1453 : Ebf,
1454 : Rb,
1455 : Rf,
1456 : nperr);
1457 : }
1458 :
1459 0 : state.dataThermalISO15099Calc->Tgap(1) = tout;
1460 0 : state.dataThermalISO15099Calc->Tgap(nlayer + 1) = tind;
1461 0 : for (i = 2; i <= nlayer; ++i) {
1462 0 : state.dataThermalISO15099Calc->Tgap(i) = (theta(2 * i - 1) + theta(2 * i - 2)) / 2;
1463 : }
1464 : //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1465 : //!!! MAIN ITERATION LOOP
1466 : //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1467 0 : while (!(iterationsFinished)) {
1468 :
1469 0 : for (i = 1; i <= 2 * nlayer; ++i) {
1470 0 : if (theta(i) < 0) {
1471 0 : theta(i) = 1.0 * i;
1472 : }
1473 : }
1474 :
1475 : // do i=1,nlayer+1
1476 : // if (i == 1) then
1477 : // Tgap(i) = tout
1478 : // else if (i == nlayer+1) then
1479 : // Tgap(i) = tind
1480 : // else
1481 : // Tgap(i) = (theta(2*i-2) + theta(2*i-1)) / 2.0d0
1482 : // end if
1483 : // end do
1484 :
1485 : // skip updating gap temperatures for shading devices. Gap temperature in that case is not simply average
1486 : // between two layer temperatures
1487 0 : for (i = 2; i <= nlayer; ++i) {
1488 0 : updateGapTemperature = false;
1489 0 : if ((!(IsShadingLayer(LayerType(i - 1)))) && (!(IsShadingLayer(LayerType(i))))) {
1490 0 : updateGapTemperature = true;
1491 : }
1492 0 : if (updateGapTemperature) {
1493 0 : state.dataThermalISO15099Calc->Tgap(i) = (theta(2 * i - 1) + theta(2 * i - 2)) / 2;
1494 : }
1495 : }
1496 :
1497 : // evaluate convective/conductive components of gap
1498 0 : hatter(state,
1499 : nlayer,
1500 : iwd,
1501 : tout,
1502 : tind,
1503 : wso,
1504 : wsi,
1505 : VacuumPressure,
1506 : VacuumMaxGapThickness,
1507 : ebsky,
1508 : tamb,
1509 : ebroom,
1510 : troom,
1511 : gap,
1512 : height,
1513 : heightt,
1514 : scon,
1515 : tilt,
1516 : theta,
1517 0 : state.dataThermalISO15099Calc->Tgap,
1518 : Radiation,
1519 : trmout,
1520 : trmin,
1521 : iprop,
1522 : frct,
1523 : presure,
1524 : nmix,
1525 : wght,
1526 : gcon,
1527 : gvis,
1528 : gcp,
1529 : gama,
1530 : SupportPillar,
1531 : PillarSpacing,
1532 : PillarRadius,
1533 0 : state.dataThermalISO15099Calc->hgas,
1534 : hcgas,
1535 : hrgas,
1536 : hcin,
1537 : hcout,
1538 : hin,
1539 : hout,
1540 : index,
1541 : ibc,
1542 : nperr,
1543 : ErrorMessage,
1544 : hrin,
1545 : hrout,
1546 : Ra,
1547 : Nu);
1548 :
1549 0 : effectiveLayerCond(state,
1550 : nlayer,
1551 : LayerType,
1552 : scon,
1553 : thick,
1554 : iprop,
1555 : frct,
1556 : nmix,
1557 : presure,
1558 : wght,
1559 : gcon,
1560 : gvis,
1561 : gcp,
1562 : EffectiveOpenness,
1563 : theta,
1564 : sconScaled,
1565 : nperr,
1566 : ErrorMessage);
1567 :
1568 : // exit on error
1569 0 : if (!(GoAhead(nperr))) {
1570 0 : return;
1571 : }
1572 :
1573 : // bi...Override hhat values near SHADING DEVICE layer(s), but only for CSM thermal model:
1574 0 : if ((ThermalMod == TARCOGThermalModel::CSM) && (SDLayerIndex > 0)) {
1575 : // adjust hhat values
1576 : // call adjusthhat(SDLayerIndex, ibc, tout, tind, nlayer, theta, wso, wsi, iwd, height, heightt, tilt, &
1577 : // & thick, gap, hout, hrout, hin, hrin, iprop, frct, presure, nmix, wght, gcon, gvis, gcp, &
1578 : // index, SDScalar, Ebf, Ebb, hgas, hhat, nperr, ErrorMessage)
1579 : // do i = 1, maxlay3
1580 : // hhatv(i) = 0.0d0
1581 : // Ebgap(i) = 0.0d0
1582 : // qv(i) = 0.0d0
1583 : // hcv(i) = 0.0d0
1584 : // end do
1585 0 : matrixQBalance(nlayer,
1586 : a,
1587 : b,
1588 : sconScaled,
1589 : hcgas,
1590 0 : state.dataThermalISO15099Calc->hcgapMod,
1591 : asol,
1592 : qv,
1593 0 : state.dataThermalISO15099Calc->hcv,
1594 : tind,
1595 : tout,
1596 : Gin,
1597 : Gout,
1598 : theta,
1599 : tir,
1600 : rir,
1601 : emis,
1602 : edgeGlCorrFac);
1603 : } else {
1604 : // bi...There are no Venetian layers, or ThermalMod is not CSM, so carry on as usual:
1605 0 : shading(state,
1606 : theta,
1607 : gap,
1608 0 : state.dataThermalISO15099Calc->hgas,
1609 : hcgas,
1610 : hrgas,
1611 : frct,
1612 : iprop,
1613 : presure,
1614 : nmix,
1615 : wght,
1616 : gcon,
1617 : gvis,
1618 : gcp,
1619 : nlayer,
1620 : width,
1621 : height,
1622 : tilt,
1623 : tout,
1624 : tind,
1625 : Atop,
1626 : Abot,
1627 : Al,
1628 : Ar,
1629 : Ah,
1630 : vvent,
1631 : tvent,
1632 : LayerType,
1633 0 : state.dataThermalISO15099Calc->Tgap,
1634 : qv,
1635 0 : state.dataThermalISO15099Calc->hcv,
1636 : nperr,
1637 : ErrorMessage,
1638 : vfreevent);
1639 :
1640 : // exit on error
1641 0 : if (!(GoAhead(nperr))) {
1642 0 : return;
1643 : }
1644 :
1645 0 : matrixQBalance(nlayer,
1646 : a,
1647 : b,
1648 : sconScaled,
1649 : hcgas,
1650 0 : state.dataThermalISO15099Calc->hcgapMod,
1651 : asol,
1652 : qv,
1653 0 : state.dataThermalISO15099Calc->hcv,
1654 : tind,
1655 : tout,
1656 : Gin,
1657 : Gout,
1658 : theta,
1659 : tir,
1660 : rir,
1661 : emis,
1662 : edgeGlCorrFac);
1663 :
1664 : } // end if
1665 :
1666 0 : FResOld = FRes;
1667 :
1668 : // Pack results in one array
1669 0 : for (i = 1; i <= nlayer; ++i) {
1670 0 : k = 4 * i - 3;
1671 0 : j = 2 * i - 1;
1672 :
1673 0 : x(k) = theta(j);
1674 0 : x(k + 1) = Radiation(j);
1675 0 : x(k + 2) = Radiation(j + 1);
1676 0 : x(k + 3) = theta(j + 1);
1677 : }
1678 :
1679 0 : CalculateFuncResults(nlayer, a, b, x, FRes);
1680 :
1681 0 : FResDiff = FRes - FResOld;
1682 :
1683 0 : LeftHandSide = a;
1684 0 : RightHandSide = b;
1685 0 : EquationsSolver(state, LeftHandSide, RightHandSide, 4 * nlayer, nperr, ErrorMessage);
1686 :
1687 : // Simon: This is much better, but also much slower convergence criteria. Think of how to make this flexible and allow
1688 : // user to change this from outside (through argument passing)
1689 : // curDifference = ABS(FRes(1))
1690 : // do i = 2, 4*nlayer
1691 : // curDifference = MAX(curDifference, ABS(FRes(i)))
1692 : // curDifference = curDifference + ABS(FRes(i))
1693 : // end do
1694 :
1695 0 : curDifference = std::abs(theta(1) - told(1));
1696 : // curDifference = ABS(FRes(1))
1697 0 : for (i = 2; i <= 2 * nlayer; ++i) {
1698 : // do i = 2, 4*nlayer
1699 0 : curDifference = max(curDifference, std::abs(theta(i) - told(i)));
1700 : // curDifference = MAX(ABS(FRes(i)), curDifference)
1701 : }
1702 :
1703 0 : for (i = 1; i <= nlayer; ++i) {
1704 0 : k = 4 * i - 3;
1705 0 : j = 2 * i - 1;
1706 : // if (TurnOnNewton) then
1707 : // theta(j) = theta(j) + Relaxation*dx(k)
1708 : // theta(j+1) = theta(j+1) + Relaxation*dx(k+1)
1709 : // Radiation(j) = Radiation(j) + Relaxation*dx(k+2)
1710 : // Radiation(j+1) = Radiation(j+1) + Relaxation*dx(k+3)
1711 : // else
1712 : // dX(k) = RightHandSide(k) - theta(j)
1713 : // dX(k+1) = RightHandSide(k + 1) - theta(j+1)
1714 : // dX(k+2) = RightHandSide(k + 2) - Radiation(j)
1715 : // dX(k+3) = RightHandSide(k + 3) - Radiation(j+1)
1716 0 : told(j) = theta(j);
1717 0 : told(j + 1) = theta(j + 1);
1718 0 : theta(j) = (1 - Relaxation) * theta(j) + Relaxation * RightHandSide(k);
1719 0 : Radiation(j) = (1 - Relaxation) * Radiation(j) + Relaxation * RightHandSide(k + 1);
1720 0 : Radiation(j + 1) = (1 - Relaxation) * Radiation(j + 1) + Relaxation * RightHandSide(k + 2);
1721 0 : theta(j + 1) = (1 - Relaxation) * theta(j + 1) + Relaxation * RightHandSide(k + 3);
1722 : // end if
1723 : }
1724 :
1725 : // it is important not to update gaps around shading layers since that is already calculated by
1726 : // shading routines
1727 0 : for (i = 1; i <= nlayer + 1; ++i) {
1728 0 : updateGapTemperature = true;
1729 0 : if ((i == 1) || (i == nlayer + 1)) {
1730 : // update gap array with interior and exterior temperature
1731 0 : updateGapTemperature = true;
1732 : } else {
1733 : // update gap temperature only if gap on both sides
1734 0 : updateGapTemperature = false;
1735 0 : if ((!(IsShadingLayer(LayerType(i - 1)))) && (!(IsShadingLayer(LayerType(i))))) {
1736 0 : updateGapTemperature = true;
1737 : }
1738 : }
1739 0 : j = 2 * (i - 1);
1740 0 : if (updateGapTemperature) {
1741 0 : if (i == 1) {
1742 0 : state.dataThermalISO15099Calc->Tgap(1) = tout;
1743 0 : } else if (i == (nlayer + 1)) {
1744 0 : state.dataThermalISO15099Calc->Tgap(i) = tind;
1745 : } else {
1746 0 : state.dataThermalISO15099Calc->Tgap(i) = (theta(j) + theta(j + 1)) / 2;
1747 : }
1748 : }
1749 : }
1750 :
1751 : // and store results during iterations
1752 0 : if (saveIterationResults) {
1753 0 : storeIterationResults(state,
1754 : files,
1755 : nlayer,
1756 : index + 1,
1757 : theta,
1758 : trmout,
1759 : tamb,
1760 : trmin,
1761 : troom,
1762 : ebsky,
1763 : ebroom,
1764 : hcin,
1765 : hcout,
1766 : hrin,
1767 : hrout,
1768 : hin,
1769 : hout,
1770 : Ebb,
1771 : Ebf,
1772 : Rb,
1773 : Rf,
1774 : nperr);
1775 : }
1776 :
1777 0 : if (!(GoAhead(nperr))) {
1778 0 : return;
1779 : }
1780 :
1781 0 : prevDifference = curDifference;
1782 :
1783 0 : if ((index == 0) || (curDifference < AchievedErrorTolerance)) {
1784 0 : AchievedErrorTolerance = curDifference;
1785 0 : currentTry = 0;
1786 0 : for (i = 1; i <= 2 * nlayer; ++i) {
1787 0 : RadiationSave(i) = Radiation(i);
1788 0 : thetaSave(i) = theta(i);
1789 : }
1790 : } else {
1791 : // This is case when program solution diverged
1792 0 : ++currentTry;
1793 0 : if (currentTry >= NumOfTries) {
1794 0 : currentTry = 0;
1795 0 : for (i = 1; i <= 2 * nlayer; ++i) {
1796 0 : Radiation(i) = RadiationSave(i);
1797 0 : theta(i) = thetaSave(i);
1798 : }
1799 : // if (.not.TurnOnNewton) then
1800 : // TurnOnNewton = .TRUE.
1801 : // else
1802 0 : Relaxation -= RelaxationDecrease;
1803 0 : TotalIndex += index;
1804 0 : index = 0;
1805 : // Start from best achieved convergence
1806 0 : if (Relaxation <= 0.0) { // cannot continue with relaxation equal to zero
1807 0 : iterationsFinished = true;
1808 : }
1809 : // TurnOnNewton = .TRUE.
1810 : // end if ! if (.not.TurnOnNewton) then
1811 : } // f (currentTry == NumOfTries) then
1812 : }
1813 :
1814 : // Check if results were found:
1815 0 : if (curDifference < ConvergenceTolerance) {
1816 0 : CalcOutcome = CalculationOutcome::OK;
1817 0 : TotalIndex += index;
1818 0 : iterationsFinished = true;
1819 : }
1820 :
1821 0 : if (index >= maxiter) {
1822 0 : Relaxation -= RelaxationDecrease;
1823 0 : TotalIndex += index;
1824 0 : index = 0;
1825 : // TurnOnNewton = .TRUE.
1826 :
1827 : // Start from best achieved convergence
1828 0 : for (i = 1; i <= 2 * nlayer; ++i) {
1829 0 : Radiation(i) = RadiationSave(i);
1830 0 : theta(i) = thetaSave(i);
1831 : }
1832 0 : if (Relaxation <= 0.0) { // cannot continue with relaxation equal to zero
1833 0 : iterationsFinished = true;
1834 : }
1835 : }
1836 :
1837 0 : ++index;
1838 : }
1839 :
1840 : // Get results from closest iteration and store it
1841 0 : if (CalcOutcome == CalculationOutcome::OK) {
1842 0 : for (i = 1; i <= 2 * nlayer; ++i) {
1843 0 : Radiation(i) = RadiationSave(i);
1844 0 : theta(i) = thetaSave(i);
1845 : }
1846 :
1847 0 : for (i = 2; i <= nlayer; ++i) {
1848 0 : updateGapTemperature = false;
1849 0 : if ((!(IsShadingLayer(LayerType(i - 1)))) && (!(IsShadingLayer(LayerType(i))))) {
1850 0 : updateGapTemperature = true;
1851 : }
1852 :
1853 0 : if (updateGapTemperature) {
1854 0 : state.dataThermalISO15099Calc->Tgap(i) = (theta(2 * i - 1) + theta(2 * i - 2)) / 2;
1855 : }
1856 : }
1857 :
1858 : // Simon: It is important to recalculate coefficients from most accurate run
1859 0 : hatter(state,
1860 : nlayer,
1861 : iwd,
1862 : tout,
1863 : tind,
1864 : wso,
1865 : wsi,
1866 : VacuumPressure,
1867 : VacuumMaxGapThickness,
1868 : ebsky,
1869 : tamb,
1870 : ebroom,
1871 : troom,
1872 : gap,
1873 : height,
1874 : heightt,
1875 : scon,
1876 : tilt,
1877 : theta,
1878 0 : state.dataThermalISO15099Calc->Tgap,
1879 : Radiation,
1880 : trmout,
1881 : trmin,
1882 : iprop,
1883 : frct,
1884 : presure,
1885 : nmix,
1886 : wght,
1887 : gcon,
1888 : gvis,
1889 : gcp,
1890 : gama,
1891 : SupportPillar,
1892 : PillarSpacing,
1893 : PillarRadius,
1894 0 : state.dataThermalISO15099Calc->hgas,
1895 : hcgas,
1896 : hrgas,
1897 : hcin,
1898 : hcout,
1899 : hin,
1900 : hout,
1901 : index,
1902 : ibc,
1903 : nperr,
1904 : ErrorMessage,
1905 : hrin,
1906 : hrout,
1907 : Ra,
1908 : Nu);
1909 :
1910 0 : shading(state,
1911 : theta,
1912 : gap,
1913 0 : state.dataThermalISO15099Calc->hgas,
1914 : hcgas,
1915 : hrgas,
1916 : frct,
1917 : iprop,
1918 : presure,
1919 : nmix,
1920 : wght,
1921 : gcon,
1922 : gvis,
1923 : gcp,
1924 : nlayer,
1925 : width,
1926 : height,
1927 : tilt,
1928 : tout,
1929 : tind,
1930 : Atop,
1931 : Abot,
1932 : Al,
1933 : Ar,
1934 : Ah,
1935 : vvent,
1936 : tvent,
1937 : LayerType,
1938 0 : state.dataThermalISO15099Calc->Tgap,
1939 : qv,
1940 0 : state.dataThermalISO15099Calc->hcv,
1941 : nperr,
1942 : ErrorMessage,
1943 : vfreevent);
1944 : }
1945 :
1946 0 : if (CalcOutcome == CalculationOutcome::Invalid) {
1947 0 : ErrorMessage = "Tarcog failed to converge";
1948 0 : nperr = 2; // error 2: failed to converge...
1949 : }
1950 :
1951 : // Get radiation results first
1952 : // if (curEquationsApproach.eq.eaQBalance) then
1953 0 : for (i = 1; i <= nlayer; ++i) {
1954 0 : k = 2 * i - 1;
1955 0 : Rf(i) = Radiation(k);
1956 0 : Rb(i) = Radiation(k + 1);
1957 0 : Ebf(i) = Constant::StefanBoltzmann * pow_4(theta(k));
1958 0 : Ebb(i) = Constant::StefanBoltzmann * pow_4(theta(k + 1));
1959 : }
1960 : // end if
1961 :
1962 : // Finishing calcs:
1963 0 : resist(nlayer, trmout, tout, trmin, tind, hcgas, hrgas, theta, q, qv, LayerType, thick, scon, ufactor, flux, qcgas, qrgas);
1964 :
1965 : // bi... Set T6-related quantities - ratios for modified epsilon, hc for modelling external SDs:
1966 : // (using non-solar pass results)
1967 0 : if ((dir == 0.0) && (nlayer > 1)) {
1968 :
1969 0 : qr_gap_out = Rf(2) - Rb(1);
1970 0 : qr_gap_in = Rf(nlayer) - Rb(nlayer - 1);
1971 :
1972 0 : if (IsShadingLayer(LayerType(1))) {
1973 0 : ShadeEmisRatioOut = qr_gap_out / (emis(3) * Constant::StefanBoltzmann * (pow_4(theta(3)) - pow_4(trmout)));
1974 : // qc_gap_out = qprim(3) - qr_gap_out
1975 : // qcgapout2 = qcgas(1)
1976 : // Hc_modified_out = (qc_gap_out / (theta(3) - tout))
1977 : // ShadeHcModifiedOut = Hc_modified_out
1978 : }
1979 :
1980 0 : if (IsShadingLayer(LayerType(nlayer))) {
1981 0 : ShadeEmisRatioIn = qr_gap_in / (emis(2 * nlayer - 2) * Constant::StefanBoltzmann * (pow_4(trmin) - pow_4(theta(2 * nlayer - 2))));
1982 0 : qc_gap_in = q(2 * nlayer - 1) - qr_gap_in;
1983 0 : hc_modified_in = (qc_gap_in / (tind - theta(2 * nlayer - 2)));
1984 0 : ShadeHcModifiedIn = hc_modified_in;
1985 : }
1986 : } // IF dir = 0
1987 0 : }
1988 :
1989 0 : void guess(Real64 const tout,
1990 : Real64 const tind,
1991 : int const nlayer,
1992 : const Array1D<Real64> &gap,
1993 : const Array1D<Real64> &thick,
1994 : Real64 &width,
1995 : Array1D<Real64> &theta,
1996 : Array1D<Real64> &Ebb,
1997 : Array1D<Real64> &Ebf,
1998 : Array1D<Real64> &Tgap)
1999 : {
2000 : //***********************************************************************
2001 : // purpose - initializes temperature distribution assuming
2002 : // a constant temperature gradient across the window
2003 : //***********************************************************************
2004 : // Input
2005 : // tout outdoor air temperature (k)
2006 : // tind indoor air temperature (k)
2007 : // nlayer number of solid layers in window output
2008 : // gap thickness of gas gaps (m)
2009 : // thick thickness of glazing layers (m)
2010 : // Output
2011 : // width total width of the glazing system
2012 : // theta array of surface temps starting from outdoor layer (k)
2013 : // Ebb vector of emissive power (?) of the back surface (# of layers)
2014 : // Ebf vector of emissive power (?) of the front surface (# of layers)
2015 : // Locals
2016 : // x Vector of running width
2017 : // delta delta T per unit length
2018 :
2019 : // Using
2020 : // Argument array dimensioning
2021 0 : EP_SIZE_CHECK(gap, MaxGap);
2022 0 : EP_SIZE_CHECK(thick, maxlay);
2023 0 : EP_SIZE_CHECK(theta, maxlay2);
2024 0 : EP_SIZE_CHECK(Ebb, maxlay);
2025 0 : EP_SIZE_CHECK(Ebf, maxlay);
2026 0 : EP_SIZE_CHECK(Tgap, maxlay1);
2027 :
2028 : // Locals
2029 0 : Array1D<Real64> x(maxlay2);
2030 : Real64 delta;
2031 : int i;
2032 : int j;
2033 : int k;
2034 :
2035 0 : x(1) = 0.001;
2036 0 : x(2) = x(1) + thick(1);
2037 :
2038 0 : for (i = 2; i <= nlayer; ++i) {
2039 0 : j = 2 * i - 1;
2040 0 : k = 2 * i;
2041 0 : x(j) = x(j - 1) + gap(i - 1);
2042 0 : x(k) = x(k - 1) + thick(i);
2043 : }
2044 :
2045 0 : width = x(nlayer * 2) + 0.01;
2046 0 : delta = (tind - tout) / width;
2047 :
2048 0 : if (delta == 0.0) {
2049 0 : delta = TemperatureQuessDiff / width;
2050 : }
2051 :
2052 0 : for (i = 1; i <= nlayer; ++i) {
2053 0 : j = 2 * i;
2054 0 : theta(j - 1) = tout + x(j - 1) * delta;
2055 0 : theta(j) = tout + x(j) * delta;
2056 0 : Ebf(i) = Constant::StefanBoltzmann * pow_4(theta(j - 1));
2057 0 : Ebb(i) = Constant::StefanBoltzmann * pow_4(theta(j));
2058 : }
2059 :
2060 0 : for (i = 1; i <= nlayer + 1; ++i) {
2061 0 : if (i == 1) {
2062 0 : Tgap(1) = tout;
2063 0 : } else if (i == (nlayer + 1)) {
2064 0 : Tgap(nlayer + 1) = tind;
2065 : } else {
2066 0 : Tgap(i) = (theta(2 * i - 1) + theta(2 * i - 2)) / 2;
2067 : }
2068 : }
2069 0 : }
2070 :
2071 0 : void solarISO15099(Real64 const totsol, Real64 const rtot, const Array1D<Real64> &rs, int const nlayer, const Array1D<Real64> &absol, Real64 &sf)
2072 : {
2073 : //***********************************************************************
2074 : // This subroutine calculates the shading coefficient for a window.
2075 : //***********************************************************************
2076 : // Inputs:
2077 : // absol array of absorbed fraction of solar radiation in lites
2078 : // totsol total solar transmittance
2079 : // rtot total thermal resistance of window
2080 : // rs array of thermal resistances of each gap and layer
2081 : // layer number of layers
2082 : // dir direct solar radiation
2083 : // Outputs:
2084 : // sf solar gain of space
2085 :
2086 : // Argument array dimensioning
2087 0 : EP_SIZE_CHECK(rs, maxlay3);
2088 0 : EP_SIZE_CHECK(absol, maxlay);
2089 :
2090 : // Locals
2091 : Real64 flowin;
2092 : Real64 fract;
2093 : int i;
2094 : int j;
2095 :
2096 0 : fract = 0.0;
2097 0 : flowin = 0.0;
2098 0 : sf = 0.0;
2099 :
2100 0 : if (rtot == 0.0) {
2101 0 : return;
2102 : }
2103 :
2104 : // evaluate inward flowing fraction of absorbed radiation:
2105 0 : flowin = (rs(1) + 0.5 * rs(2)) / rtot;
2106 0 : fract = absol(1) * flowin;
2107 :
2108 0 : for (i = 2; i <= nlayer; ++i) {
2109 0 : j = 2 * i;
2110 0 : flowin += (0.5 * (rs(j - 2) + rs(j)) + rs(j - 1)) / rtot;
2111 0 : fract += absol(i) * flowin;
2112 : }
2113 0 : sf = totsol + fract; // add inward fraction to directly transmitted fraction
2114 : }
2115 :
2116 0 : void resist(int const nlayer,
2117 : Real64 const trmout,
2118 : Real64 const Tout,
2119 : Real64 const trmin,
2120 : Real64 const tind,
2121 : const Array1D<Real64> &hcgas,
2122 : const Array1D<Real64> &hrgas,
2123 : Array1D<Real64> const &Theta,
2124 : Array1D<Real64> &qlayer,
2125 : const Array1D<Real64> &qv,
2126 : const Array1D<TARCOGLayerType> &LayerType,
2127 : const Array1D<Real64> &thick,
2128 : const Array1D<Real64> &scon,
2129 : Real64 &ufactor,
2130 : Real64 &flux,
2131 : Array1D<Real64> &qcgas,
2132 : Array1D<Real64> &qrgas)
2133 : {
2134 : //***********************************************************************
2135 : // subroutine to calculate total thermal resistance of the glazing system
2136 : //***********************************************************************
2137 :
2138 : // Locals
2139 : int i;
2140 :
2141 : // Simon: calculation of heat flow through gaps and layers as well as ventilation speed and heat flow
2142 : // are kept just for reporting purposes. U-factor calculation is performed by calculating heat flow transfer
2143 : // at indoor layer
2144 :
2145 : // calculate heat flow for external and internal environments and gaps
2146 0 : for (i = 1; i <= nlayer + 1; ++i) {
2147 0 : if (i == 1) {
2148 0 : qcgas(i) = hcgas(i) * (Theta(2 * i - 1) - Tout);
2149 0 : qrgas(i) = hrgas(i) * (Theta(2 * i - 1) - trmout);
2150 0 : qlayer(2 * i - 1) = qcgas(i) + qrgas(i);
2151 : // rs(2*i-1) = 1/hgas(i)
2152 0 : } else if (i == (nlayer + 1)) {
2153 0 : qcgas(i) = hcgas(i) * (tind - Theta(2 * i - 2));
2154 0 : qrgas(i) = hrgas(i) * (trmin - Theta(2 * i - 2));
2155 0 : qlayer(2 * i - 1) = qcgas(i) + qrgas(i);
2156 : // rs(2*i-1) = 1/hgas(i)
2157 : } else {
2158 0 : qcgas(i) = hcgas(i) * (Theta(2 * i - 1) - Theta(2 * i - 2));
2159 0 : qrgas(i) = hrgas(i) * (Theta(2 * i - 1) - Theta(2 * i - 2));
2160 0 : qlayer(2 * i - 1) = qcgas(i) + qrgas(i);
2161 : // rs(2*i-1) = 1/hgas(i)
2162 : }
2163 : }
2164 :
2165 : //.....Calculate thermal resistances for glazing layers:
2166 0 : for (i = 1; i <= nlayer; ++i) {
2167 : // rs(2*i) = thick(i)/scon(i)
2168 0 : qlayer(2 * i) = scon(i) / thick(i) * (Theta(2 * i) - Theta(2 * i - 1));
2169 : }
2170 :
2171 0 : flux = qlayer(2 * nlayer + 1);
2172 0 : if (IsShadingLayer(LayerType(nlayer))) {
2173 0 : flux += qv(nlayer);
2174 : }
2175 :
2176 0 : ufactor = 0.0;
2177 0 : if (tind != Tout) {
2178 0 : ufactor = flux / (tind - Tout);
2179 : }
2180 0 : }
2181 :
2182 0 : void hatter(EnergyPlusData &state,
2183 : int const nlayer,
2184 : int const iwd,
2185 : Real64 const tout,
2186 : Real64 const tind,
2187 : Real64 const wso,
2188 : Real64 const wsi,
2189 : Real64 const VacuumPressure,
2190 : Real64 const VacuumMaxGapThickness,
2191 : Real64 const ebsky,
2192 : Real64 &tamb,
2193 : Real64 const ebroom,
2194 : Real64 &troom,
2195 : const Array1D<Real64> &gap,
2196 : Real64 const height,
2197 : Real64 const heightt,
2198 : const Array1D<Real64> &scon,
2199 : Real64 const tilt,
2200 : Array1D<Real64> const &theta,
2201 : const Array1D<Real64> &Tgap,
2202 : Array1D<Real64> const &Radiation,
2203 : Real64 const trmout,
2204 : Real64 const trmin,
2205 : Array2_int const &iprop,
2206 : Array2<Real64> const &frct,
2207 : const Array1D<Real64> &presure,
2208 : const Array1D_int &nmix,
2209 : const Array1D<Real64> &wght,
2210 : Array2<Real64> const &gcon,
2211 : Array2<Real64> const &gvis,
2212 : Array2<Real64> const &gcp,
2213 : const Array1D<Real64> &gama,
2214 : const Array1D_int &SupportPillar,
2215 : const Array1D<Real64> &PillarSpacing,
2216 : const Array1D<Real64> &PillarRadius,
2217 : Array1D<Real64> &hgas,
2218 : Array1D<Real64> &hcgas,
2219 : Array1D<Real64> &hrgas,
2220 : Real64 &hcin,
2221 : Real64 &hcout,
2222 : Real64 const hin,
2223 : Real64 const hout,
2224 : int const index,
2225 : const Array1D_int &ibc,
2226 : int &nperr,
2227 : std::string &ErrorMessage,
2228 : Real64 &hrin,
2229 : Real64 &hrout,
2230 : Array1D<Real64> &Ra,
2231 : Array1D<Real64> &Nu)
2232 : {
2233 : //***********************************************************************
2234 : // This subroutine calculates the array of conductances/film coefficients used to model convection. The conductances/film
2235 : // coefficients are calculated as functions of temperature defined with the usual variable h and THEN are converted into an
2236 : // equivalent value interms of the black body emittance based on the surface
2237 : //***********************************************************************
2238 : // Inputs
2239 : // nlayer number of solid layers
2240 : // iwd wind direction
2241 : // tout outside temp in k
2242 : // tind inside temp in k
2243 : // wso wind speed in m/s
2244 : // wsi inside forced air speed m/s
2245 : // Ebsky ir flux from outside
2246 : // Ebroom ir flux from room
2247 : // Gout radiosity (ir flux) of the combined environment (sky+ground)
2248 : // Gin
2249 : // gap vector of gap widths in meters
2250 : // height IGU cavity height
2251 : // heightt
2252 : // thick glazing layer thickness
2253 : // scon Vector of conductivities of each glazing layer
2254 : // tilt Window tilt (in degrees)
2255 : // theta Vector of average temperatures
2256 : // Ebb
2257 : // Ebf
2258 : // iprop array of gap mixtures
2259 : // frct vector of mixture fractions
2260 : // presure
2261 : // hin Indoor Indoor combined film coefficient (if non-zero)
2262 : // hout Outdoor combined film coefficient (if non-zero)
2263 : // nmix vector of number of gasses in a mixture for each gap
2264 : // Ouputs
2265 : // hhat vector of all film coefficients (maxlay3)
2266 : // hgas vector of gap 'film' coeff.
2267 : // hcin Indoor convective surface heat transfer coefficient
2268 : // hcout Outdoor convective heat transfer coeff
2269 : // hrin Indoor radiative surface heat transfer coefficient
2270 : // hrout Outdoor radiative surface heat transfer coefficient
2271 : // hin Indoor combined film coefficient
2272 : // hout Outdoor combined film coefficient
2273 : // index iteration step
2274 : // ibc
2275 : // Inactives**
2276 : // wa - window azimuth (degrees, clockwise from south)
2277 :
2278 : // Locals
2279 : int i;
2280 : int k;
2281 : int nface;
2282 :
2283 : // evaluate convective/conductive components of gap grashof number, thermal conductivity and their derivatives:
2284 0 : nface = 2 * nlayer;
2285 :
2286 0 : filmg(state,
2287 : tilt,
2288 : theta,
2289 : Tgap,
2290 : nlayer,
2291 : height,
2292 : gap,
2293 : iprop,
2294 : frct,
2295 : VacuumPressure,
2296 : presure,
2297 : nmix,
2298 : wght,
2299 : gcon,
2300 : gvis,
2301 : gcp,
2302 : gama,
2303 : hcgas,
2304 : Ra,
2305 : Nu,
2306 : nperr,
2307 : ErrorMessage);
2308 :
2309 0 : if (!(GoAhead(nperr))) {
2310 0 : return;
2311 : }
2312 :
2313 : // this is adding influence of pillar to hgas
2314 0 : filmPillar(state, SupportPillar, scon, PillarSpacing, PillarRadius, nlayer, gap, hcgas, VacuumMaxGapThickness, nperr, ErrorMessage);
2315 :
2316 0 : if (!(GoAhead(nperr))) {
2317 0 : return;
2318 : }
2319 :
2320 : // adjust radiation coefficients
2321 : // hrgas = 0.0d0
2322 0 : for (i = 2; i <= nlayer; ++i) {
2323 0 : k = 2 * i - 1;
2324 : // if ((theta(k)-theta(k-1)) == 0) then
2325 : // theta(k-1) = theta(k-1) + tempCorrection
2326 : // end if
2327 0 : if ((theta(k) - theta(k - 1)) != 0) {
2328 0 : hrgas(i) = (Radiation(k) - Radiation(k - 1)) / (theta(k) - theta(k - 1));
2329 : }
2330 :
2331 0 : hgas(i) = hcgas(i) + hrgas(i);
2332 : }
2333 :
2334 : // convective indoor film coeff:
2335 0 : if (ibc(2) <= 0) {
2336 0 : filmi(state,
2337 : tind,
2338 0 : theta(nface),
2339 : nlayer,
2340 : tilt,
2341 : wsi,
2342 : heightt,
2343 : iprop,
2344 : frct,
2345 : presure,
2346 : nmix,
2347 : wght,
2348 : gcon,
2349 : gvis,
2350 : gcp,
2351 : hcin,
2352 0 : ibc(2),
2353 : nperr,
2354 : ErrorMessage);
2355 0 : } else if (ibc(2) == 1) {
2356 0 : hcin = hin - hrin;
2357 : // Simon: First iteration is with index = 0 and that means it should reenter iteration with whatever is provided as input
2358 : // else if (ibc(2).eq.2.and.index.eq.1) then
2359 0 : } else if ((ibc(2) == 2) && (index == 0)) {
2360 0 : hcin = hin;
2361 : }
2362 0 : if (hcin < 0) {
2363 0 : nperr = 8;
2364 0 : ErrorMessage = "Hcin is less then zero.";
2365 0 : return;
2366 : }
2367 :
2368 0 : hcgas(nlayer + 1) = hcin;
2369 : // hrin = 0.95d0*(Ebroom - Radiation(2*nlayer))/(Trmin-theta(2*nlayer))+0.05d0*hrin
2370 0 : hrin = (ebroom - Radiation(2 * nlayer)) / (trmin - theta(2 * nlayer));
2371 : // if ((Theta(2*nlayer) - Trmin).ne.0) then
2372 : // hrin = sigma * emis(2*nlayer) * (Theta(2*nlayer)**4 - Trmin**4)/(Theta(2*nlayer) - Trmin)
2373 : // else
2374 : // Theta(2*nlayer) = Theta(2*nlayer) + tempCorrection
2375 : // hrin = sigma * emis(2*nlayer) * (Theta(2*nlayer)**4 - Trmin**4)/(Theta(2*nlayer) - Trmin)
2376 : // end if
2377 0 : hrgas(nlayer + 1) = hrin;
2378 : // hgas(nlayer+1) = hcgas(nlayer+1) + hrgas(nlayer+1)
2379 0 : troom = (hcin * tind + hrin * trmin) / (hcin + hrin);
2380 :
2381 : // convective outdoor film coeff:
2382 0 : if (ibc(1) <= 0) {
2383 0 : film(tout, theta(1), wso, iwd, hcout, ibc(1));
2384 0 : } else if (ibc(1) == 1) {
2385 0 : hcout = hout - hrout;
2386 : // Simon: First iteration is with index = 0 and that means it should reenter iteration with whatever is provided as input
2387 : // else if (ibc(1).eq.2.and.index.eq.1) then
2388 0 : } else if ((ibc(1) == 2) && (index == 0)) {
2389 0 : hcout = hout;
2390 : }
2391 0 : if (hcout < 0) {
2392 0 : nperr = 9;
2393 0 : ErrorMessage = "Hcout is less than zero.";
2394 0 : return;
2395 : }
2396 :
2397 0 : hcgas(1) = hcout;
2398 0 : hrout = (Radiation(1) - ebsky) / (theta(1) - trmout);
2399 : // if ((Theta(1) - Trmout).ne.0) then
2400 : // hrout = sigma * emis(1) * (Theta(1)**4 - Trmout**4)/(Theta(1) - Trmout)
2401 : // else
2402 : // Theta(1) = Theta(1) + tempCorrection
2403 : // hrout = sigma * emis(1) * (Theta(1)**4 - Trmout**4)/(Theta(1) - Trmout)
2404 : // end if
2405 0 : hrgas(1) = hrout;
2406 : // hgas(1) = hrout + hcout
2407 0 : tamb = (hcout * tout + hrout * trmout) / (hcout + hrout);
2408 : }
2409 :
2410 0 : void effectiveLayerCond(EnergyPlusData &state,
2411 : int const nlayer,
2412 : const Array1D<TARCOGLayerType> &LayerType, // Layer type
2413 : const Array1D<Real64> &scon, // Layer thermal conductivity
2414 : const Array1D<Real64> &thick, // Layer thickness
2415 : Array2A_int const iprop, // Gas type in gaps
2416 : Array2A<Real64> const frct, // Fraction of gas
2417 : const Array1D_int &nmix, // Gas mixture
2418 : const Array1D<Real64> &pressure, // Gas pressure [Pa]
2419 : const Array1D<Real64> &wght, // Molecular weight
2420 : Array2A<Real64> const gcon, // Gas specific conductivity
2421 : Array2A<Real64> const gvis, // Gas specific viscosity
2422 : Array2A<Real64> const gcp, // Gas specific heat
2423 : const Array1D<Real64> &EffectiveOpenness, // Layer effective openness [m2]
2424 : Array1D<Real64> const &theta, // Layer surface temperatures [K]
2425 : Array1D<Real64> &sconScaled, // Layer conductivity divided by thickness
2426 : int &nperr, // Error message flag
2427 : std::string &ErrorMessage // Error message
2428 : )
2429 : {
2430 0 : for (int i = 1; i <= nlayer; ++i) {
2431 0 : if (LayerType(i) != TARCOGLayerType::SPECULAR) {
2432 0 : Real64 tLayer = (theta(2 * i - 1) + theta(2 * i)) / 2;
2433 0 : Real64 nmix1 = nmix(i);
2434 0 : Real64 press1 = (pressure(i) + pressure(i + 1)) / 2.0;
2435 0 : for (int j = 1; j <= maxgas; ++j) {
2436 0 : state.dataThermalISO15099Calc->iprop1(j) = iprop(j, i);
2437 0 : state.dataThermalISO15099Calc->frct1(j) = frct(j, i);
2438 : }
2439 :
2440 : Real64 con;
2441 : Real64 visc;
2442 : Real64 dens;
2443 : Real64 cp;
2444 : Real64 pr;
2445 0 : GASSES90(state,
2446 : tLayer,
2447 0 : state.dataThermalISO15099Calc->iprop1,
2448 0 : state.dataThermalISO15099Calc->frct1,
2449 : press1,
2450 : nmix1,
2451 : wght,
2452 : gcon,
2453 : gvis,
2454 : gcp,
2455 : con,
2456 : visc,
2457 : dens,
2458 : cp,
2459 : pr,
2460 : TARCOGGassesParams::Stdrd::ISO15099,
2461 : nperr,
2462 : ErrorMessage);
2463 0 : sconScaled(i) = (EffectiveOpenness(i) * con + (1 - EffectiveOpenness(i)) * scon(i)) / thick(i);
2464 : } else {
2465 0 : sconScaled(i) = scon(i) / thick(i);
2466 : }
2467 : }
2468 0 : }
2469 :
2470 0 : void filmi(EnergyPlusData &state,
2471 : Real64 const tair,
2472 : Real64 const t,
2473 : int const nlayer,
2474 : Real64 const tilt,
2475 : Real64 const wsi,
2476 : Real64 const height,
2477 : Array2A_int const iprop,
2478 : Array2A<Real64> const frct,
2479 : const Array1D<Real64> &presure,
2480 : const Array1D_int &nmix,
2481 : const Array1D<Real64> &wght,
2482 : Array2A<Real64> const gcon,
2483 : Array2A<Real64> const gvis,
2484 : Array2A<Real64> const gcp,
2485 : Real64 &hcin,
2486 : int const ibc,
2487 : int &nperr,
2488 : std::string &ErrorMessage)
2489 : {
2490 : //***********************************************************************
2491 : // purpose to evaluate heat flux at indoor surface of window using still air correlations (Curcija and Goss 1993)
2492 : // found in SPC142 equations 5.43 - 5.48.
2493 : //***********************************************************************
2494 : // Input
2495 : // tair - room air temperature
2496 : // t - inside surface temperature
2497 : // nlayer number of glazing layers
2498 : // tilt - the tilt of the glazing in degrees
2499 : // wsi - room wind speed (m/s)
2500 : // height - window height
2501 : // iprop
2502 : // frct
2503 : // presure
2504 : // nmix vector of number of gasses in a mixture for each gap
2505 : // Output
2506 : // hcin - indoor convective heat transfer coeff
2507 :
2508 : // If there is forced air in the room than use SPC142 correlation 5.49 to calculate the room side film coefficient.
2509 :
2510 : // Using
2511 : // Argument array dimensioning
2512 0 : iprop.dim(maxgas, maxlay1);
2513 0 : frct.dim(maxgas, maxlay1);
2514 0 : EP_SIZE_CHECK(presure, maxlay1);
2515 0 : EP_SIZE_CHECK(nmix, maxlay1);
2516 0 : EP_SIZE_CHECK(wght, maxgas);
2517 0 : gcon.dim(3, maxgas);
2518 0 : gvis.dim(3, maxgas);
2519 0 : gcp.dim(3, maxgas);
2520 :
2521 : // Locals
2522 : int j;
2523 : Real64 tiltr;
2524 : Real64 tmean;
2525 : Real64 delt;
2526 : Real64 con;
2527 : Real64 visc;
2528 : Real64 dens;
2529 : Real64 cp;
2530 : Real64 pr;
2531 : Real64 gr;
2532 : Real64 RaCrit;
2533 : Real64 RaL;
2534 0 : Real64 Gnui(0.0);
2535 :
2536 0 : if (wsi > 0.0) { // main IF
2537 0 : switch (ibc) {
2538 0 : case 0: {
2539 0 : hcin = 4.0 + 4.0 * wsi;
2540 0 : } break;
2541 0 : case -1: {
2542 0 : hcin = 5.6 + 3.8 * wsi; // SPC142 correlation
2543 0 : return;
2544 : } break;
2545 0 : default:
2546 0 : break;
2547 : }
2548 : } else { // main IF - else
2549 0 : tiltr = tilt * 2.0 * Constant::Pi / 360.0; // convert tilt in degrees to radians
2550 0 : tmean = tair + 0.25 * (t - tair);
2551 0 : delt = std::abs(tair - t);
2552 :
2553 0 : for (j = 1; j <= nmix(nlayer + 1); ++j) {
2554 0 : state.dataThermalISO15099Calc->ipropi(j) = iprop(j, nlayer + 1);
2555 0 : state.dataThermalISO15099Calc->frcti(j) = frct(j, nlayer + 1);
2556 : }
2557 :
2558 0 : GASSES90(state,
2559 : tmean,
2560 0 : state.dataThermalISO15099Calc->ipropi,
2561 0 : state.dataThermalISO15099Calc->frcti,
2562 0 : presure(nlayer + 1),
2563 0 : nmix(nlayer + 1),
2564 : wght,
2565 : gcon,
2566 : gvis,
2567 : gcp,
2568 : con,
2569 : visc,
2570 : dens,
2571 : cp,
2572 : pr,
2573 : TARCOGGassesParams::Stdrd::ISO15099,
2574 : nperr,
2575 : ErrorMessage);
2576 :
2577 : // Calculate grashoff number:
2578 : // The grashoff number is the Rayleigh Number (equation 5.29) in SPC142 divided by the Prandtl Number (prand):
2579 0 : gr = Constant::Gravity * pow_3(height) * delt * pow_2(dens) / (tmean * pow_2(visc));
2580 :
2581 0 : RaL = gr * pr;
2582 : // write(*,*)' RaCrit,RaL,gr,pr '
2583 : // write(*,*) RaCrit,RaL,gr,pr
2584 :
2585 0 : if ((0.0 <= tilt) && (tilt < 15.0)) { // IF no. 1
2586 0 : Gnui = 0.13 * std::pow(RaL, 1.0 / 3.0);
2587 0 : } else if ((15.0 <= tilt) && (tilt <= 90.0)) {
2588 : // if the room air is still THEN use equations 5.43 - 5.48:
2589 0 : RaCrit = 2.5e5 * std::pow(std::exp(0.72 * tilt) / std::sin(tiltr), 0.2);
2590 0 : if (RaL <= RaCrit) { // IF no. 2
2591 0 : Gnui = 0.56 * root_4(RaL * std::sin(tiltr));
2592 : // write(*,*) ' Nu ', Gnui
2593 : } else {
2594 : // Gnui = 0.13d0*(RaL**0.3333d0 - RaCrit**0.3333d0) + 0.56d0*(RaCrit*sin(tiltr))**0.25d0
2595 0 : Gnui = 0.13 * (std::pow(RaL, 1.0 / 3.0) - std::pow(RaCrit, 1.0 / 3.0)) + 0.56 * root_4(RaCrit * std::sin(tiltr));
2596 : } // end if no. 2
2597 0 : } else if ((90.0 < tilt) && (tilt <= 179.0)) {
2598 0 : Gnui = 0.56 * root_4(RaL * std::sin(tiltr));
2599 0 : } else if ((179.0 < tilt) && (tilt <= 180.0)) {
2600 0 : Gnui = 0.58 * std::pow(RaL, 1 / 3.0);
2601 : } else {
2602 0 : assert(false);
2603 : } // end if no. 1
2604 : // write(*,*) ' RaL ', RaL, ' RaCrit', RaCrit
2605 : // write(*,*)' Nusselt Number ',Gnui
2606 :
2607 0 : hcin = Gnui * (con / height);
2608 : // hin = 1.77d0*(ABS(t-tair))**0.25d0
2609 :
2610 : } // end main IF
2611 : }
2612 :
2613 0 : void filmg(EnergyPlusData &state,
2614 : Real64 const tilt,
2615 : const Array1D<Real64> &theta,
2616 : const Array1D<Real64> &Tgap,
2617 : int const nlayer,
2618 : Real64 const height,
2619 : const Array1D<Real64> &gap,
2620 : Array2A_int const iprop,
2621 : Array2A<Real64> const frct,
2622 : Real64 const VacuumPressure,
2623 : const Array1D<Real64> &presure,
2624 : const Array1D_int &nmix,
2625 : const Array1D<Real64> &wght,
2626 : Array2A<Real64> const gcon,
2627 : Array2A<Real64> const gvis,
2628 : Array2A<Real64> const gcp,
2629 : const Array1D<Real64> &gama,
2630 : Array1D<Real64> &hcgas,
2631 : Array1D<Real64> &Rayleigh,
2632 : Array1D<Real64> &Nu,
2633 : int &nperr,
2634 : std::string &ErrorMessage)
2635 : {
2636 : //***********************************************************************
2637 : // subroutine to calculate effective conductance of gaps
2638 : //***********************************************************************
2639 : // Inputs:
2640 : // tilt window angle (deg)
2641 : // theta vector of surface temperatures [K]
2642 : // nlayer total number of glazing layers
2643 : // height glazing cavity height
2644 : // gap vector of gap widths [m]
2645 : // iprop
2646 : // frct
2647 : // presure
2648 : // nmix vector of number of gasses in a mixture for each gap
2649 : // Output:
2650 : // hgas vector of gap coefficients
2651 : // nperr error code
2652 : // Locals:
2653 : // gr gap grashof number
2654 : // con gap gas conductivity
2655 : // visc dynamic viscosity @ mean temperature [g/m*s]
2656 : // dens density @ mean temperature [kg/m^3]
2657 : // cp specific heat @ mean temperature [J/g*K]
2658 : // pr gap gas Prandtl number
2659 : // tmean average film temperature
2660 : // delt temperature difference
2661 :
2662 : // Using
2663 : // Argument array dimensioning
2664 0 : EP_SIZE_CHECK(theta, maxlay2);
2665 0 : EP_SIZE_CHECK(Tgap, maxlay1);
2666 0 : EP_SIZE_CHECK(gap, MaxGap);
2667 0 : iprop.dim(maxgas, maxlay1);
2668 0 : frct.dim(maxgas, maxlay1);
2669 0 : EP_SIZE_CHECK(presure, maxlay1);
2670 0 : EP_SIZE_CHECK(nmix, maxlay1);
2671 0 : EP_SIZE_CHECK(wght, maxgas);
2672 0 : gcon.dim(3, maxgas);
2673 0 : gvis.dim(3, maxgas);
2674 0 : gcp.dim(3, maxgas);
2675 0 : EP_SIZE_CHECK(gama, maxgas);
2676 0 : EP_SIZE_CHECK(hcgas, maxlay1);
2677 0 : EP_SIZE_CHECK(Rayleigh, maxlay);
2678 0 : EP_SIZE_CHECK(Nu, maxlay);
2679 :
2680 : // Locals
2681 : Real64 con;
2682 : Real64 visc;
2683 : Real64 dens;
2684 : Real64 cp;
2685 : Real64 pr;
2686 : Real64 delt;
2687 : Real64 tmean;
2688 : Real64 ra;
2689 : Real64 asp;
2690 : Real64 gnu;
2691 : int i;
2692 : int j;
2693 : int k;
2694 : int l;
2695 :
2696 0 : hcgas = 0.0;
2697 :
2698 0 : for (i = 1; i <= nlayer - 1; ++i) {
2699 0 : j = 2 * i;
2700 0 : k = j + 1;
2701 : // determine the gas properties of each gap:
2702 : // tmean = (theta(j)+theta(k))/2.0d0
2703 0 : tmean = Tgap(i + 1); // Tgap(1) is exterior environment
2704 0 : delt = std::abs(theta(j) - theta(k));
2705 : // Temperatures should not be equal. This can happen in initial temperature guess before iterations started
2706 0 : if (delt == 0.0) {
2707 0 : delt = 1.0e-6;
2708 : }
2709 0 : for (l = 1; l <= nmix(i + 1); ++l) {
2710 0 : state.dataThermalISO15099Calc->ipropg(l) = iprop(l, i + 1);
2711 0 : state.dataThermalISO15099Calc->frctg(l) = frct(l, i + 1);
2712 : }
2713 :
2714 0 : if (presure(i + 1) > VacuumPressure) {
2715 0 : GASSES90(state,
2716 : tmean,
2717 0 : state.dataThermalISO15099Calc->ipropg,
2718 0 : state.dataThermalISO15099Calc->frctg,
2719 0 : presure(i + 1),
2720 0 : nmix(i + 1),
2721 : wght,
2722 : gcon,
2723 : gvis,
2724 : gcp,
2725 : con,
2726 : visc,
2727 : dens,
2728 : cp,
2729 : pr,
2730 : TARCOGGassesParams::Stdrd::ISO15099,
2731 : nperr,
2732 : ErrorMessage);
2733 :
2734 : // Calculate grashoff number:
2735 : // The grashoff number is the Rayleigh Number (equation 5.29) in SPC142 divided by the Prandtl Number (prand):
2736 0 : ra = Constant::Gravity * pow_3(gap(i)) * delt * cp * pow_2(dens) / (tmean * visc * con);
2737 0 : Rayleigh(i) = ra;
2738 : // write(*,*) 'height,gap(i),asp',height,gap(i),asp
2739 : // asp = 1
2740 : // if (gap(i).ne.0) then
2741 0 : asp = height / gap(i);
2742 : // end if
2743 : // determine the Nusselt number:
2744 0 : nusselt(tilt, ra, asp, gnu, nperr, ErrorMessage);
2745 :
2746 0 : Nu(i) = gnu;
2747 : // calculate effective conductance of the gap
2748 0 : hcgas(i + 1) = con / gap(i) * gnu;
2749 :
2750 : // write(*,*)'theta(j),theta(k),j,k',j,theta(j),k,theta(k)
2751 : // write(*,*)'Nusselt,Rayleigh,Prandtl,hgas(k),k'
2752 : // write(*,*) gnu,gr*pr,pr,hgas(k),k
2753 : } else { // low pressure calculations
2754 0 : GassesLow(tmean, wght(iprop(1, i + 1)), presure(i + 1), gama(iprop(1, i + 1)), con, nperr, ErrorMessage);
2755 0 : hcgas(i + 1) = con;
2756 : } // if (pressure(i+1).gt.VacuumPressure) then
2757 : }
2758 0 : }
2759 :
2760 0 : void filmPillar(EnergyPlusData &state,
2761 : const Array1D_int &SupportPillar, // Shows whether or not gap have support pillar
2762 : const Array1D<Real64> &scon, // Conductivity of glass layers
2763 : const Array1D<Real64> &PillarSpacing, // Pillar spacing for each gap (used in case there is support pillar)
2764 : const Array1D<Real64> &PillarRadius, // Pillar radius for each gap (used in case there is support pillar)
2765 : int const nlayer,
2766 : const Array1D<Real64> &gap,
2767 : Array1D<Real64> &hcgas,
2768 : [[maybe_unused]] Real64 const VacuumMaxGapThickness,
2769 : [[maybe_unused]] int &nperr,
2770 : [[maybe_unused]] std::string &ErrorMessage)
2771 : {
2772 : //***********************************************************************
2773 : // subroutine to calculate effective conductance of support pillars
2774 : //***********************************************************************
2775 :
2776 : // Using
2777 : // Argument array dimensioning
2778 0 : EP_SIZE_CHECK(SupportPillar, maxlay);
2779 0 : EP_SIZE_CHECK(scon, maxlay);
2780 0 : EP_SIZE_CHECK(PillarSpacing, maxlay);
2781 0 : EP_SIZE_CHECK(PillarRadius, maxlay);
2782 0 : EP_SIZE_CHECK(gap, MaxGap);
2783 0 : EP_SIZE_CHECK(hcgas, maxlay1);
2784 :
2785 : // Locals
2786 : // 0 - does not have support pillar
2787 : // 1 - have support pillar
2788 :
2789 0 : for (state.dataThermalISO15099Calc->iFP = 1; state.dataThermalISO15099Calc->iFP <= nlayer - 1; ++state.dataThermalISO15099Calc->iFP) {
2790 0 : state.dataThermalISO15099Calc->kFP = 2 * state.dataThermalISO15099Calc->iFP + 1;
2791 0 : if (SupportPillar(state.dataThermalISO15099Calc->iFP) == YES_SupportPillar) {
2792 :
2793 : // Average glass conductivity is taken as average from both glass surrounding gap
2794 0 : state.dataThermalISO15099Calc->aveGlassConductivity =
2795 0 : (scon(state.dataThermalISO15099Calc->iFP) + scon(state.dataThermalISO15099Calc->iFP + 1)) / 2;
2796 :
2797 0 : state.dataThermalISO15099Calc->cpa =
2798 0 : 2.0 * state.dataThermalISO15099Calc->aveGlassConductivity * PillarRadius(state.dataThermalISO15099Calc->iFP) /
2799 0 : (pow_2(PillarSpacing(state.dataThermalISO15099Calc->iFP)) *
2800 0 : (1.0 + 2.0 * gap(state.dataThermalISO15099Calc->iFP) / (Constant::Pi * PillarRadius(state.dataThermalISO15099Calc->iFP))));
2801 :
2802 : // It is important to add on previous values calculated for gas
2803 0 : hcgas(state.dataThermalISO15099Calc->iFP + 1) += state.dataThermalISO15099Calc->cpa;
2804 : } // if (SupportPillar(i).eq.YES_SupportPillar) then
2805 : }
2806 0 : }
2807 :
2808 0 : void nusselt(Real64 const tilt, Real64 const ra, Real64 const asp, Real64 &gnu, int &nperr, std::string &ErrorMessage)
2809 : {
2810 : //***********************************************************************
2811 : // purpose to calculate nusselt modulus for air gaps (ISO15099)
2812 : //***********************************************************************
2813 : // Input
2814 : // tilt tilt in degrees
2815 : // ra rayleigh number
2816 : // asp Aspect ratio
2817 : // Output
2818 : // gnu nusselt number
2819 : // nperr
2820 :
2821 : // Using
2822 : // Locals
2823 : Real64 subNu1;
2824 : Real64 subNu2;
2825 : Real64 subNu3;
2826 : Real64 Nu1;
2827 : Real64 Nu2;
2828 : Real64 G;
2829 : Real64 Nu60;
2830 : Real64 Nu90;
2831 : Real64 tiltr;
2832 :
2833 0 : subNu1 = 0.0;
2834 0 : subNu2 = 0.0;
2835 0 : subNu3 = 0.0;
2836 0 : Nu1 = 0.0;
2837 0 : Nu2 = 0.0;
2838 0 : Nu90 = 0.0;
2839 0 : Nu60 = 0.0;
2840 0 : G = 0.0;
2841 0 : tiltr = tilt * 2.0 * Constant::Pi / 360.0; // convert tilt in degrees to radians
2842 0 : if ((tilt >= 0.0) && (tilt < 60.0)) { // ISO/DIS 15099 - chapter 5.3.3.1
2843 0 : subNu1 = 1.0 - 1708.0 / (ra * std::cos(tiltr));
2844 0 : subNu1 = pos(subNu1);
2845 0 : subNu2 = 1.0 - (1708.0 * std::pow(std::sin(1.8 * tiltr), 1.6)) / (ra * std::cos(tiltr));
2846 0 : subNu3 = std::pow(ra * std::cos(tiltr) / 5830.0, 1.0 / 3.0) - 1.0;
2847 0 : subNu3 = pos(subNu3);
2848 0 : gnu = 1.0 + 1.44 * subNu1 * subNu2 + subNu3; // equation 42
2849 0 : if (ra >= 1.0e5) {
2850 0 : nperr = 1001; // Rayleigh number is out of range
2851 0 : ErrorMessage = "Rayleigh number out of range in Nusselt num. calc. for gaps (angle between 0 and 60 deg).";
2852 : }
2853 0 : if (asp <= 20.0) {
2854 0 : nperr = 1002; // Aspect Ratio is out of range
2855 0 : ErrorMessage = "Aspect Ratio out of range in Nusselt num. calc. for gaps (angle between 0 and 60 deg).";
2856 : }
2857 0 : } else if (tilt == 60.0) { // ISO/DIS 15099 - chapter 5.3.3.2
2858 0 : G = 0.5 / std::pow(1.0 + std::pow(ra / 3160.0, 20.6), 0.1); // equation 47
2859 0 : Nu1 = std::pow(1.0 + pow_7((0.0936 * std::pow(ra, 0.314)) / (1.0 + G)), 0.1428571); // equation 45
2860 0 : Nu2 = (0.104 + 0.175 / asp) * std::pow(ra, 0.283); // equation 46
2861 0 : gnu = max(Nu1, Nu2); // equation 44
2862 0 : } else if ((tilt > 60.0) && (tilt < 90.0)) { // ISO/DIS 15099 - chapter 5.3.3.3
2863 0 : if ((ra > 100.0) && (ra < 2.0e7) && (asp > 5.0) && (asp < 100.0)) {
2864 0 : G = 0.5 / std::pow(1.0 + std::pow(ra / 3160.0, 20.6), 0.1); // equation 47
2865 0 : Nu1 = std::pow(1.0 + pow_7((0.0936 * std::pow(ra, 0.314)) / (1.0 + G)), 0.1428571); // equation 45
2866 0 : Nu2 = (0.104 + 0.175 / asp) * std::pow(ra, 0.283); // equation 46
2867 0 : Nu60 = max(Nu1, Nu2); // equation 44
2868 0 : Nu2 = 0.242 * std::pow(ra / asp, 0.272); // equation 52
2869 0 : if (ra > 5.0e4) {
2870 0 : Nu1 = 0.0673838 * std::pow(ra, 1.0 / 3.0); // equation 49
2871 0 : } else if ((ra > 1.0e4) && (ra <= 5.0e4)) {
2872 0 : Nu1 = 0.028154 * std::pow(ra, 0.4134); // equation 50
2873 0 : } else if (ra <= 1.0e4) {
2874 0 : Nu1 = 1.0 + 1.7596678e-10 * std::pow(ra, 2.2984755); // equation 51
2875 : }
2876 0 : } else if (ra <= 100.0) {
2877 0 : G = 0.5 / std::pow(1.0 + std::pow(ra / 3160.0, 20.6), 0.1); // equation 47
2878 0 : Nu1 = std::pow(1.0 + pow_7((0.0936 * std::pow(ra, 0.314)) / (1.0 + G)), 0.1428571); // equation 45
2879 0 : Nu2 = (0.104 + 0.175 / asp) * std::pow(ra, 0.283); // equation 46
2880 0 : Nu60 = max(Nu1, Nu2); // equation 44
2881 0 : Nu2 = 0.242 * std::pow(ra / asp, 0.272); // equation 52
2882 0 : Nu1 = 1.0 + 1.7596678e-10 * std::pow(ra, 2.2984755); // equation 51
2883 0 : nperr = 1003; // Rayleigh number is less than 100
2884 0 : ErrorMessage = "Rayleigh number is less than 100 in Nusselt number calculations for gaps (angle between 60 and 90 degrees).";
2885 0 : } else if (ra > 2.0e7) {
2886 0 : G = 0.5 / std::pow(1.0 + std::pow(ra / 3160.0, 20.6), 0.1); // equation 47
2887 0 : Nu1 = std::pow(1.0 + pow_7((0.0936 * std::pow(ra, 0.314)) / (1.0 + G)), 0.1428571); // equation 45
2888 0 : Nu2 = (0.104 + 0.175 / asp) * std::pow(ra, 0.283); // equation 46
2889 0 : Nu60 = max(Nu1, Nu2); // equation 44
2890 0 : Nu2 = 0.242 * std::pow(ra / asp, 0.272); // equation 52
2891 0 : Nu1 = 0.0673838 * std::pow(ra, 1.0 / 3.0); // equation 49
2892 0 : nperr = 1004; // Rayleigh number is great from 2e7
2893 0 : ErrorMessage = "Rayleigh number is greater than 2e7 in Nusselt number calculations for gaps (angle between 60 and 90 degrees).";
2894 0 : } else if ((asp <= 5.0) || (asp >= 100.0)) {
2895 0 : G = 0.5 / std::pow(1.0 + std::pow(ra / 3160.0, 20.6), 0.1); // equation 47
2896 0 : Nu1 = std::pow(1.0 + pow_7((0.0936 * std::pow(ra, 0.314)) / (1.0 + G)), 0.1428571); // equation 45
2897 0 : Nu2 = (0.104 + 0.175 / asp) * std::pow(ra, 0.283); // equation 46
2898 0 : Nu60 = max(Nu1, Nu2); // equation 44
2899 0 : Nu2 = 0.242 * std::pow(ra / asp, 0.272); // equation 52
2900 0 : if (ra > 5.0e4) {
2901 0 : Nu1 = 0.0673838 * std::pow(ra, 1.0 / 3.0); // equation 49
2902 0 : } else if ((ra > 1.0e4) && (ra <= 5.0e4)) {
2903 0 : Nu1 = 0.028154 * std::pow(ra, 0.4134); // equation 50
2904 0 : } else if (ra <= 1.0e4) {
2905 0 : Nu1 = 1.0 + 1.7596678e-10 * std::pow(ra, 2.2984755); // equation 51
2906 : }
2907 0 : nperr = 1005; // Aspect Ratio is out of range
2908 0 : ErrorMessage = "Aspect Ratio is out of range in Nusselt number calculations for gaps (angle between 60 and 90 degrees).";
2909 : }
2910 0 : Nu90 = max(Nu1, Nu2); // equation 48
2911 0 : gnu = ((Nu90 - Nu60) / (90.0 - 60.0)) * (tilt - 60.0) + Nu60; // linear interpolation between 60 and 90 degrees
2912 0 : } else if (tilt == 90.0) { // ISO/DIS 15099 - chapter 5.3.3.4
2913 0 : Nu2 = 0.242 * std::pow(ra / asp, 0.272); // equation 52
2914 0 : if (ra > 5.0e4) {
2915 0 : Nu1 = 0.0673838 * std::pow(ra, 1.0 / 3.0); // equation 49
2916 0 : } else if ((ra > 1.0e4) && (ra <= 5.0e4)) {
2917 0 : Nu1 = 0.028154 * std::pow(ra, 0.4134); // equation 50
2918 : // Nu1 = 0.028154d0 * ra ** 0.414d0 !equation 50 - DISCONTINUITY CORRECTED
2919 0 : } else if (ra <= 1.0e4) {
2920 0 : Nu1 = 1.0 + 1.7596678e-10 * std::pow(ra, 2.2984755); // equation 51
2921 : }
2922 0 : gnu = max(Nu1, Nu2); // equation 48
2923 0 : } else if ((tilt > 90.0) && (tilt <= 180.0)) {
2924 0 : Nu2 = 0.242 * std::pow(ra / asp, 0.272); // equation 52
2925 0 : if (ra > 5.0e4) {
2926 0 : Nu1 = 0.0673838 * std::pow(ra, 1.0 / 3.0); // equation 49
2927 0 : } else if ((ra > 1.0e4) && (ra <= 5.0e4)) {
2928 0 : Nu1 = 0.028154 * std::pow(ra, 0.4134); // equation 50
2929 0 : } else if (ra <= 1.0e4) {
2930 0 : Nu1 = 1.0 + 1.7596678e-10 * std::pow(ra, 2.2984755); // equation 51
2931 : }
2932 0 : gnu = max(Nu1, Nu2); // equation 48
2933 0 : gnu = 1.0 + (gnu - 1.0) * std::sin(tiltr); // equation 53
2934 : } else {
2935 0 : nperr = 10; // error flag: angle is out of range
2936 0 : ErrorMessage = "Window tilt angle is out of range.";
2937 0 : return;
2938 : }
2939 : }
2940 :
2941 0 : void storeIterationResults(EnergyPlusData &state,
2942 : Files &files,
2943 : int const nlayer,
2944 : int const index,
2945 : const Array1D<Real64> &theta,
2946 : Real64 const trmout,
2947 : Real64 const tamb,
2948 : Real64 const trmin,
2949 : Real64 const troom,
2950 : Real64 const ebsky,
2951 : Real64 const ebroom,
2952 : Real64 const hcin,
2953 : Real64 const hcout,
2954 : Real64 const hrin,
2955 : Real64 const hrout,
2956 : Real64 const hin,
2957 : Real64 const hout,
2958 : const Array1D<Real64> &Ebb,
2959 : const Array1D<Real64> &Ebf,
2960 : const Array1D<Real64> &Rb,
2961 : const Array1D<Real64> &Rf,
2962 : int &)
2963 : {
2964 0 : auto &dynFormat = state.dataThermalISO15099Calc->dynFormat;
2965 : int i;
2966 :
2967 : // write(a,1000) index
2968 0 : print(files.TarcogIterationsFile, "*************************************************************************************************\n");
2969 0 : print(files.TarcogIterationsFile, "Iteration number: {:5}\n", index);
2970 :
2971 0 : print(files.TarcogIterationsFile, "Trmin = {:8.4F}\n", trmin - Constant::Kelvin);
2972 0 : print(files.TarcogIterationsFile, "Troom = {:12.6F}\n", troom - Constant::Kelvin);
2973 0 : print(files.TarcogIterationsFile, "Trmout = {:8.4F}\n", trmout - Constant::Kelvin);
2974 0 : print(files.TarcogIterationsFile, "Tamb = {:12.6F}\n", tamb - Constant::Kelvin);
2975 :
2976 0 : print(files.TarcogIterationsFile, "Ebsky = {:8.4F}\n", ebsky);
2977 0 : print(files.TarcogIterationsFile, "Ebroom = {:8.4F}\n", ebroom);
2978 :
2979 0 : print(files.TarcogIterationsFile, "hcin = {:8.4F}\n", hcin);
2980 0 : print(files.TarcogIterationsFile, "hcout = {:8.4F}\n", hcout);
2981 0 : print(files.TarcogIterationsFile, "hrin = {:8.4F}\n", hrin);
2982 0 : print(files.TarcogIterationsFile, "hrout = {:8.4F}\n", hrout);
2983 0 : print(files.TarcogIterationsFile, "hin = {:8.4F}\n", hin);
2984 0 : print(files.TarcogIterationsFile, "hout = {:8.4F}\n", hout);
2985 :
2986 : // Write headers for Ebb and Ebf
2987 0 : for (i = 1; i <= 2 * nlayer; ++i) {
2988 0 : if (i == 1) {
2989 0 : dynFormat = "";
2990 : }
2991 0 : if (mod(i, 2) == 1) {
2992 0 : dynFormat += fmt::format("Ebf({:3})", (i + 1) / 2);
2993 : } else {
2994 0 : dynFormat += fmt::format("Ebb({:3})", (i + 1) / 2);
2995 : }
2996 0 : if (i != 2 * nlayer) {
2997 0 : dynFormat += "===";
2998 : }
2999 : }
3000 0 : print(files.TarcogIterationsFile, dynFormat);
3001 0 : print(files.TarcogIterationsFile, "\n");
3002 :
3003 : // write Ebb and Ebf
3004 0 : print(files.TarcogIterationsFile, "{:16.8F} {:16.8F}", Ebf(1), Ebb(1));
3005 0 : for (i = 2; i <= nlayer; ++i) {
3006 0 : print(files.TarcogIterationsFile, " {:16.8F} {:16.8F}", Ebf(i), Ebb(i));
3007 : }
3008 0 : print(files.TarcogIterationsFile, "\n");
3009 :
3010 : // Write headers for Rb and Rf
3011 0 : for (i = 1; i <= 2 * nlayer; ++i) {
3012 0 : const std::string a = fmt::format("{:3}", (i + 1) / 2); // this is just to simulate correct integer in brackets
3013 0 : if (i == 1) {
3014 0 : dynFormat = "";
3015 : }
3016 0 : if (mod(i, 2) == 1) {
3017 0 : dynFormat += "Rf(" + a + ')';
3018 : } else {
3019 0 : dynFormat += "Rb(" + a + ')';
3020 : }
3021 0 : if (i != 2 * nlayer) {
3022 0 : dynFormat += "===";
3023 : }
3024 0 : }
3025 0 : print(files.TarcogIterationsFile, dynFormat);
3026 0 : print(files.TarcogIterationsFile, "\n");
3027 : // write Rb and Rf
3028 0 : print(files.TarcogIterationsFile, "{:16.8F} {:16.8F}", Rf(1), Rb(1));
3029 0 : for (i = 1; i <= nlayer; ++i) {
3030 0 : print(files.TarcogIterationsFile, " {:16.8F} {:16.8F}", Rf(i), Rb(i));
3031 : }
3032 0 : print(files.TarcogIterationsFile, "\n");
3033 :
3034 : // Write header for temperatures
3035 0 : for (i = 1; i <= 2 * nlayer; ++i) {
3036 0 : const std::string a = fmt::format("{:3}", i);
3037 0 : if (i == 1) {
3038 0 : dynFormat = "";
3039 : }
3040 0 : dynFormat += "theta(" + a + ')';
3041 0 : if (i != (2 * nlayer)) {
3042 0 : dynFormat += "==";
3043 : }
3044 0 : }
3045 0 : print(files.TarcogIterationsFile, dynFormat);
3046 0 : print(files.TarcogIterationsFile, "\n");
3047 :
3048 : // write temperatures
3049 0 : print(files.TarcogIterationsFile, "{:16.8F} \n", theta(1) - Constant::Kelvin);
3050 0 : for (i = 2; i <= 2 * nlayer; ++i) {
3051 0 : print(files.TarcogIterationsFile, " {:16.8F} \n", theta(i) - Constant::Kelvin);
3052 : }
3053 0 : print(files.TarcogIterationsFile, "\n");
3054 :
3055 : // close(TarcogIterationsFileNumber)
3056 :
3057 : // write results in csv file
3058 0 : if (index == 0) {
3059 0 : dynFormat = " ";
3060 0 : for (i = 1; i <= 2 * nlayer; ++i) {
3061 0 : const std::string a = fmt::format("{:3}", i);
3062 0 : if (i != 2 * nlayer) {
3063 0 : dynFormat += "theta(" + a + "),";
3064 : } else {
3065 0 : dynFormat += "theta(" + a + ')';
3066 : }
3067 0 : }
3068 0 : print(files.IterationCSVFile, dynFormat);
3069 0 : print(files.IterationCSVFile, "\n");
3070 : }
3071 0 : print(files.IterationCSVFile, "{:16.8F} \n", theta(1) - Constant::Kelvin);
3072 0 : for (i = 2; i <= 2 * nlayer; ++i) {
3073 0 : print(files.IterationCSVFile, " {:16.8F} \n", theta(i) - Constant::Kelvin);
3074 : }
3075 0 : print(files.IterationCSVFile, "\n");
3076 :
3077 : // close(IterationCSVFileNumber)
3078 0 : }
3079 :
3080 0 : void CalculateFuncResults(int const nlayer, Array2<Real64> const &a, const Array1D<Real64> &b, const Array1D<Real64> &x, Array1D<Real64> &FRes)
3081 : {
3082 : // Tuned Rewritten to traverse a in unit stride order
3083 0 : int const nlayer4(4 * nlayer);
3084 0 : for (int i = 1; i <= nlayer4; ++i) {
3085 0 : FRes(i) = -b(i);
3086 : }
3087 0 : for (int j = 1; j <= nlayer4; ++j) {
3088 0 : Real64 const x_j(x(j));
3089 0 : for (int i = 1; i <= nlayer4; ++i) {
3090 0 : FRes(i) += a(j, i) * x_j;
3091 : }
3092 : }
3093 0 : }
3094 :
3095 : } // namespace EnergyPlus::ThermalISO15099Calc
|
