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