LBL logo

Buildings.Fluid.HeatExchangers.Radiators

Package with radiators models for hydronic space heating systems

Information

This package contains models for radiators that are typically found in hydronic heating systems.

Extends from Modelica.Icons.VariantsPackage (Icon for package containing variants).

Package Content

NameDescription
Buildings.Fluid.HeatExchangers.Radiators.RadiatorEN442_2 RadiatorEN442_2 Dynamic radiator for space heating
Buildings.Fluid.HeatExchangers.Radiators.Examples Examples Collection of models that illustrate model use and test models

Buildings.Fluid.HeatExchangers.Radiators.RadiatorEN442_2 Buildings.Fluid.HeatExchangers.Radiators.RadiatorEN442_2

Dynamic radiator for space heating

Buildings.Fluid.HeatExchangers.Radiators.RadiatorEN442_2

Information

This is a model of a radiator that can be used as a dynamic or steady-state model. The required parameters are data that are typically available from manufacturers that follow the European Norm EN 442-2.

However, to allow for varying mass flow rates, the transferred heat is computed using a discretization along the water flow path, and heat is exchanged between each compartment and a uniform room air and radiation temperature. This discretization is different from the computation in EN 442-2, which may yield water outlet temperatures that are below the room temperature at low mass flow rates. Furthermore, rather than using only one room temperature, this model uses a room air and room radiation temperature.

The transferred heat is modeled as follows: Let N denote the number of elements used to discretize the radiator model. For each element i ∈ {1, … , N}, the convective and radiative heat transfer Qic and Qir from the radiator to the room is

Qic = sign(Ti-Ta) (1-fr) UA ⁄ N |Ti-Ta|n

Qir = sign(Ti-Tr) fr UA ⁄ N |Ti-Tr|n

where Ti is the water temperature of the element, Ta is the temperature of the room air, Tr is the radiative temperature, 0 < fr < 1 is the fraction of radiant to total heat transfer, UA is the UA-value of the radiator, and n is an exponent for the heat transfer. The model computes the UA-value by numerically solving the above equations for given nominal heating power, nominal temperatures, fraction radiant to total heat transfer and exponent for heat transfer.

The parameter energyDynamics (in the Assumptions tab), determines whether the model computes the dynamic or the steady-state response. For the transient response, heat storage is computed using a finite volume approach for the water and the metal mass, which are both assumed to be at the same temperature.

The default parameters for the heat capacities are valid for a flat plate radiator without fins, with one plate of water carying fluid, and a height of 0.42 meters.

Extends from Fluid.Interfaces.PartialTwoPortInterface (Partial model transporting fluid between two ports without storing mass or energy), Buildings.Fluid.Interfaces.LumpedVolumeDeclarations (Declarations for lumped volumes).

Parameters

TypeNameDefaultDescription
replaceable package MediumPartialMediumMedium in the component
IntegernEle5Number of elements used in the discretization
RealfraRad0.35Fraction radiant heat transfer
Realn1.24Exponent for heat transfer
Nominal condition
MassFlowRatem_flow_nominalabs(Q_flow_nominal/cp_nomina...Nominal mass flow rate [kg/s]
PowerQ_flow_nominal Nominal heating power (positive for heating) [W]
TemperatureT_a_nominal Water inlet temperature at nominal condition [K]
TemperatureT_b_nominal Water outlet temperature at nominal condition [K]
TemperatureTAir_nominal293.15Air temperature at nominal condition [K]
TemperatureTRad_nominalTAir_nominalRadiative temperature at nominal condition [K]
Initialization
MassFlowRatem_flow.start0Mass flow rate from port_a to port_b (m_flow > 0 is design flow direction) [kg/s]
Pressuredp.start0Pressure difference between port_a and port_b [Pa]
Assumptions
BooleanallowFlowReversalsystem.allowFlowReversal= true to allow flow reversal, false restricts to design direction (port_a -> port_b)
Advanced
MassFlowRatem_flow_small1E-4*abs(m_flow_nominal)Small mass flow rate for regularization of zero flow [kg/s]
BooleanhomotopyInitializationtrue= true, use homotopy method
Diagnostics
Booleanshow_Ttrue= true, if actual temperature at port is computed
Dynamics
Equations
DynamicsenergyDynamicsModelica.Fluid.Types.Dynamic...Formulation of energy balance
DynamicsmassDynamicsenergyDynamicsFormulation of mass balance
VolumeVWat5.8E-6*abs(Q_flow_nominal)Water volume of radiator [m3]
MassmDry0.0263*abs(Q_flow_nominal)Dry mass of radiator that will be lumped to water heat capacity [kg]
Initialization
AbsolutePressurep_startMedium.p_defaultStart value of pressure [Pa]
TemperatureT_startMedium.T_defaultStart value of temperature [K]
MassFractionX_start[Medium.nX]Medium.X_defaultStart value of mass fractions m_i/m [kg/kg]
ExtraPropertyC_start[Medium.nC]fill(0, Medium.nC)Start value of trace substances
ExtraPropertyC_nominal[Medium.nC]fill(1E-2, Medium.nC)Nominal value of trace substances. (Set to typical order of magnitude.)

Connectors

TypeNameDescription
FluidPort_aport_aFluid connector a (positive design flow direction is from port_a to port_b)
FluidPort_bport_bFluid connector b (positive design flow direction is from port_a to port_b)
HeatPort_aheatPortConHeat port for convective heat transfer with room air temperature
HeatPort_aheatPortRadHeat port for radiative heat transfer with room radiation temperature

Modelica definition

model RadiatorEN442_2 "Dynamic radiator for space heating"
   extends Fluid.Interfaces.PartialTwoPortInterface(
   showDesignFlowDirection = false,
   show_T=true,
   m_flow_nominal=abs(Q_flow_nominal/cp_nominal/(T_a_nominal-T_b_nominal)));
   extends Buildings.Fluid.Interfaces.LumpedVolumeDeclarations(
     final X_start = Medium.X_default,
     final C_start = fill(0, Medium.nC),
     final C_nominal = fill(1E-2, Medium.nC));

  parameter Integer nEle(min=1) = 5 
    "Number of elements used in the discretization";
  parameter Real fraRad(min=0, max=1) = 0.35 "Fraction radiant heat transfer";
  // Assumptions

  parameter Modelica.SIunits.Power Q_flow_nominal 
    "Nominal heating power (positive for heating)";
  parameter Modelica.SIunits.Temperature T_a_nominal 
    "Water inlet temperature at nominal condition";
  parameter Modelica.SIunits.Temperature T_b_nominal 
    "Water outlet temperature at nominal condition";
  parameter Modelica.SIunits.Temperature TAir_nominal = 293.15 
    "Air temperature at nominal condition";
  parameter Modelica.SIunits.Temperature TRad_nominal = TAir_nominal 
    "Radiative temperature at nominal condition";

  parameter Real n = 1.24 "Exponent for heat transfer";
  parameter Modelica.SIunits.Volume VWat = 5.8E-6*abs(Q_flow_nominal) 
    "Water volume of radiator";
  parameter Modelica.SIunits.Mass mDry = 0.0263*abs(Q_flow_nominal) 
    "Dry mass of radiator that will be lumped to water heat capacity";
  parameter Boolean homotopyInitialization = true "= true, use homotopy method";

  // Heat flow rates
  Modelica.SIunits.HeatFlowRate QCon_flow = heatPortCon.Q_flow 
    "Heat input into the water due to convective heat transfer with room air";
  Modelica.SIunits.HeatFlowRate QRad_flow = heatPortRad.Q_flow 
    "Heat input into the water due to radiative heat transfer with room";
  Modelica.SIunits.HeatFlowRate Q_flow = QCon_flow + QRad_flow 
    "Heat input into the water";

  // Heat ports
  Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a heatPortCon 
    "Heat port for convective heat transfer with room air temperature";
  Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a heatPortRad 
    "Heat port for radiative heat transfer with room radiation temperature";

  Fluid.MixingVolumes.MixingVolume[nEle] vol(
    redeclare each package Medium = Medium,
    each nPorts = 2,
    each V=VWat/nEle,
    each final m_flow_nominal = m_flow_nominal,
    each final energyDynamics=energyDynamics,
    each final massDynamics=energyDynamics,
    each final p_start=p_start,
    each final T_start=T_start,
    each final X_start=X_start,
    each final C_start=C_start) "Volume for fluid stream";
protected 
   parameter Modelica.SIunits.SpecificHeatCapacity cp_nominal=
      Medium.specificHeatCapacityCp(
        Medium.setState_pTX(Medium.p_default, T_a_nominal, Medium.X_default)) 
    "Specific heat capacity at nominal conditions";
   parameter Modelica.SIunits.HeatFlowRate QEle_flow_nominal[nEle](
      each fixed=false, each start=Q_flow_nominal/nEle) 
    "Nominal heating power of each element";
   parameter Modelica.SIunits.Temperature TWat_nominal[nEle](
      each fixed=false,
      start={T_a_nominal - i/nEle * (T_a_nominal-T_b_nominal) for i in 1:nEle}) 
    "Water temperature in each element at nominal conditions";
   parameter Modelica.SIunits.TemperatureDifference[nEle] dTRad_nominal(
    each fixed=false, start={T_a_nominal - i/nEle * (T_a_nominal-T_b_nominal) - TRad_nominal
    for i in 1:nEle}) 
    "Temperature difference for radiative heat transfer at nominal conditions";
   parameter Modelica.SIunits.TemperatureDifference[nEle] dTCon_nominal(
    each fixed=false, start={T_a_nominal - i/nEle * (T_a_nominal-T_b_nominal) - TAir_nominal
    for i in 1:nEle}) 
    "Temperature difference for convective heat transfer at nominal conditions";

   parameter Modelica.SIunits.ThermalConductance UAEle(fixed=false, min=0,
     start=Q_flow_nominal/((T_a_nominal+T_b_nominal)/2-((1-fraRad)*TAir_nominal+fraRad*TRad_nominal))/nEle) 
    "UA value at nominal condition for each element";

   final parameter Real k = if T_b_nominal > TAir_nominal then 1 else -1 
    "Parameter that is used to compute QEle_flow_nominal for heating or cooling mode";

   Modelica.Thermal.HeatTransfer.Components.HeatCapacitor[nEle] heaCap(
     each C=500*mDry/nEle,
     each T(start=T_start)) if 
       not (energyDynamics == Modelica.Fluid.Types.Dynamics.SteadyState) 
    "heat capacity of radiator metal";

   Buildings.HeatTransfer.Sources.PrescribedHeatFlow[nEle] preCon 
    "Heat input into radiator from convective heat transfer";
   Buildings.HeatTransfer.Sources.PrescribedHeatFlow[nEle] preRad 
    "Heat input into radiator from radiative heat transfer";

   Modelica.SIunits.TemperatureDifference dTCon[nEle] = heatPortCon.T .- vol.T 
    "Temperature difference for convective heat transfer";
   Modelica.SIunits.TemperatureDifference dTRad[nEle] = heatPortRad.T .- vol.T 
    "Temperature difference for radiative heat transfer";

  Modelica.Blocks.Sources.RealExpression QCon[nEle](y=if homotopyInitialization
         then homotopy(actual=(1 - fraRad) .* UAEle .* dTCon .*
        Buildings.Utilities.Math.Functions.regNonZeroPower(
        x=dTCon,
        n=n - 1,
        delta=0.05), simplified=(1 - fraRad) .* UAEle .* abs(dTCon_nominal) .^ (
        n - 1) .* dTCon) else (1 - fraRad) .* UAEle .* dTCon .*
        Buildings.Utilities.Math.Functions.regNonZeroPower(
        x=dTCon,
        n=n - 1,
        delta=0.05)) "Convective heat flow rate";

  Modelica.Blocks.Sources.RealExpression QRad[nEle](y=if homotopyInitialization
         then homotopy(actual=fraRad .* UAEle .* dTRad .*
        Buildings.Utilities.Math.Functions.regNonZeroPower(
        x=dTRad,
        n=n - 1,
        delta=0.05), simplified=fraRad .* UAEle .* abs(dTRad_nominal) .^ (n - 1)
         .* dTRad) else fraRad .* UAEle .* dTRad .*
        Buildings.Utilities.Math.Functions.regNonZeroPower(
        x=dTRad,
        n=n - 1,
        delta=0.05)) "Convective heat flow rate";

  Buildings.HeatTransfer.Sources.PrescribedHeatFlow preSumCon 
    "Heat input into radiator from convective heat transfer";
  Modelica.Blocks.Math.Sum sumCon(nin=nEle, k=-ones(nEle)) 
    "Sum of convective heat flow rate";
  Modelica.Blocks.Math.Sum sumRad(nin=nEle, k=-ones(nEle)) 
    "Sum of radiative heat flow rate";
  Buildings.HeatTransfer.Sources.PrescribedHeatFlow preSumRad 
    "Heat input into radiator from convective heat transfer";
initial equation 
  if T_b_nominal > TAir_nominal then
     assert(T_a_nominal > T_b_nominal,
       "In RadiatorEN442_2, T_a_nominal must be higher than T_b_nominal");
     assert(Q_flow_nominal > 0,
       "In RadiatorEN442_2, nominal power must be bigger than zero if T_b_nominal > TAir_nominal");
  else
     assert(T_a_nominal < T_b_nominal,
       "In RadiatorEN442_2, T_a_nominal must be lower than T_b_nominal");
     assert(Q_flow_nominal < 0,
       "In RadiatorEN442_2, nominal power must be smaller than zero if T_b_nominal < TAir_nominal");
  end if;
  TWat_nominal[1] = T_a_nominal - QEle_flow_nominal[1]/m_flow_nominal/
  Medium.specificHeatCapacityCp(
        Medium.setState_pTX(Medium.p_default, T_a_nominal, Medium.X_default));
  for i in 2:nEle loop
    TWat_nominal[i] = TWat_nominal[i-1] - QEle_flow_nominal[i]/m_flow_nominal/
    Medium.specificHeatCapacityCp(
        Medium.setState_pTX(Medium.p_default, TWat_nominal[i-1], Medium.X_default));
  end for;
  dTRad_nominal = TWat_nominal .- TRad_nominal;
  dTCon_nominal = TWat_nominal .- TAir_nominal;
  Q_flow_nominal = sum(QEle_flow_nominal);

  for i in 1:nEle loop
    QEle_flow_nominal[i] = k * UAEle * ((1-fraRad) *
                   Buildings.Utilities.Math.Functions.powerLinearized(x=k*dTRad_nominal[i],
                   n=n,
                   x0=0.1*k*(T_b_nominal-TRad_nominal))
                   + fraRad *
                   Buildings.Utilities.Math.Functions.powerLinearized(x=k*dTCon_nominal[i],
                   n=n,
                   x0=0.1*k*(T_b_nominal-TAir_nominal)));
   end for;


equation 
  connect(preCon.port, vol.heatPort);
  connect(preRad.port, vol.heatPort);
  connect(heaCap.port, vol.heatPort);
  connect(port_a, vol[1].ports[1]);
  connect(vol[nEle].ports[2], port_b);
  for i in 1:nEle-1 loop
    connect(vol[i].ports[2], vol[i+1].ports[1]);
  end for;
  connect(QCon.y, preCon.Q_flow);
  connect(sumCon.u, QCon.y);
  connect(sumCon.y, preSumCon.Q_flow);
  connect(preSumCon.port, heatPortCon);
  connect(QRad.y, preRad.Q_flow);
  connect(QRad.y, sumRad.u);
  connect(sumRad.y, preSumRad.Q_flow);
  connect(preSumRad.port, heatPortRad);
end RadiatorEN442_2;

Automatically generated Thu Oct 24 15:09:17 2013.