Direct Expansion Cycle Solver

In the cycle solver, the goal is to use all the physical component models together in order to obtain the performance of the combined cycle. The process is slightly different for each of the 4 standard configurations of ACHP, but there are many common themes among the configurations. ACHP has been designed to investigate direct expansion (DX) systems as well as systems employing a secondary working loop, in both cooling and heating modes.

A number of simplifying assumptions are employed at the cycle level:

  • Pressure drops in each component are calculated, but not used directly(see Pressure-Drop Correction for explanation)
  • There is no charge inventory in the compressor shell
  • The evaporator outlet refrigerant superheat is imposed (implies an expansion device with perfect superheat control)
  • No pressure drops except for in the components considered

Direct Expansion Cooling Mode

The most common and straightforward system is the direct expansion cooling mode system. The schematic of this cycle is shown in this figure:

(Source code, png, hires.png, pdf)

../_images/DXACSchematic.png

In addition to the geometry of each component, required inputs are the superheat at the outlet of the evaporator \Delta T_{sh,evap}, and either the total refrigerant charge m_r or the refrigerant subcooling at the outlet of the condenser \Delta T_{sc,cond}.

The primary independent variables in the cycle solver are the saturation temperatures of the refrigerant at evaporation and condensation, given by T_{evap} and T_{cond} respectively. For pure refrigerants the saturation temperature for a given pressure is the same for saturated liquid and saturated vapor. For pseudo-pure fluids and blends, the refrigerant dew temperature is used which corresponds to saturated vapor.

Before the full cycle-solver is run, the preconditioner in section Cycle Solver Preconditioner is used to get a good guess for the temperatures T_{evap} and T_{cond} using extremely simple cycle models.

Once the preconditioner has been run, preliminary values for T_{evap} and T_{cond} are known, and the first iteration of the cycle model may begin. Execution of the cycle model follows in the same direction as the refrigerant flow. The pressure drops through low-pressure and high-pressure parts of the system are assumed to be zero initially.

The cycle analysis begins with the vapor line which returns refrigerant vapor from the evaporator that is typically in the air duct back to the condensing unit outdoors. In the first iteration of the cycle solver, the mass flow rate of refrigerant in the vapor line is not known, but the compressor map is used to provide a reasonable guess value for the mass flow rate of refrigerant. The Line Set model is used to calculate the process from point 1 to point 2. The model is run with known flow rates and an inlet temperature of T_{evap}+\Delta T_{sh} in order to calculate the state at the inlet to the compressor.

The compressor compresses refrigerant from state point 2 to state point 3. The Compressor model is used which is based on an empirical correlation with superheat correction and it yields the outlet state 3 as well as the compressor mass flow rate, compressor electrical power, etc.. The compressor model requires the temperatures T_{evap} and T_{cond} as well as the compressor inlet superheat.

The condenser takes the superheated refrigerant at state point 3 and condenses it to a subcooled refrigerant at state point 4. The Condenser model is used to calculate the process, and this condenser model is based on a moving-boundary model. The condenser model requires T_{cond} as an input, among others.

After the condenser, the subcooled refrigerant passes through the liquid-line that takes refrigerant from the condenser to the indoor coil that is inside the ductwork in the home. The Line Set model is used to model the flow that passes from state point 4 to state point 5.

The expansion device then expands the refrigerant from the high-pressure side of the system to the low-pressure side of the system. In the current iteration of ACHP, no expansion device is included. Thus the expansion device is just a constant-enthalpy throttling device that takes refrigerant from the high-side pressure at state point 5 to the low-side pressure at state point 6.

At the outlet of the expansion device, the refrigerant passes into the evaporator at some two-phase quality. The Evaporator model is used to model the performance of the evaporator. In the evaporator, refrigerant is heated by the air-stream (which cools the air stream) from state point 6 back to state point 1’, which should be superheated (note: not necessarily so at intermediate iterations of cycle model) close to state point 1. If the values of the cycle independent variables T_{evap} and T_{cond} have been well selected, the state points 1 and 1’ will be coincident.

Once the cycle model has been run around the loop from state point 1 to state point 1’, the residuals can be calculated. The residuals are terms that when the cycle model has converged should all be equal to zero.

One of the residuals that is always active is an energy balance over the cycle. You left from state point 1, and you should hopefully arrive back there if energy is conserved in the cycle. Thus, the first residual is given by

(1)\vec \Delta_1=\dot m_r(h_1-h_{1'})

where h_1 and h_{1'} are the enthalpies of the refrigerant at state points 1 and 1’ respectively.

Since there are two independent variables, there must be a second constraint, which in this case is a charge-level constraint. Either the charge is constrained directly with an imposed charge level, or indirectly with an imposed refrigerant subcooling. Thus the second residual can be given by

(2)\vec \Delta_2=\left \lbrace \begin{array}{cc} m_r-m_{r,target} & (\mbox{Charge imposed}) \\ \Delta T_{sc,cond}-\Delta T_{sc,cond,target} & (\mbox{Subcooling imposed})\end{array} \right .

and a numerical solver is used to drive the residual vector \vec \Delta to sufficiently close to zero by altering T_{evap} and T_{cond} (see Multi-dimensional solver).

Pressure-Drop Correction

After the cycle model has iterated to convergence, the pressure drops are then considered. In reality, the pressure drop in each component results in a lower pressure at the inlet of the next component in the refrigerant loop. In terms of modeling, coupling pressure drop and the component models causes great numerical difficulties. The compromise that is employed instead is to run all the component models without pressure drop, but calculate the high-side pressure drop as the pressure drop of the condenser and liquid-line:

(3)\Delta p_{high}=\Delta p_{cond}+\Delta p_{liquid-line}

and similarly, the low-side pressure drop is defined by

(4)\Delta p_{low}=\Delta p_{evap}+\Delta p_{vapor-line}

These pressure drops are then employed to shift the saturation temperatures used in the compressor map in order to yield less refrigerant mass flow rate, and a higher compressor power.

The new effective compressor suction and discharge pressures are

(5)p_{evap}^*=p(T_{evap})-\Delta p_{low}

p_{cond}^*=p(T_{cond})+\Delta p_{high}

and the new, effective compressor dew temperatures are

(6)T_{evap}^*=T(p=p_{evap},x=1)

T_{cond}^*=T(p=p_{cond},x=1)

This calculated pressure drop is used until the model reaches convergence again, at which point the pressure drop tems are updated, and the model is run again. This process continues until the imposed low- and high-side pressure drops are equal to the pressure drop terms calculated from the converged cycle model.

Cycle Performance Parameters

The common metrics of system efficiency are

(7)COP=\frac{\dot Q_{evap}}{\dot W_{comp}}

COSP=\frac{\dot Q_{evap}-\dot W_{fan,evap}}{\dot W_{comp}+\dot W_{fan,evap}+\dot W_{fan,cond}}

from Cycle import DXCycleClass

#Instantiate the cycle class
Cycle=DXCycleClass()

#--------------------------------------
#       Cycle parameters
#--------------------------------------
Cycle.Verbosity = 0 #the idea here is to have different levels of debug output
Cycle.ImposedVariable = 'Subcooling'                                                           
Cycle.DT_sc_target = 7.0
Cycle.Mode='AC' 
Cycle.Ref='R410A'

#--------------------------------------
#       Compressor parameters
#--------------------------------------

#A 3 ton cooling capacity compressor map 
M=[217.3163128,5.094492028,-0.593170311,4.38E-02,-2.14E-02,
    1.04E-02,7.90E-05,-5.73E-05,1.79E-04,-8.08E-05]
P=[-561.3615705,-15.62601841,46.92506685,-0.217949552,
    0.435062616,-0.442400826,2.25E-04,2.37E-03,-3.32E-03,2.50E-03]

params={
    'M':M,
    'P':P,
    'Ref':Cycle.Ref,  #Refrigerant
    'fp':0.0, #Fraction of electrical power lost as heat to ambient
    'Vdot_ratio': 1.0, #Displacement Scale factor                  
    'Verbosity': 0, # How verbose should the debugging be [0-10]
  }
Cycle.Compressor.Update(**params)                                          

#--------------------------------------
#      Condenser parameters 
#--------------------------------------
Cycle.Condenser.Fins.Tubes.NTubes_per_bank=24  #number of tubes per bank=row
Cycle.Condenser.Fins.Tubes.Nbank=1             #number of banks/rows
Cycle.Condenser.Fins.Tubes.Ncircuits=3
Cycle.Condenser.Fins.Tubes.Ltube=2.252
Cycle.Condenser.Fins.Tubes.OD=0.00913
Cycle.Condenser.Fins.Tubes.ID=0.00849
Cycle.Condenser.Fins.Tubes.Pl=0.0191  #distance between center of tubes in flow direction                                                
Cycle.Condenser.Fins.Tubes.Pt=0.0254  #distance between center of tubes orthogonal to flow direction

Cycle.Condenser.Fins.Fins.FPI=25      #Number of fins per inch
Cycle.Condenser.Fins.Fins.Pd=0.001    #2* amplitude of wavy fin
Cycle.Condenser.Fins.Fins.xf=0.001    #1/2 period of fin
Cycle.Condenser.Fins.Fins.t=0.00011   #Thickness of fin material
Cycle.Condenser.Fins.Fins.k_fin=237   #Thermal conductivity of fin material

Cycle.Condenser.Fins.Air.Vdot_ha=1.7934 #rated volumetric flowrate
Cycle.Condenser.Fins.Air.Tmean=308.15   
Cycle.Condenser.Fins.Air.Tdb=308.15     #Dry Bulb Temperature
Cycle.Condenser.Fins.Air.p=101.325      #Air pressure
Cycle.Condenser.Fins.Air.RH=0.51        #Relative Humidity
Cycle.Condenser.Fins.Air.RHmean=0.51
Cycle.Condenser.Fins.Air.FanPower=260

Cycle.Condenser.Ref=Cycle.Ref
Cycle.Condenser.Verbosity=0
    
#--------------------------------------
#     Evaporator Parameters
#--------------------------------------
Cycle.Evaporator.Fins.Tubes.NTubes_per_bank=32
Cycle.Evaporator.Fins.Tubes.Nbank=3
Cycle.Evaporator.Fins.Tubes.Ltube=0.452
Cycle.Evaporator.Fins.Tubes.OD=0.00913
Cycle.Evaporator.Fins.Tubes.ID=0.00849
Cycle.Evaporator.Fins.Tubes.Pl=0.0191
Cycle.Evaporator.Fins.Tubes.Pt=0.0254
Cycle.Evaporator.Fins.Tubes.Ncircuits=5

Cycle.Evaporator.Fins.Fins.FPI=14.5
Cycle.Evaporator.Fins.Fins.Pd=0.001
Cycle.Evaporator.Fins.Fins.xf=0.001
Cycle.Evaporator.Fins.Fins.t=0.00011
Cycle.Evaporator.Fins.Fins.k_fin=237

Cycle.Evaporator.Fins.Air.Vdot_ha=0.56319
Cycle.Evaporator.Fins.Air.Tmean=297.039
Cycle.Evaporator.Fins.Air.Tdb=297.039
Cycle.Evaporator.Fins.Air.p=101.325
Cycle.Evaporator.Fins.Air.RH=0.5
Cycle.Evaporator.Fins.Air.RHmean=0.5
Cycle.Evaporator.Fins.Air.FanPower=438
       
Cycle.Evaporator.Ref=Cycle.Ref
Cycle.Evaporator.Verbosity=0
Cycle.Evaporator.DT_sh=5

# ----------------------------------
#       Line Set Parameters
# ----------------------------------
params={
    'L':7.6,
    'k_tube':0.19,
    't_insul':0.02,
    'k_insul':0.036,
    'T_air':297,
    'Ref': Cycle.Ref,
    'h_air':0.0000000001,
}

Cycle.LineSetSupply.Update(**params)
Cycle.LineSetReturn.Update(**params)
Cycle.LineSetSupply.OD=0.009525
Cycle.LineSetSupply.ID=0.007986
Cycle.LineSetReturn.OD=0.01905
Cycle.LineSetReturn.ID=0.017526

#Now solve
from time import time
t1=time()
Cycle.PreconditionedSolve()
print 'Took '+str(time()-t1)+' seconds to run Cycle model'
print 'Cycle coefficient of system performance is '+str(Cycle.COSP)
print 'Cycle refrigerant charge is '+str(Cycle.Charge)+' kg'

which should yield the output, when run, of

Took 2.01799988747 seconds to run Cycle model
Cycle coefficient of system performance is 3.22662342401
Cycle refrigerant charge is 2.03698681412 kg

Direct Expansion Heating Mode

The heat pump configuration of the system is as shown here:

(Source code, png, hires.png, pdf)

../_images/DXHPSchematic.png

Physically, reversing valves are used to switch the mode of the system and the directions of the flows. What was the condenser of the air conditioning system becomes the evaporator of the heat pump and vice versa, and the line sets are configured in a slightly different way. Other than that, the analysis of the heat pump is directly analogous to that of the Direct Expansion Cooling Mode system.

A Preconditioner is used to get approximate values for T_{evap} and T_{cond}, and using these values (which are iteratively modified using numerical methods), the solution for the cycle performance is found.

As with the cooling mode, the cycle analysis follows the refrigerant flow path around the loop.

Beginning at the outlet of the evaporator, state point 1 is known because T_{evap} and \Delta T_{sh} are known. Thus the compressor model is used directly to calculate the electrical power, refrigerant mass flow rate and state point 2.

The Line Set model is then applied to the flow from the outlet of the compressor at state point 2 to the inlet of the condenser at state point 3.

The Condenser model is used to model the condensing process from state point 3 to a subcooled state at state point 4.

The Line Set model is used to model the flow of subcooled refrigerant at state point 4 back to the expansion device at state point 5.

As in cooling mode, the expansion device is assumed to be an ideal expansion device, which means that the working process is a constant-enthalpy expansion from state point 5 to state point 6.

The Evaporator model is then used to calculate the evaporation process of refrigerant from state point 6 to state point 1’.

As in cooling mode, the residual vector is given by

(8)\vec \Delta=[\vec \Delta_1, \vec \Delta_2]

with

(9)\vec \Delta_1=\dot m_r(h_1-h_{1'})

\vec \Delta_2=\left \lbrace \begin{array}{cc} m_r-m_{r,target} & (\mbox{Charge imposed}) \\ \Delta T_{sc}-\Delta T_{sc,target} & (\mbox{Subcooling imposed})\end{array} \right .

The set of T_{evap} and T_{cond} which solve the residual equations are obtained by a multi-dimensional solver (see Multi-dimensional solver).

Cycle Performance Parameters

The common metrics of system efficiency are

(10)COP=\frac{\dot Q_{cond}}{\dot W_{comp}}

COSP=\frac{\dot Q_{cond}+\dot W_{fan,cond}}{\dot W_{comp}+\dot W_{fan,evap}+\dot W_{fan,cond}}

Minimal working example

from Cycle import DXCycleClass,F2K

#Instantiate the class
Cycle=DXCycleClass()

#--------------------------------------
#       Cycle parameters
#--------------------------------------
Cycle.Verbosity = 0 #the idea here is to have different levels of debug output
Cycle.ImposedVariable = 'Subcooling' #or it could be 'Charge'
Cycle.DT_sc_target = 7.0
#Cycle.Charge_target = 3.3 #uncomment for use with imposed charge
Cycle.Mode='HP' 
Cycle.Ref='R410A'

#--------------------------------------
#       Compressor parameters
#--------------------------------------

#A few 3 ton cooling capacity compressor maps 
M=[217.3163128,5.094492028,-0.593170311,4.38E-02,
    -2.14E-02,1.04E-02,7.90E-05,-5.73E-05,1.79E-04,-8.08E-05]
P=[-561.3615705,-15.62601841,46.92506685,-0.217949552,
    0.435062616,-0.442400826,2.25E-04,2.37E-03,-3.32E-03,2.50E-03]

params={
'M':M,
'P':P,
'Ref':Cycle.Ref,                                                              #refrigerant
'fp':0.0, #Fraction of electrical power lost as heat to ambient            #shell heat loss
'Vdot_ratio': 1.0, #Displacement Scale factor                               #up- or downsize compressor (1=original)
'Verbosity': 0, # How verbose should the debugging statements be [0 to 10]                        
}
Cycle.Compressor.Update(**params)                                          

#--------------------------------------
#      Condenser parameters 
#--------------------------------------
Cycle.Condenser.Fins.Tubes.NTubes_per_bank=32
Cycle.Condenser.Fins.Tubes.Nbank=3
Cycle.Condenser.Fins.Tubes.Ncircuits=6
Cycle.Condenser.Fins.Tubes.Ltube=0.452
Cycle.Condenser.Fins.Tubes.OD=0.009525
Cycle.Condenser.Fins.Tubes.ID=0.0089154
Cycle.Condenser.Fins.Tubes.Pl=0.0254
Cycle.Condenser.Fins.Tubes.Pt=0.0219964

Cycle.Condenser.Fins.Fins.FPI=14.5
Cycle.Condenser.Fins.Fins.Pd=0.001
Cycle.Condenser.Fins.Fins.xf=0.001
Cycle.Condenser.Fins.Fins.t=0.00011
Cycle.Condenser.Fins.Fins.k_fin=237

Cycle.Condenser.Fins.Air.Vdot_ha=0.5663
Cycle.Condenser.Fins.Air.Tmean=F2K(70)
Cycle.Condenser.Fins.Air.Tdb=F2K(70)
Cycle.Condenser.Fins.Air.p=101.325
Cycle.Condenser.Fins.Air.RH=0.51
Cycle.Condenser.Fins.Air.RHmean=0.51
Cycle.Condenser.Fins.Air.FanPower=438

Cycle.Condenser.Ref=Cycle.Ref
Cycle.Condenser.Verbosity=0

#--------------------------------------
#           Evaporator parameters
#--------------------------------------
Cycle.Evaporator.Fins.Tubes.NTubes_per_bank=41  #number of tubes per bank=row
Cycle.Evaporator.Fins.Tubes.Nbank=1             #number of banks/rows
Cycle.Evaporator.Fins.Tubes.Ncircuits=5
Cycle.Evaporator.Fins.Tubes.Ltube=2.286
Cycle.Evaporator.Fins.Tubes.OD=0.007
Cycle.Evaporator.Fins.Tubes.ID=0.0063904
Cycle.Evaporator.Fins.Tubes.Pl=0.0191  #distance between center of tubes in flow direction                                                
Cycle.Evaporator.Fins.Tubes.Pt=0.0222  #distance between center of tubes orthogonal to flow direction

Cycle.Evaporator.Fins.Fins.FPI=25      #Number of fins per inch
Cycle.Evaporator.Fins.Fins.Pd=0.001    #2* amplitude of wavy fin
Cycle.Evaporator.Fins.Fins.xf=0.001    #1/2 period of fin
Cycle.Evaporator.Fins.Fins.t=0.00011   #Thickness of fin material
Cycle.Evaporator.Fins.Fins.k_fin=237   #Thermal conductivity of fin material

Cycle.Evaporator.Fins.Air.Vdot_ha=1.7934 #rated volumetric flowrate
Cycle.Evaporator.Fins.Air.Tmean=F2K(47)   
Cycle.Evaporator.Fins.Air.Tdb=F2K(47)     #Dry Bulb Temperature
Cycle.Evaporator.Fins.Air.p=101.325      #Air pressure
Cycle.Evaporator.Fins.Air.RH=0.51        #Relative Humidity
Cycle.Evaporator.Fins.Air.RHmean=0.51
Cycle.Evaporator.Fins.Air.FanPower=160

Cycle.Evaporator.Ref=Cycle.Ref
Cycle.Evaporator.Verbosity=0
Cycle.Evaporator.DT_sh=5

# ----------------------------------
#       Line Set parameters
# ----------------------------------
params={
'L':7.6,
'k_tube':0.19,
't_insul':0.02,
'k_insul':0.036,
'T_air':297,
'Ref': Cycle.Ref,
'h_air':6,
}

Cycle.LineSetSupply.Update(**params)
Cycle.LineSetReturn.Update(**params)
Cycle.LineSetSupply.OD=0.01905 
Cycle.LineSetSupply.ID=0.017526
Cycle.LineSetReturn.OD=0.009525
Cycle.LineSetReturn.ID=0.007986

#Now solve
from time import time
t1=time()
Cycle.PreconditionedSolve()

#Outputs
print 'Took '+str(time()-t1)+' seconds to run Cycle model'
print 'Cycle coefficient of system performance is '+str(Cycle.COSP)
print 'Cycle refrigerant charge is '+str(Cycle.Charge)+' kg'

which should yield the output, when run, of

Took 2.96899986267 seconds to run Cycle model
Cycle coefficient of system performance is 3.662196962
Cycle refrigerant charge is 1.72330956768 kg

Cycle Solver Code Documentation

class Cycle.DXCycleClass[source]
Calculate(DT_evap, DT_cond)[source]

Inputs are differences in temperature [K] between HX air inlet temperature and the dew temperature for the heat exchanger.

Required Inputs:
DT_evap:
Difference in temperature [K] between evaporator air inlet temperature and refrigerant dew temperature
DT_cond:
Difference in temperature [K] between condenser air inlet temperature and refrigeant dew temperature
OutputList()[source]

Return a list of parameters for this component for further output

It is a list of tuples, and each tuple is formed of items:
[0] Description of value [1] Units of value [2] The value itself
PreconditionedSolve()[source]

Solver that will precondition by trying a range of DeltaT until the model can solve, then will kick into 2-D Newton Raphson solve

The two input variables for the system solver are the differences in temperature between the inlet air temperature of the heat exchanger and the dew temperature of the refrigerant. This is important for refrigerant blends with temperature glide during constant-pressure evaporation or condensation. Good examples of common working fluid with glide would be R404A or R410A.

Nomenclature

Variable Description
COP Coefficient of Performance [-]
COSP Coefficient of System Performance [-]
h_1 Enthalpy at outlet of evaporator [J/kg]
h_{1'} Enthalpy after going around the cycle [J/kg]
\dot m_r Refrigerant mass flow rate [kg/s]
m_r Model-predicted charge [kg]
m_{r,target} Target refrigerant charge in system [kg]
T_{evap} Evaporating (dewpoint) temperature [K]
T_{cond} Condensing (dewpoint) temperature [K]
T_{evap}^* Effective evaporating temperature [K]
T_{cond}^* Effective condensing temperature [K]
\Delta T_{sh} Evaporator outlet superheat [K]
\Delta T_{sc,cond} Condenser outlet subcooling [K]
\Delta T_{sc,cond,target} Condenser outlet subcooling target [K]
\Delta p_{cond} Pressure drop in condenser [kPa]
\Delta p_{evap} Pressure drop in evaporator [kPa]
\Delta p_{high} Pressure drop on high-pressure side of system [kPa]
\Delta p_{liquid-line} Pressure drop in liquid line [kPa]
\Delta p_{low} Pressure drop on low-pressure side of system [kPa]
\Delta p_{vapor-line} Pressure drop in vapor line [kPa]
p_{evap}^* Effective evaporation saturation pressure [kPa]
p_{cond}^* Effective condensing saturation pressure [kPa]
\dot Q_{evap} Evaporator heat transfer rate [W]
\dot W_{fan,evap} Evaporator fan power [W]
\dot W_{fan,cond} Condenser fan power [W]
\dot W_{comp} Compressor power input [W]
\vec \Delta_1 Residual vector [varied]