model ControlledFlowMachine "Partial model for fan or pump with ideally controlled mass flow rate or head as input signal" extends Buildings.Fluid.Movers.BaseClasses.PowerInterface( final use_powerCharacteristic = false, final rho_default = Medium.density(sta_default)); extends Buildings.Fluid.Movers.BaseClasses.PartialFlowMachine( preSou(final control_m_flow=control_m_flow)); import cha = Buildings.Fluid.Movers.BaseClasses.Characteristics; // parameter Modelica.SIunits.MassFlowRate m_flow_nominal // "Nominal mass flow rate, used as flow rate if control_m_flow"; // parameter Modelica.SIunits.MassFlowRate m_flow_max = m_flow_nominal // "Maximum mass flow rate (at zero head)"; // what to control constant Boolean control_m_flow "= false to control head instead of m_flow"; Real r_V(start=1) "Ratio V_flow/V_flow_max = V_flow/V_flow(dp=0, N=N_nominal)"; protected final parameter Medium.AbsolutePressure p_a_default(displayUnit="Pa") = Medium.p_default "Nominal inlet pressure for predefined fan or pump characteristics"; parameter Medium.ThermodynamicState sta_default = Medium.setState_pTX( T=Medium.T_default, p=Medium.p_default, X=Medium.X_default[1:Medium.nXi]) "Default medium state"; Modelica.Blocks.Sources.RealExpression PToMedium_flow(y=Q_flow + WFlo) if addPowerToMedium "Heat and work input into medium"; initial equation V_flow_max=m_flow_nominal/rho_default; equation r_V = VMachine_flow/V_flow_max; etaHyd = cha.efficiency(data=hydraulicEfficiency, r_V=r_V, d=hydDer); etaMot = cha.efficiency(data=motorEfficiency, r_V=r_V, d=motDer); dpMachine = -dp; VMachine_flow = -port_b.m_flow/rho_in; // To compute the electrical power, we set a lower bound for eta to avoid // a division by zero. P = WFlo / Buildings.Utilities.Math.Functions.smoothMax(x1=eta, x2=1E-5, deltaX=1E-6); connect(PToMedium_flow.y, prePow.Q_flow); end ControlledFlowMachine;

partial model FlowMachineInterface "Partial model with performance curves for fans or pumps" extends Buildings.Fluid.Movers.BaseClasses.PowerInterface( VMachine_flow(nominal=V_flow_nominal, start=V_flow_nominal), V_flow_max(nominal=V_flow_nominal, start=V_flow_nominal)); import Modelica.Constants; import cha = Buildings.Fluid.Movers.BaseClasses.Characteristics; parameter Modelica.SIunits.Conversions.NonSIunits.AngularVelocity_rpm N_nominal = 1500 "Nominal rotational speed for flow characteristic"; final parameter Modelica.SIunits.VolumeFlowRate V_flow_nominal= pressure.V_flow[size(pressure.V_flow,1)] "Nominal volume flow rate, used for homotopy"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.flowParameters pressure "Volume flow rate vs. total pressure rise"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.powerParameters power "Volume flow rate vs. electrical power consumption"; parameter Boolean homotopyInitialization = true "= true, use homotopy method"; // Classes used to implement the filtered speed parameter Boolean filteredSpeed=true "= true, if speed is filtered with a 2nd order CriticalDamping filter"; parameter Modelica.SIunits.Time riseTime=30 "Rise time of the filter (time to reach 99.6 % of the speed)"; parameter Modelica.Blocks.Types.Init init=Modelica.Blocks.Types.Init.InitialOutput "Type of initialization (no init/steady state/initial state/initial output)"; parameter Real N_start=0 "Initial value of speed"; // Speed Modelica.Blocks.Interfaces.RealOutput N_actual(min=0, max=N_nominal, final quantity="AngularVelocity", final unit="1/min", nominal=N_nominal); // "Shaft rotational speed in rpm"; Real r_N(min=0, start=N_start/N_nominal, unit="1") "Ratio N_actual/N_nominal"; Real r_V(start=1, unit="1") "Ratio V_flow/V_flow_max"; protected Modelica.Blocks.Interfaces.RealOutput N_filtered(min=0, start=N_start, max=N_nominal) if filteredSpeed "Filtered speed in the range 0..N_nominal"; Modelica.Blocks.Continuous.Filter filter( order=2, f_cut=5/(2*Modelica.Constants.pi*riseTime), final init=init, final y_start=N_start, x(each stateSelect=StateSelect.always), u_nominal=N_nominal, u(final quantity="AngularVelocity", final unit="1/min", nominal=N_nominal), y(final quantity="AngularVelocity", final unit="1/min", nominal=N_nominal), final analogFilter=Modelica.Blocks.Types.AnalogFilter.CriticalDamping, final filterType=Modelica.Blocks.Types.FilterType.LowPass) if filteredSpeed "Second order filter to approximate valve opening time, and to improve numerics"; parameter Modelica.SIunits.VolumeFlowRate VDelta_flow(fixed=false, start=delta*V_flow_nominal) "Small volume flow rate"; parameter Modelica.SIunits.Pressure dpDelta(fixed=false, start=100) "Small pressure"; parameter Real delta = 0.05 "Small value used to transition to other fan curve"; parameter Real cBar[2](each fixed=false) "Coefficients for linear approximation of pressure vs. flow rate"; parameter Modelica.SIunits.Pressure dpMax(min=0, fixed=false); parameter Real kRes(min=0, unit="kg/(s.m4)", fixed=false) "Coefficient for internal pressure drop of fan or pump"; parameter Integer curve(min=1, max=3, fixed=false) "Flag, used to pick the right representatio of the fan or pump pressure curve"; parameter Integer nOri = size(pressure.V_flow,1) "Number of data points for pressure curve"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.flowParametersInternal pCur1( final n = nOri, V_flow(each fixed=false), dp(each fixed=false)) "Volume flow rate vs. total pressure rise with correction for pump resistance added"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.flowParametersInternal pCur2( final n = nOri + 1, V_flow(each fixed=false), dp(each fixed=false)) "Volume flow rate vs. total pressure rise with correction for pump resistance added"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.flowParametersInternal pCur3( final n = nOri + 2, V_flow(each fixed=false), dp(each fixed=false)) "Volume flow rate vs. total pressure rise with correction for pump resistance added"; parameter Real preDer1[nOri](each fixed=false) "Derivatives of flow rate vs. pressure at the support points"; parameter Real preDer2[nOri+1](each fixed=false) "Derivatives of flow rate vs. pressure at the support points"; parameter Real preDer3[nOri+2](each fixed=false) "Derivatives of flow rate vs. pressure at the support points"; parameter Real powDer[size(power.V_flow,1)]= if use_powerCharacteristic then Buildings.Utilities.Math.Functions.splineDerivatives( x=power.V_flow, y=power.P, ensureMonotonicity=Buildings.Utilities.Math.Functions.isMonotonic(x=power.P, strict=false)) else zeros(size(power.V_flow,1)) "Coefficients for polynomial of pressure vs. flow rate"; parameter Boolean haveMinimumDecrease(fixed=false) "Flag used for reporting"; parameter Boolean haveDPMax(fixed=false) "Flag, true if user specified data that contain dpMax"; parameter Boolean haveVMax(fixed=false) "Flag, true if user specified data that contain V_flow_max"; // Variables Modelica.SIunits.Density rho "Medium density"; function getPerformanceDataAsString input Buildings.Fluid.Movers.BaseClasses.Characteristics.flowParameters pressure "Performance data"; input Real derivative[:](unit="kg/(s.m4)") "Derivative"; input Integer minimumLength = 6 "Minimum width of result"; input Integer significantDigits = 6 "Number of significant digits"; output String str "String representation"; algorithm str :=""; for i in 1:size(derivative, 1) loop str :=str + " V_flow[" + String(i) + "]=" + String( pressure.V_flow[i], minimumLength=minimumLength, significantDigits=significantDigits) + "\t" + "dp[" + String(i) + "]=" + String( pressure.dp[i], minimumLength=minimumLength, significantDigits=significantDigits) + "\tResulting derivative dp/dV_flow = " + String( derivative[i], minimumLength=minimumLength, significantDigits=significantDigits) + "\n"; end for; end getPerformanceDataAsString; function getArrayAsString input Real array[:] "Array to be printed"; input String varName "Variable name"; input Integer minimumLength = 6 "Minimum width of result"; input Integer significantDigits = 6 "Number of significant digits"; output String str "String representation"; algorithm str :=""; for i in 1:size(array, 1) loop str :=str + " " + varName + "[" + String(i) + "]=" + String( array[i], minimumLength=minimumLength, significantDigits=significantDigits) + "\n"; end for; end getArrayAsString; initial algorithm // Check validity of data assert(size(pressure.V_flow, 1) > 1, "Must have at least two data points for pressure.V_flow."); assert(Buildings.Utilities.Math.Functions.isMonotonic(x=pressure.V_flow, strict=true) and pressure.V_flow[1] > -Modelica.Constants.eps, "The volume flow rate for the fan pressure rise must be a strictly decreasing sequence with the first element being non-zero. The following performance data have been entered: " + getArrayAsString(pressure.V_flow, "pressure.V_flow")); // Check if V_flow_max or dpMax are provided by user haveVMax :=(abs(pressure.dp[nOri]) < Modelica.Constants.eps); haveDPMax :=(abs(pressure.V_flow[1]) < Modelica.Constants.eps); // Assign V_flow_max and dpMax if haveVMax then V_flow_max :=pressure.V_flow[nOri]; else assert((pressure.V_flow[nOri]-pressure.V_flow[nOri-1])/((pressure.dp[nOri]-pressure.dp[nOri-1]))<0, "The last two pressure points for the fan or pump performance curve must be decreasing. You need to set more reasonable parameters. Received " + getArrayAsString(pressure.dp, "dp")); V_flow_max :=pressure.V_flow[nOri] - (pressure.V_flow[nOri] - pressure.V_flow[ nOri - 1])/((pressure.dp[nOri] - pressure.dp[nOri - 1]))*pressure.dp[nOri]; end if; if haveDPMax then dpMax :=pressure.dp[1]; else dpMax :=pressure.dp[1] - ((pressure.dp[2] - pressure.dp[1])/(pressure.V_flow[ 2] - pressure.V_flow[1]))*pressure.V_flow[1]; end if; // Check if minimum decrease condition is satisfied haveMinimumDecrease :=true; kRes :=dpMax/V_flow_max*delta^2/10; for i in 1:nOri-1 loop if ((pressure.dp[i+1]-pressure.dp[i])/(pressure.V_flow[i+1]-pressure.V_flow[i]) >= -kRes) then haveMinimumDecrease :=false; end if; end for; // Write warning if the volumetric flow rate versus pressure curve does not satisfy // the minimum decrease condition if (not haveMinimumDecrease) then Modelica.Utilities.Streams.print(" Warning: ======== It is recommended that the volume flow rate versus pressure relation of the fan or pump satisfies the minimum decrease condition (pressure.dp[i+1]-pressure.dp[i]) d[i] = ----------------------------------------- < " + String(-kRes) + " (pressure.V_flow[i+1]-pressure.V_flow[i]) is " + getArrayAsString({(pressure.dp[i+1]-pressure.dp[i])/(pressure.V_flow[i+1]-pressure.V_flow[i]) for i in 1:nOri-1}, "d") + " Otherwise, a solution to the equations may not exist if the fan or pump speed is reduced. In this situation, the solver will fail due to non-convergence and the simulation stops."); end if; // Correction for flow resistance of pump or fan // Case 1: if (haveVMax and haveDPMax) or (nOri == 2) then // ----- Curve 1 curve :=1; // V_flow_max and dpMax are provided by the user, or we only have two data points for i in 1:nOri loop pCur1.dp[i] :=pressure.dp[i] + pressure.V_flow[i] * kRes; pCur1.V_flow[i] := pressure.V_flow[i]; end for; pCur2.V_flow := zeros(nOri + 1); pCur2.dp := zeros(nOri + 1); pCur3.V_flow := zeros(nOri + 2); pCur3.dp := zeros(nOri + 2); preDer1:=Buildings.Utilities.Math.Functions.splineDerivatives(x=pCur1.V_flow, y=pCur1.dp); preDer2:=zeros(nOri+1); preDer3:=zeros(nOri+2); // Equation to compute dpDelta dpDelta :=cha.pressure( data=pCur1, V_flow=0, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer1, cBar=zeros(2), kRes= kRes); // Equation to compute VDelta_flow. By the affinity laws, the volume flow rate is proportional to the speed. VDelta_flow :=V_flow_max*delta; // Linear equations to determine cBar // Conditions for r_N=delta, V_flow = VDelta_flow // Conditions for r_N=delta, V_flow = 0 cBar[1] :=cha.pressure( data=pCur1, V_flow=0, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer1, cBar=zeros(2), kRes= kRes) * (1-delta)/delta^2; cBar[2] :=((cha.pressure( data=pCur1, V_flow=VDelta_flow, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer1, cBar=zeros(2), kRes= kRes) - delta*dpDelta)/delta^2 - cBar[1])/VDelta_flow; elseif haveVMax or haveDPMax then // ----- Curve 2 curve :=2; // V_flow_max or dpMax is provided by the user, but not both if haveVMax then pCur2.V_flow[1] := 0; pCur2.dp[1] := dpMax; for i in 1:nOri loop pCur2.dp[i+1] :=pressure.dp[i] + pressure.V_flow[i] * kRes; pCur2.V_flow[i+1] := pressure.dp[i]; end for; else for i in 1:nOri loop pCur2.dp[i] :=pressure.dp[i] + pressure.V_flow[i] * kRes; pCur2.V_flow[i] := pressure.V_flow[i]; end for; pCur2.V_flow[nOri+1] := V_flow_max; pCur2.dp[nOri+1] := 0; end if; pCur1.V_flow := zeros(nOri); pCur1.dp := zeros(nOri); pCur3.V_flow := zeros(nOri + 2); pCur3.dp := zeros(nOri + 2); preDer1:=zeros(nOri); preDer2:=Buildings.Utilities.Math.Functions.splineDerivatives(x=pCur2.V_flow, y=pCur2.dp); preDer3:=zeros(nOri+2); // Equation to compute dpDelta dpDelta :=cha.pressure( data=pCur2, V_flow=0, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer2, cBar=zeros(2), kRes= kRes); // Equation to compute VDelta_flow. By the affinity laws, the volume flow rate is proportional to the speed. VDelta_flow :=V_flow_max*delta; // Linear equations to determine cBar // Conditions for r_N=delta, V_flow = VDelta_flow // Conditions for r_N=delta, V_flow = 0 cBar[1] :=cha.pressure( data=pCur2, V_flow=0, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer2, cBar=zeros(2), kRes= kRes) * (1-delta)/delta^2; cBar[2] :=((cha.pressure( data=pCur2, V_flow=VDelta_flow, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer2, cBar=zeros(2), kRes= kRes) - delta*dpDelta)/delta^2 - cBar[1])/VDelta_flow; else // ----- Curve 3 curve :=3; // Neither V_flow_max nor dpMax are provided by the user pCur3.V_flow[1] := 0; pCur3.dp[1] := dpMax; for i in 1:nOri loop pCur3.dp[i+1] :=pressure.dp[i] + pressure.V_flow[i] * kRes; pCur3.V_flow[i+1] := pressure.V_flow[i]; end for; pCur3.V_flow[nOri+2] := V_flow_max; pCur3.dp[nOri+2] := 0; pCur1.V_flow := zeros(nOri); pCur1.dp := zeros(nOri); pCur2.V_flow := zeros(nOri + 1); pCur2.dp := zeros(nOri + 1); preDer1:=zeros(nOri); preDer2:=zeros(nOri+1); preDer3:=Buildings.Utilities.Math.Functions.splineDerivatives(x=pCur3.V_flow, y=pCur3.dp); // Equation to compute dpDelta dpDelta :=cha.pressure( data=pCur3, V_flow=0, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer3, cBar=zeros(2), kRes= kRes); // Equation to compute VDelta_flow. By the affinity laws, the volume flow rate is proportional to the speed. VDelta_flow :=V_flow_max*delta; // Linear equations to determine cBar // Conditions for r_N=delta, V_flow = VDelta_flow // Conditions for r_N=delta, V_flow = 0 cBar[1] :=cha.pressure( data=pCur3, V_flow=0, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer3, cBar=zeros(2), kRes= kRes) * (1-delta)/delta^2; cBar[2] :=((cha.pressure( data=pCur3, V_flow=VDelta_flow, r_N=delta, VDelta_flow=0, dpDelta=0, V_flow_max=Modelica.Constants.eps, dpMax=0, delta=0, d=preDer3, cBar=zeros(2), kRes= kRes) - delta*dpDelta)/delta^2 - cBar[1])/VDelta_flow; end if; equation // Hydraulic equations r_N = N_actual/N_nominal; r_V = VMachine_flow/V_flow_max; // For the homotopy method, we approximate dpMachine by an equation // that is linear in VMachine_flow, and that goes linearly to 0 as r_N goes to 0. // The three branches below are identical, except that we pass either // pCur1, pCur2 or pCur3, and preDer1, preDer2 or preDer3 if (curve == 1) then if homotopyInitialization then dpMachine = homotopy(actual=cha.pressure(data=pCur1, V_flow=VMachine_flow, r_N=r_N, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer1, cBar=cBar, kRes=kRes), simplified=r_N* (cha.pressure(data=pCur1, V_flow=V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer1, cBar=cBar, kRes=kRes) +(VMachine_flow-V_flow_nominal)* (cha.pressure(data=pCur1, V_flow=(1+delta)*V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer1, cBar=cBar, kRes=kRes) -cha.pressure(data=pCur1, V_flow=(1-delta)*V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer1, cBar=cBar, kRes=kRes)) /(2*delta*V_flow_nominal))); else dpMachine = cha.pressure(data=pCur1, V_flow=VMachine_flow, r_N=r_N, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer1, cBar=cBar, kRes=kRes); end if; // end of computation for this branch elseif (curve == 2) then if homotopyInitialization then dpMachine = homotopy(actual=cha.pressure(data=pCur2, V_flow=VMachine_flow, r_N=r_N, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer2, cBar=cBar, kRes=kRes), simplified=r_N* (cha.pressure(data=pCur2, V_flow=V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer2, cBar=cBar, kRes=kRes) +(VMachine_flow-V_flow_nominal)* (cha.pressure(data=pCur2, V_flow=(1+delta)*V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer2, cBar=cBar, kRes=kRes) -cha.pressure(data=pCur2, V_flow=(1-delta)*V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer2, cBar=cBar, kRes=kRes)) /(2*delta*V_flow_nominal))); else dpMachine = cha.pressure(data=pCur2, V_flow=VMachine_flow, r_N=r_N, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer2, cBar=cBar, kRes=kRes); end if; // end of computation for this branch else if homotopyInitialization then dpMachine = homotopy(actual=cha.pressure(data=pCur3, V_flow=VMachine_flow, r_N=r_N, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer3, cBar=cBar, kRes=kRes), simplified=r_N* (cha.pressure(data=pCur3, V_flow=V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer3, cBar=cBar, kRes=kRes) +(VMachine_flow-V_flow_nominal)* (cha.pressure(data=pCur3, V_flow=(1+delta)*V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer3, cBar=cBar, kRes=kRes) -cha.pressure(data=pCur3, V_flow=(1-delta)*V_flow_nominal, r_N=1, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer3, cBar=cBar, kRes=kRes)) /(2*delta*V_flow_nominal))); else dpMachine = cha.pressure(data=pCur3, V_flow=VMachine_flow, r_N=r_N, VDelta_flow=VDelta_flow, dpDelta=dpDelta, V_flow_max=V_flow_max, dpMax=dpMax, delta=delta, d=preDer3, cBar=cBar, kRes=kRes); end if; // end of computation for this branch end if; // Power consumption if use_powerCharacteristic then // For the homotopy, we want P/V_flow to be bounded as V_flow -> 0 to avoid a very high medium // temperature near zero flow. if homotopyInitialization then P = homotopy(actual=cha.power(data=power, V_flow=VMachine_flow, r_N=r_N, d=powDer), simplified=VMachine_flow/V_flow_nominal* cha.power(data=power, V_flow=V_flow_nominal, r_N=1, d=powDer)); else P = (rho/rho_default)*cha.power(data=power, V_flow=VMachine_flow, r_N=r_N, d=powDer); end if; // To compute the efficiency, we set a lower bound on the electricity consumption. // This is needed because WFlo can be close to zero when P is zero, thereby // causing a division by zero. // Earlier versions of the model computed WFlo = eta * P, but this caused // a division by zero. eta = WFlo / Buildings.Utilities.Math.Functions.smoothMax(x1=P, x2=1E-5, deltaX=1E-6); // In this configuration, we only now the total power consumption. // Because nothing is known about etaMot versus etaHyd, we set etaHyd=1. This will // cause etaMot=eta, because eta=etaHyd*etaMot. // Earlier versions used etaMot=sqrt(eta), but as eta->0, this function has // and infinite derivative. etaHyd = 1; else if homotopyInitialization then etaHyd = homotopy(actual=cha.efficiency(data=hydraulicEfficiency, r_V=r_V, d=hydDer), simplified=cha.efficiency(data=hydraulicEfficiency, r_V=1, d=hydDer)); etaMot = homotopy(actual=cha.efficiency(data=motorEfficiency, r_V=r_V, d=motDer), simplified=cha.efficiency(data=motorEfficiency, r_V=1, d=motDer)); else etaHyd = cha.efficiency(data=hydraulicEfficiency, r_V=r_V, d=hydDer); etaMot = cha.efficiency(data=motorEfficiency, r_V=r_V, d=motDer); end if; // To compute the electrical power, we set a lower bound for eta to avoid // a division by zero. P = WFlo / Buildings.Utilities.Math.Functions.smoothMax(x1=eta, x2=1E-5, deltaX=1E-6); end if; end FlowMachineInterface;

model IdealSource "Base class for pressure and mass flow source with optional power input" extends Modelica.Fluid.Interfaces.PartialTwoPortTransport(show_T=false); // what to control parameter Boolean control_m_flow "= false to control dp instead of m_flow"; Modelica.Blocks.Interfaces.RealInput m_flow_in if control_m_flow "Prescribed mass flow rate"; Modelica.Blocks.Interfaces.RealInput dp_in if not control_m_flow "Prescribed outlet pressure"; protected Modelica.Blocks.Interfaces.RealInput m_flow_internal "Needed to connect to conditional connector"; Modelica.Blocks.Interfaces.RealInput dp_internal "Needed to connect to conditional connector"; equation // Ideal control if control_m_flow then m_flow = m_flow_internal; dp_internal = 0; else dp = dp_internal; m_flow_internal = 0; end if; connect(dp_internal, dp_in); connect(m_flow_internal, m_flow_in); // Energy balance (no storage) port_a.h_outflow = inStream(port_b.h_outflow); port_b.h_outflow = inStream(port_a.h_outflow); end IdealSource;

partial model PartialFlowMachine "Partial model to interface fan or pump models with the medium" extends Buildings.Fluid.Interfaces.LumpedVolumeDeclarations; import Modelica.Constants; extends Buildings.Fluid.Interfaces.PartialTwoPortInterface(show_T=false, port_a( h_outflow(start=h_outflow_start), final m_flow(min = if allowFlowReversal then -Constants.inf else 0)), port_b( h_outflow(start=h_outflow_start), p(start=p_start), final m_flow(max = if allowFlowReversal then +Constants.inf else 0)), final showDesignFlowDirection=false); Delays.DelayFirstOrder vol( redeclare package Medium = Medium, tau=tau, energyDynamics=if dynamicBalance then energyDynamics else Modelica.Fluid.Types.Dynamics.SteadyState, massDynamics=if dynamicBalance then massDynamics else Modelica.Fluid.Types.Dynamics.SteadyState, T_start=T_start, X_start=X_start, C_start=C_start, m_flow_nominal=m_flow_nominal, p_start=p_start, prescribedHeatFlowRate=true, allowFlowReversal=allowFlowReversal, nPorts=2) "Fluid volume for dynamic model"; parameter Boolean dynamicBalance = true "Set to true to use a dynamic balance, which often leads to smaller systems of equations"; parameter Boolean addPowerToMedium=true "Set to false to avoid any power (=heat and flow work) being added to medium (may give simpler equations)"; parameter Modelica.SIunits.Time tau=1 "Time constant of fluid volume for nominal flow, used if dynamicBalance=true"; // Models Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a heatPort "Heat dissipation to environment"; protected Modelica.SIunits.Density rho_in "Density of inflowing fluid"; Buildings.Fluid.Movers.BaseClasses.IdealSource preSou( redeclare package Medium = Medium, allowFlowReversal=allowFlowReversal) "Pressure source"; Buildings.HeatTransfer.Sources.PrescribedHeatFlow prePow if addPowerToMedium "Prescribed power (=heat and flow work) flow for dynamic model"; parameter Medium.ThermodynamicState sta_start=Medium.setState_pTX( T=T_start, p=p_start, X=X_start) "Medium state at start values"; parameter Modelica.SIunits.SpecificEnthalpy h_outflow_start = Medium.specificEnthalpy(sta_start) "Start value for outflowing enthalpy"; equation // For computing the density, we assume that the fan operates in the design flow direction. rho_in = Medium.density( Medium.setState_phX(port_a.p, inStream(port_a.h_outflow), inStream(port_a.Xi_outflow))); connect(prePow.port, vol.heatPort); connect(vol.heatPort, heatPort); connect(port_a, vol.ports[1]); connect(vol.ports[2], preSou.port_a); connect(preSou.port_b, port_b); end PartialFlowMachine;

partial model PowerInterface "Partial model to compute power draw and heat dissipation of fans and pumps" import Modelica.Constants; parameter Boolean use_powerCharacteristic = false "Use powerCharacteristic (vs. efficiencyCharacteristic)"; parameter Boolean motorCooledByFluid = true "If true, then motor heat is added to fluid stream"; parameter Boolean homotopyInitialization = true "= true, use homotopy method"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.efficiencyParameters motorEfficiency(r_V={1}, eta={0.7}) "Normalized volume flow rate vs. efficiency"; parameter Buildings.Fluid.Movers.BaseClasses.Characteristics.efficiencyParameters hydraulicEfficiency(r_V={1}, eta={0.7}) "Normalized volume flow rate vs. efficiency"; parameter Modelica.SIunits.Density rho_default "Fluid density at medium default state"; Modelica.Blocks.Interfaces.RealOutput P(quantity="Modelica.SIunits.Power", unit="W") "Electrical power consumed"; Modelica.SIunits.Power WHyd "Hydraulic power input (converted to flow work and heat)"; Modelica.SIunits.Power WFlo "Flow work"; Modelica.SIunits.HeatFlowRate Q_flow "Heat input from fan or pump to medium"; Real eta(min=0, max=1) "Global efficiency"; Real etaHyd(min=0, max=1) "Hydraulic efficiency"; Real etaMot(min=0, max=1) "Motor efficiency"; Modelica.SIunits.Pressure dpMachine(displayUnit="Pa") "Pressure increase"; Modelica.SIunits.VolumeFlowRate VMachine_flow "Volume flow rate"; protected parameter Modelica.SIunits.VolumeFlowRate V_flow_max(fixed=false) "Maximum volume flow rate, used for smoothing"; //Modelica.SIunits.HeatFlowRate QThe_flow "Heat input into the medium"; parameter Modelica.SIunits.VolumeFlowRate delta_V_flow = 1E-3*V_flow_max "Factor used for setting heat input into medium to zero at very small flows"; final parameter Real motDer[size(motorEfficiency.r_V, 1)](each fixed=false) "Coefficients for polynomial of pressure vs. flow rate"; final parameter Real hydDer[size(hydraulicEfficiency.r_V,1)](each fixed=false) "Coefficients for polynomial of pressure vs. flow rate"; Modelica.SIunits.HeatFlowRate QThe_flow "Heat input from fan or pump to medium"; initial algorithm // Compute derivatives for cubic spline motDer := if use_powerCharacteristic then zeros(size(motorEfficiency.r_V, 1)) elseif ( size(motorEfficiency.r_V, 1) == 1) then {0} else Buildings.Utilities.Math.Functions.splineDerivatives( x=motorEfficiency.r_V, y=motorEfficiency.eta, ensureMonotonicity=Buildings.Utilities.Math.Functions.isMonotonic(x=motorEfficiency.eta, strict=false)); hydDer := if use_powerCharacteristic then zeros(size(hydraulicEfficiency.r_V, 1)) elseif ( size(hydraulicEfficiency.r_V, 1) == 1) then {0} else Buildings.Utilities.Math.Functions.splineDerivatives( x=hydraulicEfficiency.r_V, y=hydraulicEfficiency.eta); equation eta = etaHyd * etaMot; // WFlo = eta * P; // Flow work WFlo = dpMachine*VMachine_flow; // Hydraulic power (transmitted by shaft), etaHyd = WFlo/WHyd etaHyd * WHyd = WFlo; // Heat input into medium QThe_flow + WFlo = if motorCooledByFluid then P else WHyd; // At m_flow = 0, the solver may still obtain positive values for QThe_flow. // The next statement sets the heat input into the medium to zero for very small flow rates. if homotopyInitialization then Q_flow = homotopy(actual=Buildings.Utilities.Math.Functions.spliceFunction(pos=QThe_flow, neg=0, x=noEvent(abs(VMachine_flow))-2*delta_V_flow, deltax=delta_V_flow), simplified=0); else Q_flow = Buildings.Utilities.Math.Functions.spliceFunction(pos=QThe_flow, neg=0, x=noEvent(abs(VMachine_flow))-2*delta_V_flow, deltax=delta_V_flow); end if; end PowerInterface;

partial model PrescribedFlowMachine "Partial model for fan or pump with speed (y or Nrpm) as input signal" extends Buildings.Fluid.Movers.BaseClasses.FlowMachineInterface( V_flow_max(start=V_flow_nominal), final rho_default = Medium.density_pTX( p=Medium.p_default, T=Medium.T_default, X=Medium.X_default)); extends Buildings.Fluid.Movers.BaseClasses.PartialFlowMachine( final m_flow_nominal = max(pressure.V_flow)*rho_default, preSou(final control_m_flow=false)); // Models protected Modelica.Blocks.Sources.RealExpression dpMac(y=-dpMachine) "Pressure drop of the pump or fan"; Modelica.Blocks.Sources.RealExpression PToMedium_flow(y=Q_flow + WFlo) if addPowerToMedium "Heat and work input into medium"; equation VMachine_flow = -port_b.m_flow/rho; rho = rho_in; connect(preSou.dp_in, dpMac.y); connect(PToMedium_flow.y, prePow.Q_flow); end PrescribedFlowMachine;