Skip to content

Commit ee78d63

Browse files
mwettermwetter
authored
Merge pull request #2053 from ibpsa/issue2043_TwoWayFlowElement
Improve computing performance of multizone air exchange models.
2 parents fd9aa08 + 9e1bec0 commit ee78d63

25 files changed

+489
-120
lines changed

‎IBPSA/Airflow/Multizone/BaseClasses/Door.mo‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -43,6 +43,9 @@ protected
4343
constant Real conTP = IBPSA.Media.Air.dStp*Modelica.Media.IdealGases.Common.SingleGasesData.Air.R_s
4444
"Conversion factor for converting temperature difference to pressure difference";
4545

46+
final parameter Real sqrt_dp_turbulent(min=0) = sqrt(dp_turbulent)
47+
"Square root of pressure difference where laminar and turbulent flow relation coincide";
48+
4649
parameter Medium.ThermodynamicState sta_default=Medium.setState_pTX(
4750
T=Medium.T_default,
4851
p=Medium.p_default,
@@ -118,6 +121,13 @@ This is a partial model for the bi-directional air flow through a door.
118121
revisions="<html>
119122
<ul>
120123
<li>
124+
September 19, 2025, by Michael Wetter:<br/>
125+
Introduced protected parameter <code>sqrt_dp_turbulent</code>,
126+
which is needed to improve computing efficiency if flow exponent is <i>0.5</i>.<br/>
127+
This is for
128+
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/2043\">IBPSA, #2043</a>.
129+
</li>
130+
<li>
121131
October 6, 2020, by Michael Wetter:<br/>
122132
First implementation for
123133
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/1353\">#1353</a>.

‎IBPSA/Airflow/Multizone/BaseClasses/DoorDiscretized.mo‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -28,6 +28,9 @@ protected
2828
parameter Modelica.Units.SI.Density rho_default=Medium.density(sta_default)
2929
"Density, used to compute fluid volume";
3030

31+
final parameter Real sqrt_dp_turbulent(min=0) = sqrt(dp_turbulent)
32+
"Square root of pressure difference where laminar and turbulent flow relation coincide";
33+
3134
input Real hAg[nCom](each unit="m2/s2")=
3235
{Modelica.Constants.g_n*(hA - (i - 0.5)*dh) for i in 1:nCom}
3336
"Product g*h_i for each compartment";
@@ -118,6 +121,13 @@ using the model for a door that can be open or closed.
118121
revisions="<html>
119122
<ul>
120123
<li>
124+
September 19, 2025, by Michael Wetter:<br/>
125+
Introduced protected parameter <code>sqrt_dp_turbulent</code>,
126+
which is needed to improve computing efficiency if flow exponent is <i>0.5</i>.<br/>
127+
This is for
128+
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/2043\">IBPSA, #2043</a>.
129+
</li>
130+
<li>
121131
October 29, 2024, by Klaas De Jonge:<br/>
122132
Unprotected <code>dh</code> and changed prefixes of <code>dh</code>,<code>hAg</code> and <code>hBg</code> to <code>input</code>.<br/>
123133
This is for <a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/1935\">#1935</a>.

‎IBPSA/Airflow/Multizone/BaseClasses/Examples/PowerLaw05.mo‎

Lines changed: 73 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,73 @@
1+
within IBPSA.Airflow.Multizone.BaseClasses.Examples;
2+
model PowerLaw05
3+
"Test model for power law function with constant exponent of 0.5"
4+
extends Modelica.Icons.Example;
5+
parameter Real C = 2/10^m "Flow coefficient, C = V_flow/ dp^m";
6+
7+
constant Real m(min=0.5, max=1) = 0.5
8+
"Flow exponent, m=0.5 for turbulent, m=1 for laminar";
9+
parameter Modelica.Units.SI.PressureDifference dp_turbulent(min=0) = 5
10+
"Pressure difference where regularization starts";
11+
final parameter Real sqrt_dp_turbulent(min=0) = sqrt(dp_turbulent)
12+
"Square root of pressure difference where regularization starts";
13+
14+
Modelica.Units.SI.PressureDifference dp "Pressure difference";
15+
Modelica.Units.SI.VolumeFlowRate V_flow
16+
"Volume flow rate computed with model powerLaw";
17+
Modelica.Units.SI.VolumeFlowRate VFixed_flow
18+
"Volume flow rate computed with model powerLawFixed";
19+
20+
constant Real gamma(min=1) = 1.5
21+
"Normalized flow rate where dphi(0)/dpi intersects phi(1)";
22+
constant Real a = gamma
23+
"Polynomial coefficient for regularized implementation of flow resistance";
24+
constant Real b = 1/8*m^2 - 3*gamma - 3/2*m + 35.0/8
25+
"Polynomial coefficient for regularized implementation of flow resistance";
26+
constant Real c = -1/4*m^2 + 3*gamma + 5/2*m - 21.0/4
27+
"Polynomial coefficient for regularized implementation of flow resistance";
28+
constant Real d = 1/8*m^2 - gamma - m + 15.0/8
29+
"Polynomial coefficient for regularized implementation of flow resistance";
30+
31+
equation
32+
dp = 10*(-1+2*time);
33+
V_flow = IBPSA.Airflow.Multizone.BaseClasses.powerLaw(
34+
dp=dp,
35+
C=C,
36+
m=m,
37+
dp_turbulent=dp_turbulent);
38+
VFixed_flow = IBPSA.Airflow.Multizone.BaseClasses.powerLaw05(
39+
C=C,
40+
dp=dp,
41+
a=a,
42+
b=b,
43+
c=c,
44+
d=d,
45+
dp_turbulent=dp_turbulent,
46+
sqrt_dp_turbulent=sqrt_dp_turbulent);
47+
assert(abs(V_flow-VFixed_flow) < 1E-10, "Error: The two implementations of the power law model need to give identical results");
48+
annotation (
49+
experiment(Tolerance=1e-6, StopTime=1.0),
50+
__Dymola_Commands(file="modelica://IBPSA/Resources/Scripts/Dymola/Airflow/Multizone/BaseClasses/Examples/PowerLaw05.mos"
51+
"Simulate and plot"), Documentation(info="<html>
52+
<p>
53+
This examples demonstrates the
54+
<a href=\"modelica://IBPSA.Airflow.Multizone.BaseClasses.powerLaw\">
55+
IBPSA.Airflow.Multizone.BaseClasses.powerLaw</a>
56+
and
57+
<a href=\"modelica://IBPSA.Airflow.Multizone.BaseClasses.powerLaw05\">
58+
IBPSA.Airflow.Multizone.BaseClasses.powerLaw05</a>
59+
functions.
60+
They need to return the same function value.
61+
This is verified by an <code>assert</code> statement.
62+
</p>
63+
</html>", revisions="<html>
64+
<ul>
65+
<li>
66+
September 19. 2025, by Michael Wetter:<br/>
67+
First implementation.<br/>
68+
This is for
69+
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/2043\">IBPSA, #2043</a>.
70+
</li>
71+
</ul>
72+
</html>"));
73+
end PowerLaw05;

‎IBPSA/Airflow/Multizone/BaseClasses/Examples/package.order‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,5 @@
11
PowerLaw
2+
PowerLaw05
23
PowerLawFixedM
34
WindPressureLowRise
45
WindPressureProfile

‎IBPSA/Airflow/Multizone/BaseClasses/PartialOneWayFlowElement.mo‎

Lines changed: 23 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -31,13 +31,15 @@ protected
3131
Medium.dynamicViscosity(sta_default)
3232
"Dynamic viscosity at the medium default properties";
3333

34+
final parameter Real sqrt_dp_turbulent(min=0) = sqrt(dp_turbulent)
35+
"Square root of pressure difference where laminar and turbulent flow relation coincide";
36+
3437
Medium.ThermodynamicState sta "State of the medium in the component";
3538
Modelica.Units.SI.DynamicViscosity dynVis "Dynamic viscosity";
36-
Real mExc(quantity="Mass", final unit="kg")
39+
Modelica.Units.SI.Mass mExc(start=0, fixed=true, unbounded=true)
3740
"Air mass exchanged (for purpose of error control only)";
3841

3942
initial equation
40-
mExc=0;
4143
assert(homotopyInitialization, "In " + getInstanceName() +
4244
": The constant homotopyInitialization has been modified from its default value. This constant will be removed in future releases.",
4345
level = AssertionLevel.warning);
@@ -110,6 +112,25 @@ not in the equation section since this model sets both
110112
revisions="<html>
111113
<ul>
112114
<li>
115+
September 22, 2025, by Michael Wetter:<br/>
116+
Set <code>unbounded=true</code> for <code>mExc_flow</code>.
117+
This is to have the same implementation as in
118+
<a href=\"modelica://IBPSA.Airflow.Multizone.BaseClasses.TwoWayFlowElement\">
119+
IBPSA.Airflow.Multizone.BaseClasses.TwoWayFlowElement</a>
120+
where this change was done to avoid spikes in <code>port_a.m_flow</code> in
121+
<a href=\"modelica://IBPSA.Airflow.Multizone.Examples.OneOpenDoor\">
122+
IBPSA.Airflow.Multizone.Examples.OneOpenDoor</a>.<br/>
123+
This is for
124+
<a href=\"https://github.com/lbl-srg/modelica-buildings/issues/4360\">Buildings, #4360</a>.
125+
</li>
126+
<li>
127+
September 19, 2025, by Michael Wetter:<br/>
128+
Introduced protected parameter <code>sqrt_dp_turbulent</code>,
129+
which is needed to improve computing efficiency if flow exponent is <i>0.5</i>.<br/>
130+
This is for
131+
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/2043\">IBPSA, #2043</a>.
132+
</li>
133+
<li>
113134
February 2, 2022, by Michael Wetter:<br/>
114135
Revised implementation.<br/>
115136
This is for

‎IBPSA/Airflow/Multizone/BaseClasses/TwoWayFlowElement.mo‎

Lines changed: 28 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -41,9 +41,9 @@ protected
4141
Modelica.Units.SI.VolumeFlowRate VZer_flow(fixed=false)
4242
"Minimum net volume flow rate to prevent zero flow";
4343

44-
Modelica.Units.SI.Mass mExcAB(start=0, fixed=true)
44+
Modelica.Units.SI.Mass mExcAB(start=0, fixed=true, unbounded=true)
4545
"Air mass exchanged (for purpose of error control only)";
46-
Modelica.Units.SI.Mass mExcBA(start=0, fixed=true)
46+
Modelica.Units.SI.Mass mExcBA(start=0, fixed=true, unbounded=true)
4747
"Air mass exchanged (for purpose of error control only)";
4848

4949
Medium.MassFraction Xi_a1_inflow[Medium1.nXi]
@@ -55,8 +55,8 @@ equation
5555
// gives higher robustness. The reason may be that for bi-directional flow,
5656
// (VAB_flow - VBA_flow) may be close to zero.
5757
if forceErrorControlOnFlow then
58-
der(mExcAB) = mAB_flow;
59-
der(mExcBA) = mBA_flow;
58+
der(mExcAB) = port_a1.m_flow;
59+
der(mExcBA) = port_a2.m_flow;
6060
else
6161
der(mExcAB) = 0;
6262
der(mExcBA) = 0;
@@ -82,14 +82,14 @@ equation
8282

8383
VZer_flow = vZer*A;
8484

85-
mAB_flow = rho_a1_inflow*VAB_flow;
86-
mBA_flow = rho_a2_inflow*VBA_flow;
85+
mAB_flow = port_a1.m_flow;
86+
mBA_flow = port_a2.m_flow;
8787
// Average velocity (using the whole orifice area)
8888
vAB = VAB_flow/A;
8989
vBA = VBA_flow/A;
9090

91-
port_a1.m_flow = mAB_flow;
92-
port_a2.m_flow = mBA_flow;
91+
port_a1.m_flow = rho_a1_inflow*VAB_flow;
92+
port_a2.m_flow = rho_a2_inflow*VBA_flow;
9393

9494
// Energy balance (no storage, no heat loss/gain)
9595
port_a1.h_outflow = inStream(port_b1.h_outflow);
@@ -135,6 +135,26 @@ for doors that can be open or closed as a function of an input signal.
135135
revisions="<html>
136136
<ul>
137137
<li>
138+
September 22, 2025, by Michael Wetter:<br/>
139+
Set <code>unbounded=true</code> for <code>mExcAB_flow</code> and <code>mExcBA_flow</code>
140+
to avoid spikes in <code>port_a.m_flow</code> in
141+
<a href=\"modelica://IBPSA.Airflow.Multizone.Examples.OneOpenDoor\">
142+
IBPSA.Airflow.Multizone.Examples.OneOpenDoor</a>.<br/>
143+
This is for
144+
<a href=\"https://github.com/lbl-srg/modelica-buildings/issues/4360\">Buildings, #4360</a>.
145+
</li>
146+
<li>
147+
September 19, 2025, by Michael Wetter:<br/>
148+
Refactored implementation of underlying function that computes the flow rate
149+
to allow function to be inlined.
150+
This leads to a 20% faster simulation of
151+
<a href=\"modelica://IBPSA.Airflow.Multizone.Examples.OneOpenDoor\">
152+
IBPSA.Airflow.Multizone.Examples.OneOpenDoor</a> compared to the previous
153+
implementation.<br/>
154+
This is for
155+
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/2043\">IBPSA, #2043</a>.
156+
</li>
157+
<li>
138158
May 12, 2020, by Michael Wetter:<br/>
139159
Changed assignment of <code>m1_flow_small</code> and
140160
<code>m2_flow_small</code> to <code>final</code>.

‎IBPSA/Airflow/Multizone/BaseClasses/package.order‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -7,6 +7,7 @@ TwoWayFlowElement
77
TwoWayFlowElementBuoyancy
88
ZonalFlow
99
powerLaw
10+
powerLaw05
1011
powerLawFixedM
1112
windPressureLowRise
1213
windPressureProfile

‎IBPSA/Airflow/Multizone/BaseClasses/powerLaw.mo‎

Lines changed: 26 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -14,32 +14,31 @@ function powerLaw "Power law used in orifice equations"
1414
protected
1515
constant Real gamma(min=1) = 1.5
1616
"Normalized flow rate where dphi(0)/dpi intersects phi(1)";
17-
Real a
18-
"Polynomial coefficient for regularized implementation of flow resistance";
19-
Real b
20-
"Polynomial coefficient for regularized implementation of flow resistance";
21-
Real c
22-
"Polynomial coefficient for regularized implementation of flow resistance";
23-
Real d
24-
"Polynomial coefficient for regularized implementation of flow resistance";
25-
Real pi "Normalized pressure";
26-
Real pi2 "Square of normalized pressure";
17+
Real pi = dp/dp_turbulent "Normalized pressure";
18+
Real pi2 = pi*pi "Square of normalized pressure";
2719
algorithm
28-
if (dp >= dp_turbulent) then
29-
V_flow :=C *dp^m;
20+
V_flow := if (dp >= dp_turbulent) then
21+
C *dp^m
3022
elseif (dp <= -dp_turbulent) then
31-
V_flow :=-C*(-dp)^m;
23+
-C*(-dp)^m
3224
else
25+
C *dp_turbulent^m * pi *
26+
( gamma + pi2 *
27+
( (1/8*m^2 - 3*gamma - 3/2*m + 35.0/8) + pi2 *
28+
( (-1/4*m^2 + 3*gamma + 5/2*m - 21.0/4) + pi2 *
29+
(1/8*m^2 - gamma - m + 15.0/8))));
30+
/*
3331
a := gamma;
3432
b := 1/8*m^2 - 3*gamma - 3/2*m + 35.0/8;
3533
c := -1/4*m^2 + 3*gamma + 5/2*m - 21.0/4;
3634
d := 1/8*m^2 - gamma - m + 15.0/8;
37-
pi := dp/dp_turbulent;
38-
pi2 := pi*pi;
3935
V_flow :=C *dp_turbulent^m * pi * ( a + pi2 * ( b + pi2 * ( c + pi2 * d)));
40-
end if;
36+
*/
4137

42-
annotation (smoothOrder=2,
38+
39+
annotation (
40+
smoothOrder=2,
41+
Inline=true,
4342
Documentation(info="<html>
4443
<p>
4544
This model describes the mass flow rate and pressure difference relation
@@ -82,6 +81,16 @@ of a model.
8281
revisions="<html>
8382
<ul>
8483
<li>
84+
September 19, 2025, by Michael Wetter:<br/>
85+
Refactored implementation to allow function to be inlined.
86+
This leads to a 20% faster simulation of
87+
<a href=\"modelica://IBPSA.Airflow.Multizone.Examples.OneOpenDoor\">
88+
IBPSA.Airflow.Multizone.Examples.OneOpenDoor</a> compared to the previous
89+
implementation.<br/>
90+
This is for
91+
<a href=\"https://github.com/ibpsa/modelica-ibpsa/issues/2043\">IBPSA, #2043</a>.
92+
</li>
93+
<li>
8594
February 8, 2022, by Michael Wetter:<br/>
8695
Changed to use <code>C</code> for volume flow coefficient (<i>C = V_flow/dp^m</i>),
8796
and <code>k</code> for mass flow coefficient (<i>k = m_flow/dp^m</i>).

0 commit comments

Comments
 (0)