# Computational psychrometric analysis as a control problem: case of cooling and dehumidification systems

Christian Ghiaus

► To cite this version:

Christian Ghiaus. Computational psychrometric analysis as a control problem: case of cooling and dehumidification systems. Journal of Building Performance Simulation, 2022, 15 (1), pp.21-38. 10.1080/19401493.2021.1995498. hal-03484064

HAL Id: hal-03484064

https://hal.science/hal-03484064

Submitted on 16 Dec 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Computational psychrometric analysis as a control problem: case of cooling and dehumidification systems

Christian Ghiaus, Univ. Lyon, CNRS, INSA Lyon, Université Claude Bernard Lyon1, CETHIL UMR5008, F-69621, Villeurbanne, France

## Abstract

Psychrometric chart is a basic tool for analysis of processes in air-conditioning systems. While psychrometric calculators and numerical psychrometric charts are widely available, there is a lack of numerical algorithms for solving the sizing problem by using psychrometric analysis. After introducing a classification of modelling problems in direct and inverse (based on the relation between physical and computational causality), this paper defines the design problem as a set of inverse problems of control and parameter optimization. A non-linear model is obtained by assembling the models of the elementary processes. The paper proposes to solve the direct and control problems by using a method similar to Newton-Raphson's, and the parameter optimization problem by using least-squares. Open-source implementation, published on Zenodo repository, and interactive, reproducible environment, available on Binder web service, accompany this paper (DOI: 10.5281/zenodo.5236450).

Keywords: inverse problem, HVAC, design, sizing, steady-state, Python

Python code and Jupyter notebook:

Jupyter interactive figures:

DOI

10.5281/zenodo.5236450

<figure><img src="image_2.png" /></figure>

launch

binder

DOI

10.5281/zenodo.5579167

<figure><img src="image_3.png" /></figure>

launch

binder
## Nomenclature

<table><tr><td>A</td><td>surface area, m<sup>2</sup></td></tr><tr><td>AHU</td><td>air handling unit</td></tr><tr><td>C</td><td>constant</td></tr><tr><td>CAV</td><td>constant air volume</td></tr><tr><td>K</td><td>controller gain</td></tr><tr><td>H&#x0307;<sub>i</sub></td><td>enthalpy rate of the air in state i, W</td></tr><tr><td>M</td><td>molar mass, kg/kmol</td></tr><tr><td>Q&#x0307;</td><td>heat flow rate, W</td></tr><tr><td>R</td><td>universal gas constant, J/K kmol</td></tr><tr><td>SHR</td><td>sensible heat rate, -</td></tr><tr><td>T</td><td>temperature, K</td></tr><tr><td>U</td><td>overall heat transfer coefficient, W/m<sup>2</sup> K</td></tr><tr><td>VAV</td><td>variable air volume</td></tr><tr><td>c</td><td>specific heat capacity of the dry air, J/kg K</td></tr><tr><td>c<sub>w</sub></td><td>specific heat capacity of the water vapour, J/kg K</td></tr><tr><td>d</td><td>derivative</td></tr><tr><td>f</td><td>function representing the saturation curve</td></tr><tr><td>h<sub>i</sub></td><td>specific enthalpy of dry air in state i, J/kg</td></tr><tr><td>l</td><td>specific latent heat for water vaporization, J/kg</td></tr><tr><td>l<sub>k</sub></td><td>length of thermal bridge k, m</td></tr><tr><td>m&#x0307;</td><td>mass flow rate of dry air, kg/s</td></tr><tr><td>n</td><td>molar fraction</td></tr><tr><td>p</td><td>pressure, Pa</td></tr><tr><td>v</td><td>unspecified variable</td></tr><tr><td>w<sub>i</sub></td><td>humidity ratio of moist air in the state i, kg<sub>vapor</sub>/kg<sub>dry air</sub></td></tr></table>

### Greek letters

<table><tr><td>Δ</td><td>difference operator</td></tr><tr><td>α</td><td>mixing ratio</td></tr><tr><td>β</td><td>bypass factor</td></tr></table>

<table><tr><td>ε</td><td>arbitrarily small number</td></tr><tr><td>θ<sub>i</sub></td><td>temperature of the air in state i, °C</td></tr><tr><td>φ</td><td>relative humidity, %</td></tr><tr><td>ψ</td><td>overall heat transfer coefficient of a thermal bridge, W/m K</td></tr></table>

### Bold letters

<table><tr><td>A</td><td>matrix of coefficients</td></tr><tr><td>b</td><td>vector of constant terms</td></tr><tr><td>x</td><td>vector of unknowns</td></tr></table>

### Subscripts

<table><tr><td>BL</td><td>building</td></tr><tr><td>CC</td><td>cooling coil</td></tr><tr><td>HC</td><td>heating coil</td></tr><tr><td>TZ</td><td>thermal zone</td></tr><tr><td>da</td><td>dry air</td></tr><tr><td>i</td><td>air infiltration</td></tr><tr><td>1 ... 5</td><td>states in the psychrometric chart</td></tr><tr><td>l</td><td>latent</td></tr><tr><td>ls</td><td>least squares</td></tr><tr><td>o</td><td>outdoor</td></tr><tr><td>s</td><td>sensible or saturation</td></tr><tr><td>sp</td><td>set-point</td></tr><tr><td>t</td><td>total</td></tr><tr><td>vs</td><td>vapor at saturation</td></tr><tr><td>w</td><td>humidity ratio</td></tr><tr><td>θ</td><td>temperature</td></tr><tr><td>θ<sub>s</sub><sup>0</sup></td><td>initial guess of the temperature corresponding to the saturation point</td></tr></table>

### Superscripts

<table><tr><td>'</td><td>derivative of a function</td></tr><tr><td>0</td><td>initial value</td></tr></table>

# 1 Introduction

Since 1911, when Carrier published his seminal paper *Rational psychrometric formulae*, psychrometric chart has been used for calculations in air conditioning industry (Simha, 2012). Psychometrics chart is used usually in two approaches:

- manual procedures of psychrometric analysis for sizing (Carrier, 1972; TRANE, 2011; TRANE, 2002);

- computer calculations of psychrometric values for simulations (EQUA Simulation AB, 2013; TRNSYS 16, 2007; EnergyPlus, 2012; Satchell, 2016; Satchell, 2017).

There is a need for computational procedures for psychrometric analysis, similar to the manual methods used for sizing. This paper proposes to define the sizing problem of the cooling and dehumidification systems as a control problem and to use linear controllers and least-squares optimization to solve the sizing problem. The paper is accompanied by a software implementation for which the inputs can be changed in a large range. The code is published as open source (Ghiaus, 2021) and it is sharable, interactive, and reproducible (Project Jupyter et al., 2018).
## 2 Types of problems: direct and inverse

When a system is modelled in the view of deterministic classical physics, the causality principle states that any modification of a physical variable is caused by other variable(s) (see *Physical causality* in Figure 1). The dependent variables (or the effects) are the *outputs* and the independent variables (or the causes) are the *inputs* of the physical system. Likewise, a computational model has *inputs* (values provided to the computational algorithm) and *outputs* (numerical results of the algorithm).

<figure><img src="image_4.png" /><figcaption>Figure 1 Direct and inverse problems (Ghiaus, 2014). For the type of problem, in green is what is known and in red what is unknown.</figcaption></figure>

The *direct* problem is when the set of inputs of the computational model is a subset of the inputs of the physical model and the set of outputs of the computational model is a subset of the outputs of the physical model (see *Direct problem* in Figure 1). For a direct problem, the model (i.e., its structure and parameters) and the inputs are known; the problem is to find the outputs of the physical model. In the case of a direct problem, the flows of information (from inputs to outputs) are the same for the physical system and for the computational model. Direct problems are called *analysis* in the case of electrical, hydraulic, or thermal circuits. The term *simulation* is also often used. The advantage of direct problems is that, in the context of deterministic classical physics, they are well posed and have a unique solution.

The *inverse* problems are of two types: control and model identification (see *Inverse problem* in Figure 1). The *control* problem consists in finding an input of the physical system that produces a desired output when its model is known (see *Control problem* in Figure 1). For simple cases, the problem can be solved by inversing the direct problem. Usually, and especially for dynamic models, the inversion is not possible. Therefore, the problem is stated as a control loop for which a set of direct problems is solved (Ghiaus, 2013).

In *model identification*, the problem is to find the model of the physical system when its inputs and outputs are given. Usually, the problem is treated as *parameter identification* by imposing the structure of the model and solving an optimisation problem in which the arguments are the parameters and the cost function depends on the type of problem (see *Parameter identification or optimization* problem in Figure 1).
## 3 Elementary psychrometric processes: steady-state models

In air conditioning, the psychrometric properties (dry-bulb temperature, $θ$, humidity ratio, $w$, specific enthalpy, $h$, etc.) of a mass flow rate of dry air, $\dot{m}$, are controlled. Since the moist air is not oversaturated, the humidity ratio, $w$, is always smaller or equal to the humidity ratio at saturation, $w_s = f(θ)$.

### 3.1 Decoupling the sensible and the latent heat

When a mass flow rate of moist air changes from state 1, characterized by the enthalpy rate

$$ \dot{H}_1 = \dot{m}h_1 = \dot{m}[c\theta_1 + (l + c_w\theta_1)w_1] \quad (1) $$

where $\dot{H}_1$ is the enthalpy rate of the air in state 1 [W];

$\dot{m}$ -mass flow rate of dry air [kg/s];

$h_1$ -specific enthalpy of dry air in state 1 [J/kg<sub>da</sub>];

$c$ -specific heat capacity of the dry air [J/kg K];

$w_1$ -humidity ratio of moist air in state 1 [kg/kg<sub>da</sub>];

$l$ -specific latent heat for water vaporization [J/kg];

$c_w$ -specific heat capacity of the water vapours [J/kg K];

$\theta_1$ -temperature of the air in state 1 [°C],

to state 2, having the enthalpy rate

$$ \dot{H}_2 = \dot{m}h_2 = \dot{m}[c_{da}\theta_2 + (l + c_w\theta_2)w_2] \quad (2) $$

where $\dot{H}_2$ is the enthalpy rate of the air in state 2 [W];

$h_2$ -specific enthalpy of the dry air in state 2 [J/kg<sub>da</sub>];

$w_2$ -humidity ratio (or the moisture content) in state 2 [kg/kg<sub>da</sub>];

$\theta_2$ -temperature of the air in state 2 [°C],

the total variation of the enthalpy rate is

$$ \dot{H}_2 - \dot{H}_1 = \dot{m}(h_2 - h_1) \quad (3) $$

The total variation of enthalpy is the sum of sensible and latent changes:

$$ \dot{H}_2 - \dot{H}_1 = \Delta \dot{H}_s + \Delta \dot{H}_l \quad (4) $$

The change of sensible heat in the process is

$$ \Delta \dot{H}_s = \dot{m} \Delta h_s = \dot{m}[c(\theta_2 - \theta_1) + c_w(w_2\theta_2 - w_1\theta_1)] \quad (5) $$

where $\dot{m}c(\theta_2 - \theta_1)$ is the sensible heat of the dry air and $\dot{m}c_w(w_2\theta_2 - w_1\theta_1)$ is the sensible heat of the water vapor. The change in latent heat during the process is

$$ \Delta \dot{H}_l = \dot{m} \Delta h_l = \dot{m}l(w_2 - w_1) \quad (6) $$

In air conditioning systems, the error introduced in the evaluation of the enthalpy variation by neglecting the variation of sensible specific enthalpy of water vapor,
$$c_w(w_2\theta_2 - w_1\theta_1) \ll c(\theta_2 - \theta_1) \quad (7)$$

is less than 5 % in relative value or 2.5 kJ/kg in absolute value. Neglecting the sensible heat of the moisture is a key assumption that allows us to decouple the sensible and the latent heat and to obtain the model as a system of linear equations (Ghiaus, 2014).

<figure><img src="image_5.png" /><figcaption>Figure 2 Elementary psychrometric processes. The arrows show the inputs and the outputs of the model of the process.</figcaption></figure>
The specific air enthalpy is then

$$ \begin{cases} h_s = c\theta \\ h_l = lw \end{cases} \qquad (8) $$

where $h_s = c\theta$ is the sensible specific enthalpy and $h_l = lw$ is the latent specific enthalpy. With this assumption, the transformation of the moist air can be decomposed in sensible, $\Delta\dot{H}_s$, and latent, $\Delta\dot{H}_l$, exchange:

$$ \begin{cases} \Delta\dot{H}_s = \dot{m}c\theta_2 - \dot{m}c\theta_1 \\ \Delta\dot{H}_l = \dot{m}lw_2 - \dot{m}lw_1 \end{cases} \qquad (9) $$

Therefore, the models that describe the elementary air-conditioning processes may be based on decoupled sensible and latent heat balance. In the equations describing the elementary processes, the dependent variables (or the outputs) are written on the left-hand side of the equation and the independent variables (or the inputs) are on the right-hand side; in a block diagram, the inputs and the outputs are indicated by arrows (Figure 2). The specific heat capacity, c, and the specific latent heat, l, are considered constants. The mass flow of dry air, $\dot{m}$, is considered known. In these conditions, the elementary processes are represented on the psychrometric chart by the input and the output points of the process. The elementary processes are adiabatic mixing and processes with (sensible and/or latent) heat exchange (Figure 2).

## 3.2 Adiabatic mixing

The mixing of two airflows is without heat exchange with the environment. There are two cases (Figure 2):

1) mixing with given ratio $\alpha$: the two flows are $\alpha\dot{m}$ and $(1 - \alpha)\dot{m}$ (Figure 2, process 1);

2) mixing with a given mass flow rate $\dot{m}_o$: the two flows are $\dot{m}_o$, the given mass flow rate, generally the outdoor air, and $\dot{m} - \dot{m}_o$, generally the recycled or recirculated air (Figure 2, process 2).

## 3.3 Processes with heat exchange

The state of the moist air can be modified by sensible, $\dot{Q}_s$, and/or latent, $\dot{Q}_l$, heat exchange (Figure 2). If both sensible and latent heat are exchanged, the temperature and the humidity ratio of the air at the end of the process are modified (Figure 2, process 3). This is the case of a thermal zone. The simplest psychrometric processes consist in changing the humidity ratio by adding water vapor (Figure 2, process 4) or the temperature by heating or dry cooling (Figure 2, process 5). If the air is cooled at a temperature lower than its dew temperature, then the moist condenses and the process evolves on the saturation curve (Figure 2, process 6) (Eppelheimer, 2002; Carrier, 1972). In adiabatic humidification or dehumidification (Figure 2, process 7), the transformation is isenthalpic: the latent heat received by the air is equal to the sensible heat that it lost by the air, $\dot{Q}_l = -\dot{Q}_s = \dot{m}_v l$. The mass flow rate of evaporated water, $\dot{m}_v = \dot{m}(w_1 - w_2)$, is an output (i.e., a dependent variable) of the physical process.

## 3.4 Thermal loads of the building in steady-state models

Heating or cooling loads can be estimated by steady-state or dynamic calculations. Many steady-state methods use $U A \Delta\theta$ formulation (ASHRAE, 2017). The loads are decomposed in sensible and latent heat (Wang, 2001; Kreider, et al., 2010; ASHRAE, 2017; ASHRAE, 2017; Spitler, 2014). Generally, the total sensible load for a thermal zone $TZ$,
$$ \dot{Q}_{STZ} = \sum_{k} U_{k} A_{k} (\theta_{o,k} - \theta_{TZ}) + \sum_{k} \psi_{k} l_{k} (\theta_{o,k} - \theta_{TZ}) + \sum_{k} \dot{m}_{i,k} c (\theta_{o,k} - \theta_{TZ}) + \sum_{k} \dot{Q}_{SBL,k} \quad (10) $$

is the sum of:

1) $\sum_k U_k A_k (\theta_{o,k} - \theta_{TZ})$ heat flow from outdoor to the thermal zone through wall k of overall surface conductivity $U_k$ and surface area $A_k$ (opaque walls, glazing, roof, floor, doors, etc.) [W];

2) $\sum_k \psi_k l_k (\theta_{o,k} - \theta_{TZ})$ heat flow from outdoor to the thermal zone through the thermal bridge k of overall linear conductivity $\psi_k$ and length $l_k$ [W];

3) $\sum_k \dot{m}_{i,k} c (\theta_{o,k} - \theta_{TZ})$ advective heat transfer due to air infiltration of mass flow rate $\dot{m}_{i,k}$ [W];

4) $\sum_k \dot{Q}_{SBL,k}$ solar and all internal sensible heat gains of the building from occupants, electrical devices, cooking, etc. [W];

where $\theta_{o,k}$ is the temperature of the air outdoor of element k of the thermal zone [°C]; $\theta_{TZ}$ - temperature of the air in the thermal zone $TZ$ [°C].

The total latent load of the thermal zone $TZ$,

$$ \dot{Q}_{LTZ} = \sum_{k} \dot{m}_{i,k} l (w_{o,k} - w_{TZ}) + \sum_{k} \dot{Q}_{IBL,k} \quad (11) $$

is the sum of:

1) $\sum_k \dot{m}_{i,k} l (w_{o,k} - w_{TZ})$ latent heat flow from outdoor to the thermal zone advected by the mass flow rate $\dot{m}_{i,k}$ [W];

2) $\sum_k \dot{Q}_{IBL,k}$ all internal latent heat gains of the building from: occupants, electrical devices, cooking, baths, etc. [W],

where $w_{o,k}$ is the humidity ratio of the air outdoor of element k [kg/kg$_{da}$]; $w_{TZ}$ - humidity ratio of the air in the thermal zone $TZ$ [kg/kg$_{da}$].

The dynamic effects of the heat storage may be important for cooling, but, generally, they are negligible for heating. Storage of latent heat is negligible for most cases (Kreider, et al., 2010; TRANE, 2011). Cooling load temperature difference (CLTD), cooling load factor (CLF) or radiant time-series (RTS) methods, which take into account the dynamic and solar effects, may be put in the form of the system of equations (10). In the following, it is considered that the thermal loads of the building are:

$$ \begin{cases} \dot{Q}_{STZ} = (UA + \dot{m}_i c)(\theta_o - \theta_{TZ}) + \dot{Q}_{SBL} \\ \dot{Q}_{LTZ} = \dot{m}_i l(w_o - w_{TZ}) + \dot{Q}_{IBL} \end{cases} \quad (12) $$

where:

1) $(UA + \dot{m}_i c)(\theta_o - \theta_{TZ})$ is the sensible load through the envelope which takes into account the heat transmission, the solar radiation and the time lag due to wall inertia [W];

2) $\dot{Q}_{SBL}$ represents the solar and all internal sensible gains of the building from: occupants, electrical devices, cooking, etc. [W];

3) $\dot{m}_i l(w_o - w_{TZ})$ is the latent load through air infiltration [W];
4) $\dot{Q}_{IBL}$ represents all internal latent gains of the building from: occupants, electrical devices, cooking, etc. [W].

If the dynamic response of the building is considered, the term $(UA + \dot{m}_i c)$ would become variable in time; it could be expressed as a state-space model or as transfer functions (Ghiaus, 2013). Therefore, the first term of right-hand term of the 1st equation of systems (12) could not be expressed in function of the difference of temperatures $(\theta_o - \theta_{TZ})$. The dynamic formulation is out of scope of this paper.

# 4 Assembling the elementary processes to obtain a direct problem

A simulation model is obtained by assembling the elementary processes presented in Figure 2. The result is a system of linear and nonlinear equations that solves a direct problem. This system may be solved iteratively by linearizing the saturation curve (Ghiaus, 2014; Ghiaus, 2016).

## 4.1 Linearisation of saturation curve

The humidity ratio (mixing ratio or moisture content) is defined as the ratio of the water vapor mass, $m_w$, to the dry air mass, $m$:

$$w = \frac{m_w}{m} = \frac{M_w}{M_{da}} \frac{n_w}{n_{da}} \qquad (13)$$

where

$M_w$ = 18.015 286 kg/kmol and $M$ = 28.966 kg/kmol are the molar masses of water vapor and dry air, respectively,

$n_w$ = $m_w/M_w$ and $n$ = $m/M$ are the molar fraction of water vapor and of dry air, respectively.

In the moist air, the water vapor temperature is equal to that of the dry air. By using the ideal gas law, it results that the humidity content on the saturation curve, $w_s$, has a non-linear dependence on temperature (Gatley & Perez-Ganlindo, 2009):

$$w_s = f(\theta) = \frac{M_w}{M_{da}} \frac{p_{vs}}{p - p_{vs}} = 0.621 945 \frac{p_{vs}}{p - p_{vs}} \qquad (14)$$

where

$p$ is the barometric pressure [Pa];

$p_{vs}$ - saturation pressure of water vapor [Pa].

For the temperature range of -100 °C to 0 °C, the saturation pressure over ice is given by (Gatley & Perez-Ganlindo, 2009):

$$p_{vs} = \exp(C_1 / T + C_2 + C_3 T + C_4 T^2 + C_5 T^3 + C_6 T^4 + C_7 \ln T) \qquad (15)$$

where

$C_1 = -5.674 535 9 \cdot 10^3$, $C_2 = 6.392 524 7$, $C_3 = -9.677 843 0 \cdot 10^{-3}$,

$C_4 = 6.221 570 1 \cdot 10^{-7}$, $C_5 = 2.074 782 5 \cdot 10^{-9}$, $C_6 = -9.484 024 0 \cdot 10^{-13}$,

$C_7 = 4.163 501 9$.
For temperature range of 0 °C to 200 °C, the saturation pressure over liquid water is given by (Gatley & Perez-Ganlindo, 2009):

$$p_{vs} = \exp(C_8 / T + C_9 + C_{10}T + C_{11}T^2 + C_{12}T^3 + C_{13} \ln T) \quad (16)$$

where

$$C_8 = -5.800\,220\,6 \cdot 10^3, C_9 = 1.391\,499\,3 \cdot 10^0, C_{10} = -4.864\,023\,9 \cdot 10^{-2}, \\ C_{11} = 4.176\,476\,8 \cdot 10^{-5}, C_{12} = -1.445\,209\,3 \cdot 10^{-8}, C_{13} = 6.545\,967\,3 \cdot 10^0$$

In both equations (15) and (16), T [K] is the absolute temperature.

Alternatively, a simpler expression gives relatively accurate results within normal ranges of temperature (Tetens, 1930; Monteith & Unsworth, 2008; Murray, 1967; Murray, 1967; Singh, et al., 2002):

$$p_{vs} = 610.78 \exp \left( \frac{17.2694 \theta}{\theta + 237.3} \right) \quad (17)$$

where θ [°C] is the air temperature.

Equation (14), expressing the humidity content on the saturation curve, can be linearized for variations around θ_s^0:

$$w_s = w_s^0 + f_{\theta_s^0}'(\theta_s - \theta_s^0) \quad (18)$$

where $f_{\theta_s^0}' \equiv f'(\theta_s^0)$ is the value of the derivative of the function given by equations expressing the humidity content on the saturation curve. If equation (17) is used, then equation (18) becomes (Ghiaus, 2014):

$$f'_{\theta_s^0} = \frac{df(\theta)}{d\theta}\bigg|_{\theta=\theta_s^0} = \frac{d}{d\theta}\left(\frac{M_w}{M_{da}}\frac{p_{vs}}{p-p_{vs}}\right) \\ = \frac{M_w}{M_{da}}p \frac{2.502994 \cdot 10^6 e^{\frac{17.2694\theta_s^0}{\theta_s^0+237.3}}}{(\theta_s^0+237.3)^2 \left(p - 610.78 e^{\frac{17.2694\theta_s^0}{\theta_s^0+237.3}}\right)^2} \quad (19)$$

The use of approximate formula (17) for the computation of the derivative does not change the value of the humidity ratio for which the more precise equation (16) is used.

## 4.2 Assembling the elementary processes in a direct model

A direct (or simulation) model may be obtained by assembling the elementary processes given in Figure 2. The procedure has been already exemplified for a system with heating, adiabatic humidification and heat recovery processes (Ghiaus, 2014; Ghiaus, 2016).

Let's consider the widely used cooling and dehumidification system for a single zone (Stanke, 2000; Murphy, 2002; Murphy, 2004; Murphy, 2010). The cooling coil has an intrinsic bypass factor, depending on its geometrical characteristics, with typical values given in Table 1, or can be bypassed by a mixture of outdoor and return air, with typical values given in Table 2.
<table><caption>Table 1 Typical values of the bypass factor for finned coils (adapted from Carrier (1972), Table 61).</caption><thead><tr><th rowspan="3">Depth of coils<br>(rows)</th><th colspan="2">Without sprays</th><th colspan="2">With sprays</th></tr><tr><th>3 fins/cm</th><th>5 fins/cm</th><th>3 fins/cm</th><th>5 fins/cm</th></tr><tr><th colspan="4">Velocity (m/s)</th></tr><tr><th></th><th>1.5 – 3.5</th><th>1.5 – 3.5</th><th>1.5 – 3.0</th><th>1.5 – 3.0</th></tr></thead><tbody><tr><td>2</td><td>0.42 – 0.55</td><td>0.22 – 0.38</td><td></td><td></td></tr><tr><td>3</td><td>0.27 – 0.40</td><td>0.10 – 0.23</td><td></td><td></td></tr><tr><td>4</td><td>0.15 – 0.28</td><td>0.05 – 0.14</td><td>0.12 – 0.22</td><td>0.04 – 0.10</td></tr><tr><td>5</td><td>0.10 – 0.22</td><td>0.03 – 0.14</td><td>0.08 – 0.16</td><td>0.02 – 0.06</td></tr><tr><td>6</td><td>0.06 – 0.15</td><td>0.01 – 0.05</td><td>0.05 – 0.11</td><td>0.01 – 0.03</td></tr><tr><td>8</td><td>0.02 – 0.08</td><td>0.00 – 0.02</td><td>0.02 – 0.06</td><td>0.00 – 0.02</td></tr></tbody></table>

<table><caption>Table 2 Typical bypass factors for various applications (adapted from Carrier (1972), Table 62).</caption><thead><tr><th>Bypass</th><th>Type of application</th><th>Example</th></tr></thead><tbody><tr><td>0.30 – 0.50</td><td>Small total load or load with a low sensible heat factor.</td><td>Residence</td></tr><tr><td>0.20 – 0.30</td><td>Typical comfort application with small total load or low sensible heat factor</td><td>Residence</td></tr><tr><td>0.10 – 0.20</td><td>Typical comfort application</td><td>Small retail shop<br>Factory<br>Dept. store<br>Bank<br>Factory</td></tr><tr><td>0.05 – 0.10</td><td>High internal sensible loads or large amount of outdoor air</td><td>Dept. store<br>Restaurant</td></tr><tr><td>0 – 0.10</td><td>All outdoor air systems</td><td>Hospital</td></tr></tbody></table>

The model of the system, obtained by using the elementary processes from Figure 2, is represented in Figure 4. The system is composed of:

- A mixing box in which a mass flow rate of outdoor air, $\dot{m}_o$, controlled by the number of occupants or by the minimum ventilation rate (and in-field applications by $CO_2$ concentration), is mixed with recirculated indoor air. It is modelled by the elementary process 2: *Mixing with a given mass flow rate*, MR.

- A cooling coil which has a bypass factor, $\beta$, controlled by a damper $D$. They are modelled by the elementary processes 6: *cooling with dehumidification*, CC, and 2: *Mixing with given ratio*, MX.

- A reheating coil, modelled by the elementary process 5: *Heating or dry cooling*, HC.

- A fan, $F$, which imposes the mass flow rate of the supply air, $\dot{m}_s$.

- A single thermal zone. The transformation of the air in the thermal zone, $TZ$, is modelled by the elementary process 3: *Thermal zone*. The thermal loads of the zone are given by the block, $BL$, in function of the characteristics of the building, according to the system of equations (12) which has

inputs:

$UA$ global conductance of the building [J/K];

$\dot{m}_i$ air infiltration mass flow rate [kg/s];

$\dot{Q}_{SBL}$ solar and auxiliary sensible gains of the building;

$\dot{Q}_{IBL}$ auxiliary latent gains of the building;

outputs:

$\dot{Q}_{STZ}$ sensible load of the thermal zone;

$\dot{Q}_{ITZ}$ latent load of the thermal zone.
The inputs of the simulation model (indicated by inward arrows in Figure 4) are:

- temperature, $\theta_o$, and humidity ratio, $w_o$, of the outside air (symbolized by the letter $o$ above the double line),

- total load of the cooling coil, $\dot{Q}_{tcc}$ (symbolized by the red arrow $\dot{Q}_t$ related to the block CC),

- sensible load of the heating coil, $\dot{Q}_{scc}$ (symbolized by the red arrow $\dot{Q}_s$ related to the block HC),

- characteristics of the building: sensible, $\dot{Q}_{sBL}$, and latent, $\dot{Q}_{lBL}$ (symbolized by the red arrow $\dot{Q}_s$ and $\dot{Q}_l$ related to the block BL), infiltration air mass flow rate, $\dot{m}_i$, and overall conductance, UA (symbolized by the blue arrows $\dot{m}_i$ and UA related to the block BL).

The mass flow rates of supply air, $\dot{m}$, and of outdoor air, $\dot{m}_o$, and the by-pass factor, $\beta$, are parameters of the model. They are indicated below the double lines in Figures 4, 5, 7, 9, 11, 13 and 15.

In Figures 4, 5, 7, 9, 11, 13 and 15, the numbers above the double lines represent the temperatures, $\theta_k, k = 1...5$, and the humidity ratios, $w_k, k = 1...5$, of the state of the air marked by the points 1...5 in the psychrometric chart. These points (i.e., temperatures and humidity ratios), along with the heat flow rates:

$\dot{Q}_{scc}$ – sensible heat transferred by the cooling coil to the air flow [J];

$\dot{Q}_{lcc}$ – latent heat transferred by the cooling coil to the air flow [J];

$\dot{Q}_{stz}$ – sensible load of the thermal zone [J];

$\dot{Q}_{itz}$ – latent load of the thermal zone [J];

form the set of unknowns of the model. Note that:

- Point 0 (i.e., outdoor temperature, $\theta_o$, and humidity, $w_o$) on the psychrometric diagram is known (the outdoor temperature and humidity are independent variables).

- Sensible and latent loads of the thermal zone, $\dot{Q}_{stz}$ and $\dot{Q}_{itz}$, are unknowns (or outputs or dependent variables): they depend on the indoor conditions ($\theta_5, w_5$), which are unknowns.

In psychrometric analysis done by hand, the sensible and latent loads of the thermal zone, $\dot{Q}_{stz}$ and $\dot{Q}_{itz}$, are considered known. They allow us to calculate the sensible heat ratio

$$ SHR = \frac{\dot{Q}_{stz}}{\dot{Q}_{stz} + \dot{Q}_{itz}} \qquad (20) $$

that is used in graphical calculations (TRANE, 2011). However, according to equations (12), the loads depend on the indoor and outdoor conditions. Therefore, they will be considered unknowns.

For the direct problem shown in Figure 4, there are 14 unknowns (Table 3) and 14 equations obtained by assembling the elementary processes (Table 4). Equation (4) in Table 4 is the linearized form of equation (14) around the point ($\theta_s^0 \equiv \theta_2^0, w_s^0 \equiv w_2^0$) obtained by using equations (17) - (19).

A direct model is useful for verifying a solution; however, the sizing of a system is an inverse problem. In sizing, some inputs or parameters of the physical system become outputs of the computational algorithm, transforming the direct problem into an inverse one: control or parameter optimization, depending on the type of the requested values.
<table><caption>Figure 3 Legend for Figures 4 – 10.</caption><thead><tr><th colspan="2">Legend</th></tr></thead><tbody><tr><td colspan="2"><strong>Blocs</strong></td></tr><tr><td colspan="2">- Elementary processes</td></tr><tr><td>MR</td><td>mixing fix-rate in mixing box (2 eqs.)</td></tr><tr><td>CC</td><td>perfect cooling coil (4 eqs.)</td></tr><tr><td>MX</td><td>mixing process in cooling coil (2 eqs.)</td></tr><tr><td>HC</td><td>heating coil (2 eqs.)</td></tr><tr><td>TZ</td><td>thermal zone (2 eqs.)</td></tr><tr><td>BL</td><td>building (2 eqs.)</td></tr><tr><td colspan="2">- Proportional controllers</td></tr><tr><td>K<sub>w</sub></td><td>indoor air humidity ratio (1 eq.)</td></tr><tr><td>K<sub>θ</sub></td><td>indoor air temperature (1 eq.)</td></tr><tr><td colspan="2">- Other blocs</td></tr><tr><td>D</td><td>by-pass damper</td></tr><tr><td>F</td><td>fan</td></tr><tr><td>Ψ</td><td>psychrometric function w = Ψ(φ)</td></tr><tr><td>Is</td><td>least-squares controller</td></tr><tr><td colspan="2"><strong>Flows</strong></td></tr><tr><td></td><td>mass flow rate, temperature and moist ratio</td></tr><tr><td></td><td>heat flow rate</td></tr><tr><td></td><td>mass flow rate</td></tr><tr><td></td><td>information flow</td></tr></tbody></table><table><thead><tr><th colspan="2">Variables</th></tr></thead><tbody><tr><td colspan="2">- Unknowns</td></tr><tr><td>1, ..., 5</td><td>air states (θ<sub>k</sub>, w<sub>k</sub>)<sub>k=1..5</sub></td></tr><tr><td>Q̇<sub>tCC</sub></td><td>total heat of cooling coil</td></tr><tr><td>Q̇<sub>sCC</sub></td><td>sensible heat of cooling coil</td></tr><tr><td>Q̇<sub>ICC</sub></td><td>latent heat of cooling coil</td></tr><tr><td>Q̇<sub>sHC</sub></td><td>sensible heat of heating coil</td></tr><tr><td>Q̇<sub>STZ</sub></td><td>sensible heat of thermal zone</td></tr><tr><td>Q̇<sub>ITZ</sub></td><td>latent heat of thermal</td></tr><tr><td colspan="2">- Inputs</td></tr><tr><td>θ<sub>o</sub>, φ<sub>o</sub></td><td>outdoor air conditions</td></tr><tr><td>θ<sub>5,sp</sub>, φ<sub>5,sp</sub></td><td>indoor air conditions set-points</td></tr><tr><td>ṁ<sub>i</sub></td><td>mass flow rate of infiltration air</td></tr><tr><td>UA</td><td>overall conductance of building</td></tr><tr><td>Q̇<sub>sBL</sub></td><td>auxiliary sensible heat of building</td></tr><tr><td>Q̇<sub>IBL</sub></td><td>auxiliary latent heat of building</td></tr><tr><td colspan="2">- Parameters</td></tr><tr><td>ṁ</td><td>mass flow rate of supply air</td></tr><tr><td>ṁ<sub>o</sub></td><td>mass flow rate of outdoor air</td></tr><tr><td>β</td><td>by-pass factor of the cooling coil</td></tr><tr><td>K<sub>θ</sub>, K<sub>w</sub></td><td>proportional controller gains for temperature and humidity</td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>Figure 4 Cooling and dehumidification: direct (or simulation) problem obtained by assembling elementary psychrometric processes (from Figure 2) and by using the building thermal loads in steady state (§3.4). Legend in Figure 3.</figcaption></figure>
<table><caption>Table 3 List of unknowns for the direct model of cooling and dehumidification shown in Figure 4.</caption><thead><tr><th>0</th><th>1</th><th>2</th><th>3</th><th>4</th><th>5</th><th>6</th><th>7</th><th>8</th><th>9</th><th>10</th><th>11</th><th>12</th><th>13</th></tr></thead><tbody><tr><td>θ<sub>1</sub></td><td>w<sub>1</sub></td><td>θ<sub>2</sub></td><td>w<sub>2</sub></td><td>θ<sub>3</sub></td><td>w<sub>3</sub></td><td>θ<sub>4</sub></td><td>w<sub>4</sub></td><td>θ<sub>5</sub></td><td>w<sub>5</sub></td><td>Q&#x0307;<sub>SCC</sub></td><td>Q&#x0307;<sub>ICC</sub></td><td>Q&#x0307;<sub>STZ</sub></td><td>Q&#x0307;<sub>LTZ</sub></td></tr></tbody></table>

<table><caption>Table 4 System of equations for the direct model of cooling and dehumidification shown in Figure 4.</caption><thead><tr><th>N°</th><th>Equation</th><th>Symbol</th><th>Process</th></tr></thead><tbody><tr><td>0.</td><td>ṁcθ₁ - (ṁ - ḿ)θ₅ = ḿcθ₀</td><td>MR</td><td>Mixing with given flow rate</td></tr><tr><td>1.</td><td>ṁlw₁ - (ṁ - ḿ)lw₅ = ḿcw₀</td><td></td><td></td></tr><tr><td>2.</td><td>(1 - β)ṁcθ₁ - (1 - β)ṁcθ₂ + Q&#x0307;<sub>SCC</sub> = 0</td><td>CC</td><td>Cooling with dehumidification</td></tr><tr><td>3.</td><td>(1 - β)ṁlw₁ - (1 - β)ṁlw₂ + Q&#x0307;<sub>ICC</sub> = 0</td><td></td><td></td></tr><tr><td>4.</td><td>f'θ₀² θ₂ - w₂ = f'θ₀² θ₂' - w₂'</td><td></td><td></td></tr><tr><td>5.</td><td>Q&#x0307;<sub>SCC</sub> + Q&#x0307;<sub>ICC</sub> = Q&#x0307;<sub>ICC</sub></td><td></td><td></td></tr><tr><td>6.</td><td>βṁcθ₁ + (1 - β)ṁcθ₂ - ḿcθ₃ = 0</td><td>MX</td><td>Mixing with given ratio</td></tr><tr><td>7.</td><td>βṁlw₁ + (1 - β)ṁlw₂ - ḿlw₃ = 0</td><td></td><td></td></tr><tr><td>8.</td><td>ṁcθ₃ - ḿcθ₄ = -Q&#x0307;<sub>SHC</sub></td><td>HC</td><td>Heating</td></tr><tr><td>9.</td><td>ṁlw₃ - ḿlw₄ = 0</td><td></td><td></td></tr><tr><td>10.</td><td>ṁcθ₄ - ḿcθ₅ + Q&#x0307;<sub>STZ</sub> = 0</td><td>TZ</td><td>Thermal zone</td></tr><tr><td>11.</td><td>ṁlw₄ - ḿlw₅ + Q&#x0307;<sub>LTZ</sub> = 0</td><td></td><td></td></tr><tr><td>12.</td><td>(UA + mᵢc)θ₅ + Q&#x0307;<sub>STZ</sub> = (UA + mᵢc)θ₀ + Q&#x0307;<sub>BL</sub></td><td>BL</td><td>Building loads</td></tr><tr><td>13.</td><td>ṁᵢlw₅ + Q&#x0307;<sub>LTZ</sub> = ḿᵢlw₀ + Q&#x0307;<sub>BL</sub></td><td></td><td></td></tr></tbody></table>

# 5 Linear control problem: constant air volume (CAV) systems

A constant air volume (CAV) mixed-air system supplies a fixed mass flow rate of air, ḿ, usually to a single space or thermal zone (Stanke, 2000; TRANE, 2002). The aim of sizing is to find the thermal load of the cooling coil and its Apparatus Dew-Point ADP (Eppelheimer, 2002), and the thermal load of the reheating coil (if supply air tempering is used).

## 5.1 Control of indoor air temperature & humidity in CAV systems with reheating

Temperature and humidity can be controlled in the thermal zone by using a cooling coil in series with a heating coil (Stanke, 2000; TRANE, 2002). The aim of sizing is to find the total load of the cooling coil, Q&#x0307;<sub>ICC</sub>, and the sensible load of the reheating coil, Q&#x0307;<sub>SHC</sub>. This is a control problem (Figure 1) since the loads, which are inputs of the physical system, become unknowns and the temperature, θ₅, and the humidity, w₅, of the thermal zone, which are outputs of the physical system, are known (Figure 4). The proposed solution is to solve the control problem by introducing two control loops: one for temperature (controller K<sub>θ</sub> in Figure 5) and one for humidity (controller K<sub>w</sub> in Figure 5). In this approach, the temperature, θ₅, and the humidity, w₅, are still outputs of the system (i.e., unknowns) but their values will be very near to their set-points, θ₅,sp and w₅,sp. The errors between the values and their setpoints tend to zero as the gains of the controllers, K<sub>θ</sub> and K<sub>w</sub>, tend to infinity. Since usually the relative indoor humidity, φ₅,sp, is indicated, a function, Ψ, is used to transform it in humidity ratio, w₅,sp (Figure 5).

By assembling the elementary processes (Figure 2), the model of the system, having the unknowns given in Table 5, is obtained as a system of linear equations, given in Table 6.
<figure><img src="image_7.png" /><figcaption>Figure 5 CAV system controlling the indoor temperature, $θ_5$, and humidity, $φ_5$, by acting on the cooling coil [CC], $\dot{Q}_{tCC}$, and the heating coil [HC], $\dot{Q}_{sHC}$, respectively. The mass flow rate of supply air, $\dot{m}$, is given and constant (CAV system). The supply air temperature, $θ_4$, is not controlled. See CAV_wd and §5.1 in Ghiaus (2021). Legend in Figure 3.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 6 Psychrometric diagram for cooling and dehumidification of a single zone modelled by the block diagram given in Figure 5. Implementation in §5.1 in Ghiaus (2021).</figcaption></figure>

<table><caption>Table 5 List of unknowns for the model of cooling shown in Figure 5</caption><tr><td>0</td><td>1</td><td>2</td><td>3</td><td>4</td><td>5</td><td>6</td><td>7</td><td>8</td><td>9</td><td>10</td><td>11</td><td>12</td><td>13</td><td>14</td><td>15</td></tr><tr><td>θ<sub>1</sub></td><td>w<sub>1</sub></td><td>θ<sub>2</sub></td><td>w<sub>2</sub></td><td>θ<sub>3</sub></td><td>w<sub>3</sub></td><td>θ<sub>4</sub></td><td>w<sub>4</sub></td><td>θ<sub>5</sub></td><td>w<sub>5</sub></td><td>Q&#x0302;<sub>tCC</sub></td><td>Q&#x0302;<sub>sCC</sub></td><td>Q&#x0302;<sub>ICC</sub></td><td>Q&#x0302;<sub>SHC</sub></td><td>Q&#x0302;<sub>STZ</sub></td><td>Q&#x0302;<sub>LTZ</sub></td></tr></table>
<table><caption>Table 6 System of equations for the model of cooling and dehumidification with reheating and linear control of temperature and humidity ratio of the thermal zone (Figure 5).</caption><thead><tr><th>N°</th><th>Equation</th><th>Symbol</th><th>Process</th></tr></thead><tbody><tr><td>0.</td><td>ṁcθ₁ - (ṁ - ḿ)θ₅ = ḿcθ₀</td><td>MR</td><td>Mixing with given flow rate</td></tr><tr><td>1.</td><td>ṁlw₁ - (ṁ - ḿ)lw₅ = ḿcw₀</td><td></td><td></td></tr><tr><td>2.</td><td>(1 - β)ṁcθ₁ - (1 - β)ṁcθ₂ + Ḍ<sub>SCC</sub> = 0</td><td>CC</td><td>Cooling with dehumidification</td></tr><tr><td>3.</td><td>(1 - β)ṁlw₁ - (1 - β)ṁlw₂ + Ḍ<sub>ICC</sub> = 0</td><td></td><td></td></tr><tr><td>4.</td><td>f'<sub>θ₀</sub>θ₂ - w₂ = f'<sub>θ₀</sub>θ₂₀ - w₂₀</td><td></td><td></td></tr><tr><td>5.</td><td>- Ḍ<sub>TCC</sub> + Ḍ<sub>SCC</sub> + Ḍ<sub>ICC</sub> = 0</td><td></td><td></td></tr><tr><td>6.</td><td>βṁcθ₁ + (1 - β)ṁcθ₂ - ḿcθ₃ = 0</td><td>MX</td><td>Mixing with given ratio</td></tr><tr><td>7.</td><td>βṁlw₁ + (1 - β)ṁlw₂ - ḿlw₃ = 0</td><td></td><td></td></tr><tr><td>8.</td><td>ṁcθ₃ - ḿcθ₄ + Ḍ<sub>SHC</sub> = 0</td><td>HC</td><td>Heating</td></tr><tr><td>9.</td><td>ṁlw₃ - ḿlw₄ = 0</td><td></td><td></td></tr><tr><td>10.</td><td>ṁcθ₄ - ḿcθ₅ + Ḍ<sub>STZ</sub> = 0</td><td>TZ</td><td>Thermal zone</td></tr><tr><td>11.</td><td>ṁlw₄ - ḿlw₅ + Ḍ<sub>LTZ</sub> = 0</td><td></td><td></td></tr><tr><td>12.</td><td>(UA + mᵢc)θ₅ + Ḍ<sub>STZ</sub> = (UA + mᵢc)θ₀ + Ḍ<sub>SBL</sub></td><td>BL</td><td>Building loads</td></tr><tr><td>13.</td><td>ṁᵢlw₅ + Ḍ<sub>LTZ</sub> = ḿᵢlw₀ + Ḍ<sub>IBL</sub></td><td></td><td></td></tr><tr><td>14.</td><td>K<sub>θ</sub>θ₅ + Ḍ<sub>TCC</sub> = K<sub>θ</sub>θ<sub>5sp</sub></td><td>K<sub>θ</sub></td><td>Temperature controller</td></tr><tr><td>15.</td><td>K<sub>w</sub>w₅ + Ḍ<sub>SHC</sub> = K<sub>w</sub>w<sub>5sp</sub></td><td>K<sub>w</sub></td><td>Humidity controller</td></tr></tbody></table>

The system of equations from Table 6 is linear; it can be written as

$$Ax = b \tag{21}$$

where **x** is the vector of unknowns (Table 5), **A** is the matrix of coefficients, and **b** is the right-hand term obtained from the inputs of the system and physical constants (Figure 5). The system of equations (21) is solved by using `np laidg`. `solve` (Harris, et al., 2020)

Algorithm 1 Solve the nonlinear system of equations modelling the problem.

```python
solve_lin(θ²)
repeat
    solve Ax = b for x
    Δw₂ = |w² - w₂|
    θ² = θ₂
    until Δw₂ < ε
return x
```

The non-linear model (21) is solved by finding iteratively the solution of the system of equations from Table 6 by guessing the initial value of the Apparatus Dew Point, θ₂₀, which is on the saturation curve, w₂₀ = f(θ₂₀). Then, the value of θ₂₀ is updated to x₂ = θ₂. The procedure is repeated until the humidity ratio x₃ = w₂ is on the saturation curve with an acceptable error Δw₂ = |w₂₀ - w₂| < ε. The procedure, equivalent to Newton-Raphson method for solving non-linear equations, is given by Algorithm 1 (see CAV_wd in Ghiaus (2021) for Python implementation). Figure 6 gives the results on the psychrometric chart. The indoor conditions are θ₅ = 24 °C, w₅ = 9.30 g/kg, φ₅ = 50 %. All the values of the problem are given in the computational implementation (Ghiaus, 2021)
## 5.2 Control of indoor air temperature in CAV systems without reheating

Typically, a constant air volume (CAV) system with mixed recirculated and outdoor air uses a single cooling coil with modulated capacity in order to control the indoor temperature in the thermal zone (Stanke, 2000). The model, obtained by assembling the blocks from Figure 2, is shown in Figure 7.

<figure><img src="image_9.png" /><figcaption>Figure 7 CAV system controlling the indoor temperature, $θ_5$, by acting on the cooling coil [CC], $\dot{Q}_{tCC}$. The heating coil is not operating, $\dot{Q}_{SHC} = 0$. Indoor humidity, $w_5$, is not controlled. The mass flow rate of supply air, $\dot{m}$, is given and constant (CAV system). The supply air temperature, $θ_4 = θ_3$, is not controlled. See CAV_wd and §4.2 in Ghiaus (2021). Legend in Figure 3.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 8 Psychrometric diagram with typical processes for cooling and dehumidification of a single zone modelled by the block diagram given in Figure 7. Implementation in §5.2 in Ghiaus (2021).</figcaption></figure>
The system of equations and the solving algorithm are the same as for cooling and reheating discussed in § 5.1 (Table 6 and Algorithm 1). The difference is in setting the gain of the humidity controller to zero, $K_w = 0$ (Figure 7). Since the load of the heating coil is zero, $\dot{Q}_{SHC} = 0$, points 3 and 4 are superposed on the psychrometric diagram (Figure 8). The temperature of the thermal zone, $\theta_5$, is controlled, but the humidity ratio, $w_5$, varies freely. For the example given, the relative humidity changed from $\varphi_5 = 50.0 \%$, when it was controlled, to $\varphi_5 = 66.8 \%$ (Figure 8). All the values of the problem are given in the computational implementation (Ghiaus, 2021).

# 6 Non-linear parameter optimization problem

The cooling and dehumidification system shown in Figure 4 may be controlled by changing the mix-air bypass, $\beta$, or by modulating the mass flow rate of the supply air, $\dot{m}$ (in which case it becomes a variable air volume, VAV, system). Since the bypass factor, $\beta$, and the supply air mass flow rate, $\dot{m}$, are parameters in the model described by the set of equations given in Table 4, the inverse problem is a parameter optimization problem (Figure 1). It is treated by considering it a non-linear control problem solved by using an optimization algorithm.

## 6.1 Mix-air bypass control in constant air volume (CAV) systems without reheating

A simple and unexpensive solution to control the temperature and the humidity in the thermal zone is to use a bypass damper to extend the dehumidification range of a constant air volume, CAV, air handling unit, AHU (Murphy, 2002; TRANE, 2002). The bypass factor, $\beta$, can be found by solving the optimization problem:

$$ \text{find } \beta \text{ so that } \varepsilon(\beta) = |v_{sp} - v| < \varepsilon_{adm} \qquad (22) $$

where $\varepsilon_{adm}$ is the admissible error between the variable $v$ and its set-point $v_{sp}$. The optimization can be achieved by least-squares (Virtanen, et al., 2020).

### 6.1.1 Mix-air bypass control of indoor air temperature & humidity in CAV systems

The block diagram of a constant air volume (CAV) system with mix-air bypass control of temperature and humidity is given in Figure 9. A non-linear controller $l_s$ changes the bypass factor, $\beta$, by acting on the damper $D$.

Let's consider that the value that we want to control is the indoor air humidity, $w_5$. By defining the error function as

$$ \varepsilon(\beta) = |w_{5,sp} - w_5| \qquad (23) $$

the problem is then to find the optimal value

$$ \varepsilon_{min} = \min_{\beta}(\varepsilon) \text{ subject to } 0 \le \beta \le 1 \qquad (24) $$

by using least squares optimization (Virtanen, et al., 2020).

If $\varepsilon_{min} \le \varepsilon$, then, the optimal solution
$$ \beta = \arg \min_{\beta} (\varepsilon) \text{ subject to } 0 \le \beta \le 1 \qquad (25) $$

is the value of the bypass factor. If the minimum value of the error is not practically zero, i.e., $\varepsilon_{min} \ge \varepsilon$, then there is no solution to the parameter optimization problem, i.e., the humidity of the thermal zone cannot be controlled by using a mix-air damper. Algorithm 2 synthesis the procedure (see VBP_wd (Ghiaus, 2021) for Python implementation).

<figure><img src="image_11.png" /><figcaption>Figure 9 CAV with mix-air bypass controlling the indoor temperature, $\theta_5$, by acting on the cooling coil [CC], $\dot{Q}_{tCC}$, and the indoor humidity, $\varphi_5$, by acting on the by-pass damper D, $\beta$. The heating coil is not operating, $\dot{Q}_{SHC} = 0$. The mass flow rate of the supply air, $\dot{m}$, is given and constant (CAV system). The supply air temperature, $\theta_4 = \theta_3$, is not controlled. See VBP_wd and §6.1.1 in Ghiaus (2021). Legend in Figure 3.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 10 Psychrometric diagram for cooling and dehumidification by a CAV system with mix-air bypass modelled by the block diagram given in Figure 9. Implementation in §6.1.1 in Ghiaus (2021).</figcaption></figure>
Algorithm 2 Find the bypass factor, $\beta$, by least squares optimization.

function $\beta_{ls}$ ($v$, $v_{sp}$):
    function $\varepsilon(\beta)$:
        solve\_lin($\theta_s^0$) (see Algorithm 1)
        return $\varepsilon(\beta) = |v_{sp} - v|$ 
    find $\varepsilon_{min} = \min_{\beta}(\varepsilon)$ subject to $0 \le \beta \le 1$
    if $\varepsilon_{min} < \varepsilon_{adm}$:
        return $\beta = \arg\min_{\beta}(\varepsilon)$ subject to $0 \le \beta \le 1$
    else print an error message

Figure 10 shows the results. As in Figure 8, the points 3 and 4 are superposed since the heating coil HC is not operating. Both temperature and humidity of the indoor air are controlled: $\theta_5 = 24$ °C, $w_5 = 9.29$ g/kg, $\phi_5 = 50\%$. They have practically the same values as in § 5.1. The difference is that the bypass factor is $\beta = 0.16$ in § 5.1 and $\beta = 0.0072$ in § 6.1.1. All values of the problem are given in the computational implementation (Ghiaus, 2021).

### 6.1.2 Control of indoor and supply air temperatures in CAV systems

Thermal comfort may be affected by the supply air temperature. Therefore, it would be interesting to control the supply air temperature by acting on the bypass factor. However, the problem has no solution. The sensible heat balance of the thermal zone, equation 10 from Table 4 or Table 6,

$$ \dot{m}c(\theta_4 - \theta_5) + \dot{Q}_{STZ} = 0 \qquad (26) $$

has no solution when all the variables are imposed:

$\dot{m}$ - the mass flow rate imposed because the system is with constant air volume (CAV);

$\theta_4$ - the supply air temperature;

$\theta_5$ - the indoor air temperature;

$\dot{Q}_{STZ}$ - the sensible load of the thermal zone.

The value of $\varepsilon_{min}$ from relation (22) will not be almost zero (save special cases). Therefore, the temperatures of the supply air and of the thermal zone cannot be controlled by the cooling coil, $\dot{Q}_{tcc}$, and the bypass factor, $\beta$.

## 6.2 Mass-flow control in variable air volume (VAV) systems

In variable air volume (VAV) systems, the mass flow rate of supply air, $\dot{m}$, is modulated (TRANE, 2002; Murphy, 2010). Since the mass flow rate, $\dot{m}$, is a parameter of the model (Figure 4, Table 4), modulating it is a parameter optimization problem (Figure 1) that can be solved by least squares (Algorithm 3).

Algorithm 3 is similar to Algorithm 2: a value, $v$, and its set-point, $v_{sp}$, are specified. A least squares optimization problem is solved (Virtanen, et al., 2020). If the minimum of the error function is almost zero, $\min_{\dot{m}}|v_{sp} - v| < \varepsilon_{adm}$, then there is a solution, i.e., a supply air mass
flow rate $\dot{m}$ which can control the variable $v$. The controlled variable can be the indoor air humidity of the supply air temperature.

<caption>Algorithm 3 Find the mass flow rate of the supply air, $\dot{m}$, by least squares optimization.</caption>

```matlab
function \dot{m}_{ls}(v, v_{sp}):
    \epsilon(\dot{m}):
        solve_lin(\theta_{2}^{0})          (see
    Algorithm 1)
    return \epsilon(\dot{m}) = |v_{sp} - v|
find \epsilon_{min} = \min_{\dot{m}}(\epsilon) \text{ subject to } 0 \le \dot{m} \le \dot{m}_{max}
if \epsilon_{min} < \epsilon_{adm}:
    \dot{m} = \arg\min_{\dot{m}}(\epsilon) \text{ subject to } 0 \le \dot{m} \le \dot{m}_{max}
else
    print an error message
end
```

### 6.2.1 Control of indoor air temperature and humidity in VAV systems without reheating

If the humidity ratio of the thermal zone is controlled, the variable $v$ in the Algorithm 3 is the relative humidity of the indoor air, $v \equiv \phi_5$. Figure 11 shows the block diagram of the model obtained with elements from Figure 2. The non-linear least-squares controller $ls$ acts on the mass flow rate of supply air, $\dot{m}$, in order to control the humidity of the indoor air, $\theta_5$. The processes are presented on the psychrometric chart on Figure 12. Since the heating coil $HC$ is not operating, points 3 and 4 are superposed. The system cannot control the indoor humidity at 50%; the optimization gives $\theta_5 = 24$ °C, $w_5 = 11.83$ g/kg (i.e., $\phi_5 = 63.4$ %). The humidity in the thermal zone may be maintained at 50% if the outdoor humidity is lower (e.g., $\phi_0 = 50$ %). All values of the problem are given in the computational implementation (Ghiaus, 2021).

<figure><img src="image_13.png" /><figcaption>Figure 11 VAV system controlling the indoor temperature, $\theta_5$, by acting on the cooling coil [CC], $\dot{Q}_{tCC}$, and the indoor air humidity, $\phi_5$, by acting on the supply air mass flow rate, $\dot{m}$. The heating coil is not operating, $\dot{Q}_{SHC} = 0$. See VAV_wd and §6.2.1 in Ghiaus (2021). Legend in Figure 3.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>Figure 12 Psychrometric diagram for cooling and dehumidification by a VAV system controlling the thermal zone temperature and humidity modelled by the block diagram given in Figure 11. Implementation in §6.2.1 in Ghiaus (2021).</figcaption></figure>

### 6.2.2 Control of supply and indoor air temperatures in VAV systems without reheating

If the humidity ratio of the thermal zone is controlled, the variable in the Algorithm 3 is the supply air temperature, $v \equiv \theta_4$. Figure 13 shows the model in which the non-linear least squares controller *ls* acts on the mass flow rate of supply air, *ṁ*, in order to control the supply air temperature $\theta_4$. Figure 14 shows the processes on the psychrometric diagram on which points 3 and 4 are identical because the heating coil is not operating. All values of the problem are given in the computational implementation (Ghiaus, 2021).
<figure><img src="image_15.png" /><figcaption>Figure 13 VAV system controlling the indoor temperature, $θ_5$, by acting on the cooling coil [CC], $\dot{Q}_{tCC}$, and supply air temperature, $θ_4$, by acting on the supply air mass flow rate, $\dot{m}$. The heating coil is not operating, $\dot{Q}_{SHC} = 0$. See VAV_wd and §6.2.2 in Ghiaus (2021). Legend in Figure 3.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Figure 14 Psychrometric diagram for cooling and dehumidification by a VAV controlling the temperatures of the supply air and the thermal zone modelled by the block diagram given in Figure 13. Implementation in §6.2.2 in Ghiaus (2021).</figcaption></figure>

### 6.2.3 Control of indoor air temperature and humidity and of supply air temperature in VAV systems with reheating

Typically, a variable air volume (VAV) system has a central air handling unit (AHU) that supplies air to several zones, each served by a VAV terminal unit. While a constant air volume
(CAV) supplies a fixed mass flow rate of air at variable temperature, the VAV terminals modulate the air mass flow rate which has a constant temperature between 7 °C and 13 °C (Murphy, 2010). Supply temperature and indoor temperature and humidity can be controlled by “tempering” (i.e., reheating) the supply air at the VAV terminal (TRANE, 2002).

<figure><img src="image_17.png" /><figcaption>Figure 15 VAV system controlling the indoor temperature, $θ_5$, by acting on the cooling coil [CC], $Q_{tCC}$, indoor humidity, $φ_5$, by acting on the heating coil [HC], $Q_{sHC}$, and supply air temperature, $θ_4$, by acting on the supply air mass flow rate, $\dot{m}$. See VAV_wd and §6.2.3 in Ghiaus (2021). Legend in Figure 3.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Figure 16 Psychrometric diagram for VAV controlling the supply air temperature and the temperature and humidity in the thermal zone modelled by the block diagram given in Figure 15. Implementation in §6.2.3 in Ghiaus (2021).</figcaption></figure>
Figure 15 shows the model. There are three control loops acting on cooling coil, $\dot{Q}_{tCC}$, heating coil, $\dot{Q}_{SHC}$, and fan, $\dot{m}$, in order to control the indoor temperature, $\theta_5$, indoor humidity, $w_5$ or $\varphi_5$, and supply air temperature, $\theta_4$. Figure 16 shows the results on a psychrometric chart. Both temperature and humidity are controlled: $\theta_5 = 24$ °C, $w_5 = 9.29$ g/kg, $\varphi_5 = 50\%$. All values of the problem are given in the computational implementation (Ghiaus, 2021).

# 7 Conclusions

This paper proposes a classification of the modelling problems into direct and inverse by using the information flow in the physical processes and the computational algorithms. Direct problems are when the physical and computational causality are the same; they are used to verify the solutions. Inverse problems are when the physical and the computational causality are different. In his paper, they are used for system sizing. The inverse problems are classified in control and parameter optimization. This paper proposes to solve the control problem by linear control (with practically infinite gain) and the parameter optimization problem by formulating it as a least squares control problem.

The classification is then used to solve the problem of sizing the cooling and dehumidification systems. It is proposed that a set of elementary psychrometric processes to be used to model the physical elements of the air handling units and terminals. By assembling the elementary blocks, the models of the air handling unit (in this case cooling coil, heating coil, mixing box) and of the thermal zone can be obtained. Most of the equations describing the processes (i.e., the direct model) are linear. After linearizing the equation of the saturation curve, the system of equations is solved iteratively; this is equivalent to using Newton-Raphson method.

The linear control problem is exemplified for a constant air volume (CAV) system with or without reheating coil. The control system has as setpoints the temperature and the humidity of the thermal zone and acts on the cooling and the heating coils. Since the outputs of the control system are inputs of the physical system, it is a linear control problem. The model consists of a system of equations which is nonlinear due to the vapor saturation curve. It is solved iteratively as a linear model of the direct problem.

The nonlinear control problem consists in changing the bypass factor of the cooling coil or the mass flow rate of the supply air. In the framework proposed in this paper, this is a parameter optimization problem, i.e., find a parameter (the bypass factor or the mass flow rate of supply air) that minimizes the error between the controlled variable (indoor air temperature and humidity, and supply air temperature) and its setpoint. The optimization algorithm used was least squares.

The advantage of the proposed modelling method is that it can be algorithmized. The solutions to a set of variable conditions can be easily obtained in the design phase. The limitation of the method is the use of steady-state models. This problem may be alleviated by corrections made to the steady-state models.
Annex: Psychrometric relations

The relations used in the psychrometric calculations (module psychro.py (Ghiaus, 2021)) are
given below.

Molar mass of dry air (Gatley & Perez-Ganlindo, 2009):

$$
M_{da} = 28.966 \text{ kg/kmol} \tag{27}
$$

Molar mass of water vapor (Gatley & Perez-Ganlindo, 2009):

$$
M_v = 18.015\ 286\ \text{kg/kmol} \tag{28}
$$

Universal gas constant:

$$
R = N_A k = 8\,314.462\,618\,153\,24\,\text{J/(K kmol)} \quad (29)
$$

where N_A = 6 022.140 76 · 10⁻²³ kmol⁻¹, Avogadro constant, and k = 1.380 649 · 10⁻²³ J/K, Boltzmann constant, have values defined in SI (BIPM, 2019).

Atmospheric pressure as a function of altitude -500 ... 10 000 m (Gatley & Perez-Ganlindo, 2009):

$$
p = 101.325 (1 - 2.225\,77 \cdot 10^{-5}Z)^{5.255\,9} \text{ Pa} \quad (30)
$$

Saturation vapor pressure as a function of temperature: see equation (16) (Gatley & Perez-Ganlindo, 2009):

$$
p_{vs} = \exp\left(\frac{C_8}{T} + C_9 + C_{10}T + C_{11}T^2 + C_{12}T^3 + C_{13}\ln T\right) \quad (31)
$$

Humidity ratio as a function of temperature and relative humidity:

$$
w = \frac{M_v}{M_{da}} \frac{\varphi p_{vs}}{p - \varphi p_{vs}} \text{ kg/kg}_{\text{da}} \quad (32)
$$

Derivative of the humidity ratio (32) at saturation for a given temperature uses the convenient
form of saturation vapor pressure over water:

$$
p_{vs} = C \exp\left(\frac{a \theta}{\theta + b}\right) \qquad (33)
$$

where a = 17.2693882, b = 273.16 – 35.86 = 237.3, and C = 610.78 (Murray, 1967). By substituting (33) in (32), the derivative of humidity ratio, w_s, at saturation, i.e., φ = 1 in (32), as a function of temperature, θ, is:
$$ \frac{dw_s}{d\theta} = \frac{dw_s}{dp_{vs}} \frac{dp_{vs}}{d\theta} = \frac{M_v}{M_{da}} \frac{p}{(p - p_{vs})^2} \frac{a b C \exp \left( \frac{a\theta}{\theta + b} \right)}{(\theta + b)^2} \\ = \frac{M_v}{M_{da}} \frac{a b}{(\theta + b)^2} \frac{p p_{vs}}{(p - p_{vs})^2} \qquad (34) $$

Specific volume as a function of temperature and humidity ratio for a given altitude (default 0 m) (Gatley & Perez-Ganlindo, 2009)is:

$$ v = \frac{R}{M_v} \left( \frac{M_v}{M_{da}} + w \right) (\theta + 273.15) \qquad (35) $$

Relative humidity as a function of temperature and humidity ratio for a given altitude (default 0 m) is:

$$ \varphi = \frac{p}{p_{vs}} \frac{w}{\frac{M_v}{M_{da}} + w} \qquad (36) $$

## 8 References

ASHRAE, 2017. 17. Residential Cooling and Heating Load Calculations. In: *Handbook Fundamentals*. Atlanta, GA: ASHRAE.

ASHRAE, 2017. 18. Nonresidential Cooling and Heating Load Calculations. In: *Handbook - Fundamentals*. Atlanta, GA: ASHRAE.

BIPM, 2019. *The International system of Units*. 9 ed. Paris: BIMP.

Carrier, 1972. *Carrier System Design Manual Part 1: Load Estimating*. ASIN: B00NN6RM54: s.n.

EnergyPlus, 2012. *EnergyPlus Engineering Reference*. University of Illinois and Lawrence Berkeley National Laboratory: s.n.

Eppelheimer, D., 2002. Cooling-Coil Heat Transfer. *Trane Engineers Newsletter*, Volume 31(1), pp. 1-4.

EQUA Simulation AB, 2013. *IDA Indoor Climate and Energy - User Manual*, Solna, Sweden: s.n.

Gatley, D. & Perez-Ganlindo, J., 2009. Chapter 1 Psychrometrics. In: *ASHRAE Handbook - Fundamentals (SI)*. Atlanta: ASHRAE, pp. 10-26.

Ghiaus, C., 2013. Causality issue in the heat balance method for calculating the design heating and cooling load. *Energy*, Volume 50, pp. 292-301.

Ghiaus, C., 2014. Linear algebra solution to psychometric analysis of air-conditioning systems. *Energy*, Volume 74, pp. 555-566.

Ghiaus, C., 2016. Analyse psychrometrique des systèmes de climatisation. *Revue générale du froid & du conditionnement d'air*, pp. 38-42.

Ghiaus, C., 2021. *PsychroAn_cool: Psychrometric analysis of cooling systems as a control problem*, doi: 10.5281/zenodo.5236450: Zenodo.

Harris, C. et al., 2020. Array programming with NumPy. *Nature*, Volume 585, pp. 357-362.

Kreider, J., Curtiss, P. & Rabl, A., 2010. 7. Heating and Cooling Loads. In: *Heating and cooling of buildings - design for efficiency*. Boca Raton, FL: CRC Press - Taylod & Francis Group, pp. 273-366.
Monteith, J. & Unsworth, M., 2008. *Principles of Environmental Physics*. 3 ed. Amsterdam: AD.

Murphy, J., 2002. Dehumidification performance of HVAC systems. *ASHRAE Journal*.

Murphy, J., 2004. Better Part-Load Dehumidification. *Trane Engineers Newsletter*, Volume 33(2), pp. 1-6.

Murphy, J., 2010. Using VAV to limit humidity at part load. *ASHRAE Journal*, pp. 18-22.

Murray, F., 1967. On the Computation of Saturation Vapor Pressure. *Journal of Applied Meteorology*, Volume 6, pp. 203-204.

Project Jupyter et al., 2018. *Binder 2.0 - Reproducible, Interactive, Sharable Environments for Science at Scale*. s.l., Proceedings of the 17th Python in Science Conference, pp. 113-120.

Satchell, B., 2016. *TAS Systems Controller Guide v1.1*, s.l.: EDSL.

Satchell, B., 2017. *TAS Systems Component Guide*, s.l.: EDSL.

Simha, R., 2012. Willis H. Carrier Father of Air Conditioning. *Resonance*, pp. 117-138.

Singh, A., Singh, H., Singh, S. & Sawhney, R., 2002. Numerical calculation of phychrometric properties on a calculator. *Building and Environment*, Volume 37, pp. 415-419.

Spitler, J., 2014. *Load Calculation Applicatins Manual*. 2 ed. Atlanta, GA: ASHRAE.

Stanke, D., 2000. Dehumidify with constant-volume systems. *TRANE Engineers Newsletter*, 29(4), pp. 4-8.

Tetens, O., 1930. Über einige meteorologische Begriffe. *Zeitschrift fur Geophysik*, Volume 6, pp. 297-309.

TRANE, 2002. *Dehumidification in HVAC Systems*. s.l.:TRANE.

TRANE, 2011. *Cooling and Heating Load Estimation*, s.l.: TRANE.

TRNSYS 16, 2007. *Mathematical Reference*. Madison, WI, USA: University of Wisconsin.

Virtanen, P. et al., 2020. SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python. *Nature Methods*, Volume 17, pp. 261-272.

Wang, S., 2001. 6. Load Calculations. In: *Handbook of air conditioning and refrigeration*. New York: McGraw-Hill.