Tube-Fin Heat Exchanger:

This tutorial provides a detailed example of modeling with a TubeFinHEX which is one of the key pieces of equipment for HVAC systems.

Overview:

Internal sub-components:

A tube-fin heat exchanger consists of:

  • A hot thermo-fluid flow stream
  • A cold thermo-fluid flow stream
  • A wall across which heat is transferred from the hot stream to the cold stream

The rest of the tutorial follows the case of an Evaporator (where the hot stream is the warm air and the cold stream is the refrigerant). The refrigerant stream passes through a closed pipe on the tube-side; the air stream is open on the fin-side.

1. Geometric Configuration, Circuitry and Parameters:

Heat exchanger vendors may provide information about its internal geometric configuration. These are passed into the model using the TubeFinHEXGeo type.

Number of refrigerant circuits:

The incoming refrigerant stream is split into a number of parallel paths which is specified in the nCir field of TubeFinHEXGeo. The figures below show the case of a single circuit and two circuits respectively. Thus, in a single circuit, the refrigerant effectively flows through one long pipe which results in a large pressure drop. Lower pressure drops are achieved by splitting the refrigerant stream into multiple circuits.

Single_CircuitDouble_Circuit

Number of tubes:

For simplicity, this tutorial considers the case where the number of tubes per circuit is the same as specified in the nTube field of TubeFinHEXGeo (the double circuit figure above show 6 tubes per circuit). Each tube is of length Lt of inner diameter Di and outer diameter Do.

Other parameters:

Other parameters such as the fin surface areas and thicknesses, thermal conductivities, specific heat capacities etc. may be provided but these are typically determined using a model calibration work flow.

2. Dynamic model:

Refrigerant-side:

A 1D physics-based dynamical model is developed within the TubeFinHEX function which returns an ODESystem. The refrigerant enters through refPort_a is split into nCir circuits by the a ref_distributor modeled using a RefrigerantSplitterMerger. One circuit is modeled using a StraightPipe which sets up a distributed dynamic mass and energy as well as a static momentum balance using a Finite Volume approach. Following this, a ref_merger is used to merge the refrigerant streams and performs the book-keeping on refrigerant mass flow rates and duties across the circuits. The refrigerant exits through refPort_b.

Air-side:

The air-side flow is also divided into a number of channels (taken here to also be the same as the number of refrigerant circuits for convenience) each of which is modeled using a HEXCrossflow. Steady-state is assumed with replaceable heat transfer and flow resistance models. A MultiComponentFluidPort is used for the air inlet and outlet. The geometric parameters are used to populate a AirChannelGeo record.

Wall-side:

A replaceable wall model is used for the heat transfer between the air and refrigerant streams with the defaults as a OneTempWall1D. The geometric parameters are used to populate aTubeWallGeoProps record.

3. Initialization:

It is essential to pass in reasonable initial guesses for certain states of by using the init kwarg of the TubeFinHEX function and making use of the TubeFinHEXInit function. This takes in initial values for a number of states such as the refrigerant side pressure, specific enthalpy, and mass flow rate, the wall and air-side temperatures etc. These initial values are propagated to the sub-components appropriately. The init_DAE function is used to determine a consistent set of initial values close to the provided guesses.

Example:

Copy-pasteable code:

using HVAC
using ModelingToolkit
using RefpropSplines
using OrdinaryDiffEq
using Plots
using Test

import ModelingToolkit: namespace_variables
import HVAC: moistair, relativeHumidity_PTX, P_amb, inspect_sol_inits

load_refrigerant("R32.yaml")
p_in, p_out = 8e5, 7.9e5
h_in, h_out = 230e3, 430e3
ref_mdot_out= 0.01

air_mdot_in = 0.01
air_Tin = 300
air_ϕin = 0.5

nTube = 4
Lt = 20
nSeg = 1

geo = TubeFinHEXGeo(;nTube = nTube,Lt=Lt, nSeg=nSeg)
init_hex = TubeFinHEXInit(;
                circuit_h_in_start=h_in,        
                circuit_h_out_start=h_out,       
                circuit_p_in_start=p_in,           
                circuit_p_out_start=p_out,       
                circuit_mdot_start=ref_mdot_out,       
                Twall_in_start=from_degC(20),    
                Twall_out_start=from_degC(20),   
                Twall_start=from_degC(20),
                air_mdot_start = air_mdot_in,
                air_Tin_start=air_Tin,       
                T_air_out_ref_port_a_start = from_degC(18),
                T_air_out_ref_port_b_start = from_degC(20),
                Q_init = 50)



@named heatExchanger = TubeFinHEX(model_structure=av_b, htc_air=50, init=init_hex, geo=geo)

defs = Dict()

#@named refSource = Boundary_Ph(P_in=p_in, h_in=h_in)

@named refSource = MassFlowSource_h(m_flow_in=ref_mdot_out, h_in=h_in)
@named refSink=Boundary_Ph(P_in=p_out, h_in=h_out)
#@named refSink = MassFlowSource_h(m_flow_in=-1*ref_mdot_out, h_in=h_out)

@named airSource = MassFlowSource_Tϕ(T_in=air_Tin, ϕ_in=air_ϕin, m_flow_in=air_mdot_in)
@named airSink = Boundary_PTϕ(; P_in = P_amb - 1e-2, T_in = from_degC(18), ϕ_in=relativeHumidity_PTX(P_amb - 1e-2, from_degC(18), moistair.Xi_ref, moistair))

eqns = [
        connect(refSource.port, heatExchanger.refPort_a)
        connect(heatExchanger.refPort_b, refSink.port1)

        connect(airSource.port, heatExchanger.airPort_a)
        connect(heatExchanger.airPort_b, airSink.port1)
]

@named sys = ODESystem(eqns, t, [], []; defaults=defs)
sysModel = compose(sys, heatExchanger, refSource, refSink, airSource, airSink)
sysRed = structural_simplify(sysModel)

@info "Reduced system has $(length(states(sysRed))) states, full system has $(length(states(sysModel))) states."


defs = ModelingToolkit.defaults(sysRed)
unset_ps = setdiff(parameters(sysRed), keys(defs))
unset_vs = setdiff(states(sysRed), keys(defs))

prob = ODEProblem(sysRed, defs, (0.0, 2.0))
prob_1 = init_DAE(prob, verbose=true)

sol = solve(prob_1, Rodas4(autodiff=false), tstops=0.1)
@test SciMLBase.successful_retcode(sol)

(; m_flows_plot, q_flows_plot, T_plot, h_plot, p_plot, ρ_plot, xq_plot, dm_flow_plot, du_flow_plot, H_flow_plot) = plot_comp(:pipe, heatExchanger.circuit, nSeg*nTube, sol)

Step-by-step explanation:

1. Import required packages:

  • HVAC: Provides the pre-built components
  • RefpropSplines: Provides the refrigerant properties
  • ModelingToolkit: Provides the acausal modeling framework
  • OrdinaryDiffEq: Provides the solver for the ODEs
using HVAC
using ModelingToolkit
using RefpropSplines
using OrdinaryDiffEq
using Plots
using Test

import ModelingToolkit: namespace_variables
import HVAC: moistair, relativeHumidity_PTX, P_amb, inspect_sol_inits

2. Load the refrigerant:

This example uses R32 and uses a fast and accurate spline-based interpolation for this case. 

load_refrigerant("R32.yaml")
Refrigerant model constructed using spline interpolations

3. Set up the refrigerant and air side initial conditions:

p_in, p_out = 8e5, 7.9e5
h_in, h_out = 230e3, 430e3
ref_mdot_out= 0.01

air_mdot_in = 0.01
air_Tin = 300
air_ϕin = 0.5
init_hex = TubeFinHEXInit(;
                circuit_h_in_start=h_in,
                circuit_h_out_start=h_out,
                circuit_p_in_start=p_in,
                circuit_p_out_start=p_out,
                circuit_mdot_start=ref_mdot_out,
                Twall_in_start=from_degC(20),
                Twall_out_start=from_degC(20),
                Twall_start=from_degC(20),
                air_mdot_start = air_mdot_in,
                air_Tin_start=air_Tin,
                T_air_out_ref_port_a_start = from_degC(18),
                T_air_out_ref_port_b_start = from_degC(20),
                Q_init = 50)
(circuit_h_in_start = 230000.0, circuit_h_out_start = 430000.0, circuit_p_in_start = 800000.0, circuit_p_out_start = 790000.0, circuit_mdot_start = 0.01, Twall_in_start = 293.15, Twall_out_start = 293.15, Twall_start = 293.15, air_mdot_start = 0.01, air_Tin_start = 300, air_Ploss_start = 1000, T_air_out_ref_port_a_start = 291.15, T_air_out_ref_port_b_start = 293.15, Q_init = 50)

4. Set up the refrigerant geometry

nTube = 4
Lt = 20
nSeg = 1

geo = TubeFinHEXGeo(;nTube = nTube,Lt=Lt, nSeg=nSeg)
(nTube = 4, nCir = 4, nHeader = 2, nRow = 1, nSeg = 1, Lt = 20, Di = 0.009449, Do = 0.010058, FPI = 12, Pt = 0.0254, Pr = 0.01905, tf = 0.00011, height_fin = 0, width_fin = 0, n_fin = 0, helix_angle = 0, apex_angle = 0, nTubePerRow = 4, Ap = 2.396482472277581, As = 30.57006460809771, Ai = 2.3747927187015967, Aface = 2.032, tubemass = 6.64320495775942, finmass = 4.53965459430251, meanCp = 592.8451527949595, meanRho = 14981.842445956869, nTubeOfCircuit = [12, 12, 12, 12], TubeNoSequenceOfCircuit = [[1, 36, 35, 34, 33, 17, 18, 19, 20, 4, 3, 2, 1, 2], [1, 40, 39, 38, 37, 21, 22, 23, 24, 8, 7, 6, 5, 2], [1, 44, 43, 42, 41, 25, 26, 27, 28, 12, 11, 10, 9, 2], [1, 48, 47, 46, 45, 29, 30, 31, 32, 16, 15, 14, 13, 2]], kt = 400, kf = 237, rhot = 8900, rhof = 2700, cpt = 385, cpf = 897)

5. Construct and connect the component model including source and sink terms:

@named heatExchanger = TubeFinHEX(model_structure=av_b, htc_air=50, init=init_hex, geo=geo)

defs = Dict()

#@named refSource = Boundary_Ph(P_in=p_in, h_in=h_in)

@named refSource = MassFlowSource_h(m_flow_in=ref_mdot_out, h_in=h_in)
@named refSink=Boundary_Ph(P_in=p_out, h_in=h_out)
#@named refSink = MassFlowSource_h(m_flow_in=-1*ref_mdot_out, h_in=h_out)

@named airSource = MassFlowSource_Tϕ(T_in=air_Tin, ϕ_in=air_ϕin, m_flow_in=air_mdot_in)
@named airSink = Boundary_PTϕ(; P_in = P_amb - 1e-2, T_in = from_degC(18), ϕ_in=relativeHumidity_PTX(P_amb - 1e-2, from_degC(18), moistair.Xi_ref, moistair))

eqns = [
        connect(refSource.port, heatExchanger.refPort_a)
        connect(heatExchanger.refPort_b, refSink.port1)

        connect(airSource.port, heatExchanger.airPort_a)
        connect(heatExchanger.airPort_b, airSink.port1)
]

@named sys = ODESystem(eqns, t, [], []; defaults=defs)
sysModel = compose(sys, heatExchanger, refSource, refSink, airSource, airSink)
sysRed = structural_simplify(sysModel)

@info "Reduced system has $(length(states(sysRed))) states, full system has $(length(states(sysModel))) states."


defs = ModelingToolkit.defaults(sysRed)
unset_ps = setdiff(parameters(sysRed), keys(defs))
unset_vs = setdiff(states(sysRed), keys(defs))
Any[]

6. Formulate, initialize and solve the DAE problem:

prob = ODEProblem(sysRed, defs, (0.0, 2.0))
prob_1 = init_DAE(prob, verbose=true)

sol = solve(prob_1, Rodas4(autodiff=false), tstops=0.1)
@test SciMLBase.successful_retcode(sol)

(; m_flows_plot, q_flows_plot, T_plot, h_plot, p_plot, ρ_plot, xq_plot, dm_flow_plot, du_flow_plot, H_flow_plot) = plot_comp(:pipe, heatExchanger.circuit, nSeg*nTube, sol)
(m_flows_plot = Plot{Plots.GRBackend() n=7}, q_flows_plot = Plot{Plots.GRBackend() n=4}, T_plot = Plot{Plots.GRBackend() n=4}, h_plot = Plot{Plots.GRBackend() n=4}, p_plot = Plot{Plots.GRBackend() n=4}, ρ_plot = Plot{Plots.GRBackend() n=4}, xq_plot = Plot{Plots.GRBackend() n=4}, dm_flow_plot = Plot{Plots.GRBackend() n=4}, du_flow_plot = Plot{Plots.GRBackend() n=4}, H_flow_plot = Plot{Plots.GRBackend() n=5})

7. Plot the results:

The plot below shows the quality of the refrigerant from approximately zero on entering the evaporator to approximately 1 on leaving the evaporator. 

plot(xq_plot)
plot(m_flows_plot)
plot(q_flows_plot)
plot(T_plot)
plot(h_plot)
plot(p_plot)
plot(ρ_plot)
plot(dm_flow_plot)
plot(du_flow_plot)
plot(H_flow_plot)