Evaluation of Thermal Comfort Inside an Office Equipped with a Fan Coil HVAC System: A CFD Approach

Giovanni Semprini; Aminhossein Jahanbin; Beatrice Pulvirenti; Paolo Guidorzi

1. Introduction

Most industrialized countries, in particular the European Union countries, are in shifting phase from a fossil-based development to a sustainable development. Approximately 40% of the energy consumption in the European Union is attributed to buildings across the private and public sectors. The European Parliament Directive 2010/31 () states that new buildings should consume very low amount of fossil-based energy (Nearly Zero Energy Buildings) by 2020. The Directive declares that all new buildings built after 31 December 2018 will have to be buildings with high energy performance and their required amount of energy should be mostly covered by renewable energy sources. Hence, energy efficiency has become a common target for all buildings, in particular historical buildings (). Since buildings with historic value are less energy-efficient than modern ones, energy demand and indoor comfort expectations are critical challenges for them.

Indeed, indoor thermal comfort in a moderate environment is considered as the human’s perception of well-being and satisfaction. Thermal comfort is defined by ASHRAE () as “that condition of mind that express satisfaction with thermal environment and is assessed by subjective evaluation”. Overall thermal comfort and the assessment of indoor environment quality depend not only on physical parameters but also on psychological and physiological response of the human’s body to the environment. Thermal comfort can be evaluated by either classic approach introduced by Fanger () or adaptive approach (). The classic approach interprets the global and local thermal comfort by means of two indexes: Predicted Mean Vote (PMV) and Predicted Percentage of Dissatisfied (PPD), adopted also in ISO 7730 (). However, acquisition of the required data to predict thermal comfort is not an easy task; it is time demanding and requires specific instruments to calculate the relevant parameters.

An interesting alternative technique to predict thermal comfort could be simulations implemented by Computational Fluid Dynamic (CFD) codes, allowing to analyse thermo-hydraulic characteristics within environments. Several studies were conducted in order to determine the thermal comfort with CFD simulations in various environments such as office room (), educational building (), theatre (), stadium (), museum () and test room (). For instance, Sevilgen and Kilic (), by means of finite volume method, performed a 3D steady-state analysis of a room heated by two-panel radiators. A virtual sitting manikin with real dimensions and physiological shape was added to the model of the room, and two different pairs of heat transfer coefficients of outer wall and window were considered. The obtained results showed that better-insulated wall and window could enhance thermal comfort while reducing the energy consumption.

A 3D numerical model was employed by Nada et al. () to investigate the airflow, temperature distribution and thermal comfort in a high-rise ceiling theatre air conditioned with UnderFloor Air Distribution (UFAD). The results showed that the numerical model can accurately predict the airflow and temperature distribution in the high-rise conditioned space and that the low air velocities of a UFAD provide the appropriate level of environmental comfort. It was recommended to use higher numbers of diffusers for achieving better thermal comfort condition.

Buratti et al. () employed a CFD tool together with the experimental data to predict indoor condition and thermal comfort inside a room, by setting only external weather conditions as input parameters. The survey and data acquisition for evaluation of thermal comfort were carried out in a classroom at the department of engineering, university of Perugia. It was concluded that CFD code could be a very useful and cost-effective tool to support the experimental campaign and to evaluate the thermal local sensation.

Cheong et al. () presented a simulation-aided approach for analysing and optimizing lighting and thermal conditions in building interior. A finite volume CFD code was adopted for prediction of air flow and ambient temperature. A methodology on a real-world case study on an office building in Singapore was demonstrated, improving lighting and thermal conditions greatly through passive daylight fixtures and the reconfiguration of climate-control air supply vents within the office space. The obtained results showed a significant reduction in energy consumption and enhancement of indoor environmental quality.

Thermal comfort in office rooms has been widely investigated both to increase people satisfaction and also to analyse solutions to reduce energy consumption (, , ). In many situations, the internal layout of HAVC terminal systems and workstations inside a single room are not properly investigated in the design phase, causing gradient in temperature and air velocity with different comfort conditions for occupants. Furthermore, in existing historical buildings with low thermal performance, asymmetry of radiant temperature occurs in winter period and too high air temperature is set inside the rooms (). CFD analysis is an important tool in design phase for determination of the internal layout, the operation mode of terminals (inlet temperature, air flow), and the position of sensors for room control system (, ): the latter plays an important role to maintain internal temperature set-point.

Accuracy of the CFD results depends mainly on the adopted boundary conditions and convection model. In application of CFD analysis to real buildings and large rooms, where many conditions affect the thermo-hydraulic behaviour of the air flow, some assumptions must be made and a comparison and/or calibration in specified conditions is recommended (, ).

In the present paper, local thermal comfort in a university office room equipped with fan coil heating system has been evaluated in typical winter conditions. The case study is an office room located inside a typical historical building, built with thick masonry walls and large single-glass windows with high thermal transmittance. Architectural constraints, often combined to economic difficulties in those public buildings, prevent the use of retrofit technical solutions (wall insulation, new performing windows, local control system, etc.), causing high energy consumption due to higher indoor air temperatures (sometimes set at 23–24°C) while at the same time maintains local discomfort at working positions.

The analysis has been performed by means of CFD simulations and the numerical model is validated by experimental results. A new methodology for sensors localization based on acoustic signals sampling and processing is presented for the experimental measurement and model validation. The spatial distribution of the air temperature and velocity are presented and values of the operative temperature for specific positions inside the room are analysed. The obtained results allow to evaluate thermal comfort condition at working places, to evidence discomfort zones that are related to the characteristics of the building and heating system under study, and to identify the best comfort areas inside the office for optimum internal layout.

2. Numerical Method

The schematic of the computational domain under consideration is illustrated in Figure 1. It is composed of a 4.73 m × 5.03 m office with 4.24 m height, equipped with a fan-coil air conditioning system, with an exit surface of the air of 0.08 m × 0.65 m = 0.052 m². The mean value of air velocity at the exit surface of the fan coil was measured equal to 0.8 m/s. The value of air temperature exiting from the fan coil, depending on the centralized water heating system, was measured in the range 33–34°C.

Figure 1

Schematic of the office under study.

The office domain consists of two external walls with windows, and two internal walls. 5 specific positions inside the room was considered to evaluate the condition of thermal comfort: 4 working desks and the central point of the room. As boundary condition, the internal walls were assumed to be isothermal with measured mean value T = 21°C, while for the external walls boundary condition of third kind was considered, with convection coefficient h = 25 W/m² k. Surface walls of the floor and ceiling were considered as adiabatic, as they are adjoining to the heated room.

The thermophysical properties of moist air were calculated and employed in numerical simulations (). The values of thermal properties of dry air at atmospheric pressure () and of the saturation properties of water () were adopted at reference temperature. The following values of the density, dynamic viscosity, specific heat capacity at constant pressure and thermal conductivity of moist air were obtained for T = 20°C: ρ_m = 1.153 kg m^–3, μ_m = 1.82 × 10^–5 kg m^–1 s^–1, c_pm = 1020 J kg^–1 K^–1, and k_m = 0.0253 W m^–1 K^–1. It should be noted that ρ_m is the reference value for the density of moist air. The fluid density term was considered to be constant except in the buoyancy force term, where it was assumed as a linear function of the temperature.

In analogy with recent literature on CFD simulations of the rooms equipped with a HVAC system (, , ), a 3D steady-state thermo-hydraulic model was considered assuming the air flow inside the room turbulent. The governing equations, i.e. the mass, momentum and energy balance equations, can be written as:

(1)

∇ ⋅ v → = 0

M1 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[\nabla \cdot \vec v = 0\] \end{document}

(2)

ρ 0 ∂ v → ∂ t + v → ⋅ ∇ v → = − ∇ p + ρ 0 g z + ρ 0 g β T − T 0 ∇ z + ∇ ⋅ τ ˜ eff

M2 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\rho _0}\left[ {\frac{{\partial \vec v}}{{\partial t}} + \left({\vec v \cdot \nabla } \right)\vec v} \right] = - \nabla \left({p + {\rho _0}g{\kern 1pt} z} \right) + {\rho _0}g\beta \left({T - {T_0}} \right)\nabla z + \nabla \cdot {\tilde \tau _{eff}} \] \end{document}

(3)

∂ ∂ t ρ e + ∇ ⋅ v → ρ e + p = ∇ ⋅ k + k t ∇ T + τ ˜ eff ⋅ v →

M3 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \frac{\partial }{{\partial t}}\left({\rho e} \right) + \nabla \cdot \left({\vec v\left({\rho e + p} \right)} \right) = \nabla \cdot \left[ {\left({k + {k_t}} \right)\nabla T + {{\tilde \tau }_{eff}} \cdot \vec v} \right] \] \end{document}

In Equations 1–3, $v →$ M16 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \vec v \] \end{document} is the velocity vector, p is pressure, β is the thermal expansion coefficient, $τ ˜ eff$ M15 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\tilde \tau _{eff}} \] \end{document} is the effective stress tensor, E = h – p/ρ + v²/2, h is specific enthalpy, v is the magnitude of $v →$ M14 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \vec v \] \end{document} , k is the thermal conductivity, k_t is the turbulent thermal conductivity and subscript 0 refers to the reference values. The effective stress tensor is given by:

(4)

τ ˜ eff = τ ˜ + τ ˜ t

M4 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\tilde \tau _{eff}} = \tilde \tau + {\tilde \tau _t} \] \end{document}

Among common turbulence models, the κ-ɛ model is the most widely used and validated model. It applies to a wide range of flows, which occur in many technical applications. The RNG (Renormalization Group) κ-ɛ model was selected for this study, due to the higher precision and stability rendering for internal flows (, ). The RNG κ-ɛ model is quite similar to the classic κ-ɛ model but includes some refinements such as: enhanced accuracy for swirling flows, an analytically derived differential formula for effective-viscosity that accounts for low-Reynolds numbers effects, and an additional term in ɛ to improve accuracy. The RNG κ-ɛ model’s governing equations are (, ):

(5)

∂ ρ к ∂ t + ∇⋅ ρ к v → = ∇ ⋅ α к μ eff ∇ к + 2 μ t E ij E ij − ρ ε

M5 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[\frac{{\partial \left( {\rho \kappa\,} \right)}}{{\partial t}} + \nabla \cdot \left( {\rho \ \kappa\,\vec v} \right) = \nabla \cdot \left( {{\alpha _\kappa} \mu _{eff} \nabla \kappa} \right) + 2{\mu _t}{E_{ij}}{E_{ij}} - \rho \varepsilon\] \end{document}

(6)

∂ ρ ε ∂ t +∇⋅ ρ ε v → = ∇ ⋅ α ε μ eff ∇ ε + C 1 ε * ε κ 2 μ t E ij E ij − C 2 ε ρ ε 2 κ

M6 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \frac{{\partial \left({\rho \,\varepsilon } \right)}}{{\partial t}}\; + \nabla \cdot \left({\rho \,\varepsilon \,\vec v} \right)\; = \nabla \cdot \left({{\alpha _\varepsilon }{\mu _{eff}}\nabla \varepsilon } \right) + C_{1\varepsilon }^*\frac{\varepsilon }{\kappa }2{\mu _t}{E_{ij}}{E_{ij}} - {C_{2\varepsilon }}\rho \frac{{{\varepsilon ^2}}}{\kappa } \] \end{document}

In Equations 5 and 6,

(7)

μ eff = μ + μ t; μ t =ρ C μ κ 2 ε

M7 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\mu _{eff}} = \mu + {\mu _t}\;;\quad {\mu _t} = \rho \,{C_\mu }\frac{{{\kappa ^2}}}{\varepsilon } \] \end{document}

where C_μ = 0.0845. In addition, α_κ = α_ɛ = 1.39; C_1ɛ = 1.42; C_2ɛ = 1.68 and $C 1 ε * = C 1 ε − η 1 − η / η 0 1 + γ η 3$ M12 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[C_{1\varepsilon }^* = {C_{1\varepsilon }} - \frac{{\eta \left( {1 - \eta /{\eta _0}} \right)}}{{1 + \gamma {\eta ^3}}}\] \end{document} where $η= 2 E ij E ij κ ε$ M13 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \eta = \sqrt {2{E_{ij}}{E_{ij}}} \frac{\kappa}{\varepsilon} \] \end{document} and η0 = 4.377; γ = 0.012.

Radiation heat transfer was considered by applying the surface-to-surface radiation model (). In this model, the air inside the room is considered as a transparent medium and all surfaces are assumed to be diffusive and grey. The following values were considered for the emissivity of surfaces: 0.88, 0.89, 0.87, and 0.9 for fan coil box, window (glass), door (glass), and walls, respectively ().

The effects of solar loading are needed in many HVAC applications, where the temperature, humidity, and velocity fields around occupants are desired. In the present study, effect of solar irradiation on windows has been taken into account by means of the solar ray tracing algorithm, for different working hours namely at 9:00, 12:00, 15:00 and 18:00. The ray tracing approach is a highly efficient and practical means of applying solar loads as heat sources in the energy equations. The resulting heat flux is coupled to the calculation as a source term in the energy equation. The heat sources are added directly to computational cells bordering each face and are assigned to adjacent cells.

The values of direct solar irradiation and diffuse solar irradiation have been adopted from the data extracted from the Italia climatic data (), on February 6^th for the Bologna, Italy, with latitude 44.4938°N and longitude 11.3387°E. Figure 2 illustrates the variation of outdoor temperature, direct and diffuse solar irradiation during a typical winter day (February 6^th) for Bologna, Italy. The values of direct solar irradiation are equal to 94.5, 325.7, 229 and 0 W/m² at 9:00, 12:00, 15:00 and 18:00 o’clock, respectively. Orientation of computational mesh with respect to the sun direction vector has been calculated by employing Solar Calculator program in the software. The values of spectral fraction and sunshine factor has been set equal to 0.5 and 0.65, respectively.

Figure 2

Outdoor temperature and solar irradiation (direct and diffuse) for a typical winter day of Bologna ().

In order to calculate the operative temperature, T_opr, the equation recommended by ASHRAE () was employed, namely:

(8)

T opr = ψ T air + 1 − ψ T rad

M8 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {T_{opr}} = \psi {T_{air}} + \left({1 - \psi } \right){T_{rad}} \] \end{document}

where T_air is the air temperature, T_rad is the mean radiant temperature, and ψ is a weight factor depending on the air velocity; this factor is usually considered between 0.4 and 0.5. In this study, it was assumed equal to 0.5. The operative temperature is calculated inside the office for 5 different positions (at 4 working desks + central point of the room), as illustrated in Figure 1, for different working hours.

The mean radiant temperature of the room surfaces with respect to a virtual person inside the room is given by:

(9)

T rad = ∑ 1 n F j T s j

M9 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[{T_{rad}} = \sum\limits_1^n {{F_{j\;}}{T_s}_j}\] \end{document}

where T_sj is the temperature of the j-th surface and F_j is the angular factor between the person and the j-th surface. In order to calculate view factors involved in Equation (9), a small cubic surface with size of 10 cm was considered as a local sensor at corresponding positions.

In general, thermal comfort is classified according to the type of environment: indoor, outdoor or semi-outdoor. In terms of indoor thermal comfort, the discussion centres mainly on two distinct approaches: classic steady-state model based on a heat balance model of human’s body, developed by Fanger (), and adaptive model which is based on field studies in naturally ventilated buildings (). In the present study, Fanger’s classic model has been employed to evaluate indoor thermal comfort condition. In fact, Fanger’s model aims to predict mean thermal sensation of a group of people and their respective percentage of dissatisfaction with the thermal environment, expressed through the indexes Predicted Mean Vote (PMV) and Predicted Percentage Dissatisfied (PPD) (). According to ASHRAE () and ISO 7730 (), there are six primary factors that must be addressed when defining conditions for thermal comfort: metabolic rate, clothing insulation, air temperature, radiant temperature, air speed and humidity. Thermal comfort can be assessed through two primary metrics, i.e. PMV and PPD. The PMV predicts the mean response of building occupants based on the ASHRAE thermal sensation scale ranging from the value of –3 to +3, corresponding to cold and hot relative to optimum comfort condition. PMV is a function of four environmental variables, namely air temperature T_air, mean radiant temperature T_rad, air velocity v, and air humidity p_air, and also of the clothing insulation I_cl and the person’s metabolic rate M. Thus, it can be written:

(10)

PMV = f T air, T rad, v, p air, I cl, M

M10 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\rm{PMV}} = f\left({{T_{air}},\,{T_{rad}},\,v,\,{p_{air}},\,{I_{cl}},\,M} \right) \] \end{document}

On the other hand, the PPD predicts the proportion of occupants dissatisfied with the specified thermal conditions. Based on experimental studies in which participants voted on their thermal sensation, an empirical relationship between PMV and PPD was developed as:

(11)

PPD = 100 − 95 exp − 0.3353 PM V 4 − 0.2179 PM V 2

M11 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\rm{PPD}} = 100 - 95{\rm exp}\left({ - 0.3353PM{V^4} - 0.2179PM{V^2}} \right) \] \end{document}

In the present work, values of PMV and PPD have been considered as main parameters to evaluate indoor thermal comfort at defined positions, namely working desks and central point of the office room. The adopted values for the metabolic rate, clothing insulation and relative humidity are 1.2 met, 1.0 clo and 60%, respectively. In addition, the local thermal discomfort caused by a vertical air temperature difference between head and ankles of the occupant and by the radiant temperature asymmetry has been assessed for the positions under study. According to ASHRAE (), allowable vertical air temperature difference between head and ankles is ΔT_air < 3°C and allowable radiant temperature asymmetry due to cool wall is ΔT_rad < 10°C.

A finite volume CFD code implemented in ANSYS Fluent was employed for a 3D fluid flow and heat transfer analysis. The governing differential equations were solved by means of the control volume approach, which converts the governing equations to a set of algebraic discretized equations. Second-order discretization method was used for the pressure, momentum and energy terms, and the well-known SIMPLE (Semi-Implicit Method for Pressure-Linked Equation) was used to solve the pressure-velocity coupling. A 3D TGrid unstructured mesh composed of tetrahedral elements was generated by employing the Gambit software. The computational domain contains 2 799 219 tetrahedral cells. The mesh independence of the results was assured by performing preliminary computations with different meshes.

3. Experimental measurements and the model validation

The determination of the temperature in the 3D measurement grid was carried out with an innovative patented system () which measures the distance of the temperature sensor (a thermocouple) from some known points, similar to the principle of the GPS. Analytically, the problem is addressed by solving a system of 4 equations with 4 unknowns through a specific trilateration algorithm, which include the four distances of the target from four known points and the three spatial coordinates associated to four known points (). The 4 distances are measured by means of acoustic waves: 4 different MLS (Maximum Length Sequence) signals orthogonal to each other are emitted simultaneously by 4 small loudspeakers placed on a grid, all in a known position (), as demonstrated in Figures 3a and 4. The microphone, fixed very close to the thermocouple as shown in Figure 3c, samples the 4 MLS signals, which, despite being emitted at the same time, can be discriminated and through the Fast-Hadamard Transform algorithm are converted into 4 distinct impulse responses, as shown in Figure 5b. From the peak of these impulse responses, the flight time of the sound from the respective loudspeaker to the single microphone is calculated. Knowing the speed of sound, corrected according to the average air temperature for more accuracy, the 4 distances from the microphone to the 4 known points (loudspeakers) are obtained and with the algorithm described above, the position in the 3D space of the microphone and thermocouple is found (). The idea of the position detection system was born from the need to correctly position microphones during measurements on road noise barriers (). The temperature from the thermocouple is detected by an Arduino system, as shown in Figure 3b, interfaced directly with the software for measuring the position on the computer, illustrated in Figure 5a. The temperature measurement sub-system allows a resolution of 0.25 degrees.

Figure 3

Measurement system: a) loudspeaker grid and microphone/thermocouple; b) detail of the Arduino system for reading the temperature from the thermocouple, interfaced with the measurement software on the computer; c) detail of the microphone and thermocouple; d) example of measurement in a high position.

Figure 4

Panoramic view of the room being measured; the measured points are visible in 3D inside the box.

Figure 5

a) Measurement system: software and hardware components (microphone preamplifier, loudspeakers power amplifiers, multi-channel soundcard); b) example measurement of the four impulse responses sampled by the microphone.

The CFD model was validated by means of experimental temperature measurement, as explained above. In order to evaluate the accuracy of numerical model, 32 positions randomly spread inside the room were considered (Figure 4). For each position, the value of temperature measured by experimental test was compared with that obtained by CFD results. The average reported error of temperature evaluation was 3.6% with reference to experimental data, with a range of variation from 0.1 to 9.5%. This can be considered a satisfying model validation.

In addition, values of the air velocity were measured and compared to the numerical results. Since the air velocity inside the room is too low, in particular for the positions far from the fan coil, 5 different points near the fan coil were selected and measured with an air velocity meter (anemometer), with accuracy ±0.015 m/s. The considered points were selected at distance between 0.6–1.1 meter from the fan coil with height between 1–2 meter from the floor. Table 1 compares values of the air velocity yielded by experiments and those obtained by numerical simulations. The mean square deviation of numerical results from experimental data is 0.038 m/s. The comparison shows a good agreement between experimental and numerical results.

Table 1

A comparison for the values of air velocity: experimental vs. numerical.

Position	distance from fan coil (m)	distance from floor (m)	v_exp (m/s)	v_num (m/s)

a	0.6	1.0	0.05	0.06
b	0.6	1.5	0.11	0.12
c	0.6	2.0	0.52	0.43
d	1.1	1.5	0.01	0.02
e	1.1	2.0	0.02	0.04

4. Results and discussion

The distribution of the air temperature and velocity, obtained by numerical simulations, allows to understand the thermal behaviour of the office room during a typical winter day whilst the heating system works in a stationary condition at fixed velocity and inlet temperature of the fan coil, as usual in centralised plants without local control systems.

Figures 6 and 7 illustrate air temperature and velocity distribution in a vertical section of the room for two different periods of a day, namely at morning and afternoon. The figures demonstrate clearly the effect of mixed convection inside the office room. In Figure 6, the vertical air stratification is evident from floor level to ceiling which is limited to an acceptable range of temperatures thanks to fan coil outlet air direction at 45°, while worst conditions would have occurred in case of vertical outlet air. The figure shows that the temperature distribution in the lower half of the vertical plane is more uniform in the afternoon, due to solar irradiation effect.

Figure 6

Temperature distribution in a vertical plane in middle of the room, at 9:00 (left) and at 15:00 (right).

Figure 7

Velocity pathline in a vertical plane in middle of the room, at 9:00 (left) and at 15:00 (right).

Figure 7 demonstrates that in vicinity of the outlet vent of fan coil, the forced condition due to air velocity is prevalent and streamlines follow the outlet direction; in the upper half of the room, the air velocity reduces to some extent and the chimney effect is more evident together with the induction effect which creates vortices and increases the air velocity on top of the room. This effect is more evident in the morning when the difference between the mean room temperature and the fan coil outlet air temperature is higher than afternoon.

The horizontal temperature map illustrated in Figure 8, at the height 1.2 m above the floor, shows a general uniform distribution inside the room, varying from 20 to 21°C at 9:00 and 22 to 23°C at 15:00, except for the zone in front of the fan coil with higher temperature. It can be noted that lateral positions of the fan coil (see also the values of position 3 in Table 2) assume the air temperature 0.5–1.2°C lower than central zone, due to reduced mixing effect.

Table 2

The values of the air, mean radiant and operative temperature, PMV, and PPD.

Position	9:00					Position	12:00

	T_air (°C)	T_rad (°C)	T_opr (°C)	PMV	PPD		T_air (°C)	T_rad (°C)	T_opr (°C)	PMV	PPD

1	21.16	13.89	17.52	–0.68	14.83	1	22.08	16.47	19.28	–0.33	7.24
2	21.13	16.74	18.93	–0.42	8.76	2	22.06	17.78	19.92	–0.21	5.89
3	20.65	15.10	17.88	–0.64	13.52	3	21.85	16.43	19.14	–0.36	7.77
4	21.29	16.52	18.90	–0.42	8.76	4	21.83	17.58	19.70	–0.26	6.36
5	21.87	16.89	19.38	–0.32	7.07	5	22.48	17.88	20.18	–0.14	5.43
Position	15:00					Position	18:00

	T_air (°C)	T_rad (°C)	T_opr (°C)	PMV	PPD		T_air (°C)	T_rad (°C)	T_opr (°C)	PMV	PPD

1	22.27	17.37	19.82	–0.25	6.28	1	21.70	15.87	18.79	–0.43	8.90
2	22.52	18.46	20.49	–0.08	5.15	2	21.66	17.42	19.54	–0.29	6.78
3	22.62	17.39	20.00	–0.30	6.85	3	21.21	16.07	18.64	–0.48	9.74
4	22.42	18.29	20.36	–0.11	5.27	4	21.84	17.20	19.52	–0.29	6.75
5	23.19	18.58	20.88	0.01	5.01	5	22.36	17.58	19.97	–0.19	5.73

Figure 8

Temperature distribution in a horizontal plane at height 1.2 m, at 9:00 (left), and at 15:00 (right).

The lack of temperature control (thermostat) in the office room causes that fan coil operates always at same conditions, and at fixed water temperature and fan speed (unless manual variation of fan speed); for this reason, the internal temperature can undergo large variation during the day, as shown in the figures.

Even if the air temperature is within usual set point for offices buildings, actual comfort conditions are strongly influenced by radiation from internal surfaces. Table 2 reports the results of operative temperature (T_opr), Predicted Mean Value (PMV) and related Predicted Percentage Dissatisfied (PPD) in 4 different hours of a day, for 5 considered positions. The table shows that due to large surface and low performance of windows, the radiant temperature reduces substantially T_opr, especially in the morning when the external temperature is lower. The results show that the comfort indexes at positions 1 and 3 are undesirable since these positions are located near the windows and have higher difference between the values of T_air and T_rad.

Table 2 shows that T_opr, in the morning (at 9:00), ranges from 17.5°C in position 1 to 19.3°C in position 5, with corresponding PMV of –0.68 and –0.32, respectively. It can be observed that in whole office room a general condition of “slightly cool” is prevalent. At other times of the day, the mean radiant temperature increases to some extent, compared with the morning. In addition, due to the solar radiation effect, the air temperature increases about 1°C, with maximum values at 15:00, namely 22.5–23.0°C. For the whole day, PMV takes negative values for all working positions. However, the comfort criterion is satisfied for positions 2, 4 and 5. It can be concluded that a minimum distance of 2 m from windows is required in order to lower the radiating effect of the windows and consequently to have an acceptable thermal comfort condition.

In order to guarantee the indoor thermal comfort condition, it is necessary to analyse local discomfort parameters such as the effects of the vertical air temperature difference and radiant asymmetry inside the room. The corresponding results of the vertical air temperature difference and radiant asymmetry are reported in Table 3.

Table 3

The values of vertical air temperature difference and radiant temperature asymmetry.

Position	9:00		12:00		15:00		18:00

	ΔT_air (°C)	ΔT_rad,max (°C)	ΔT_air (°C)	ΔT_rad,max (°C)	ΔT_air (°C)	ΔT_rad,max (°C)	ΔT_air (°C)	ΔT_rad,max (°C)

1	2.44	10.16	2.13	7.58	2.26	6.34	2.41	8.86
2	2.61	7.02	2.51	5.29	2.52	4.33	2.59	6.01
3	2.39	6.22	2.08	5.03	2.18	4.25	2.28	5.49
4	2.52	8.26	2.20	6.47	1.80	5.42	2.81	7.28
5	2.73	5.83	2.37	4.59	2.82	3.76	2.81	5.02

Vertical air temperature difference, ΔT_air, evaluated between head and ankles of a sitting person (0.1 to 1.1 m) is always below 3°C, as recommended by ISO 7730 in order to have PPD <5%. This difference takes into account only the effect of air flow due to the fan coil operation: in many real cases windows infiltrations and opening doors near colder zones causes higher ΔT_air.

Maximum horizontal radiant temperature asymmetry, ΔT_{rad, max}, caused by cool windows, gives the highest values in position 1 for the whole day. In particular, this difference exceeds 10°C giving a PPD of almost 15% in morning hours (see value at 9:00). As reported in ISO 7730, PPD increases rapidly with ΔT_rad leading to higher discomfort condition in winter season.

5. Conclusions

A 3D Computational Fluid Dynamic (CFD) code has been employed to evaluate thermo-hydraulic comportment inside a historic university office, built with thick masonry walls and large glass windows using fan coil as an air conditioning system. By means of finite volume method, turbulent air flow and local heat transfer characteristics of the room have been analysed for a typical winter day and thermal comfort indexes for working places inside the office have been determined. An innovative experimental set-up based on acoustic sources has been used to validate the numerical model and to obtain required initial and boundary conditions.

The obtained results show that even when the air temperature is quite uniform inside the room, the presence of low insulated envelope such as large glass windows leads to lower radiant temperature and as a consequence, to have higher possibility of local discomfort. Furthermore, the lack of local control system (thermostat) causes large variation in comfort conditions during working day which demands higher energy consumption. The results indicate that internal layout of working desks should be rearranged in the way that enough space between the working desks and the windows be taken into account in order to have a local comfort condition.

It can be concluded that CFD is an important tool allowing evaluation of local comfort sensation and optimization of thermal comfort in a specific workstation. Moreover, it could be used to evaluate innovative approaches to improve comfort indexes. The possibility of thermal comfort optimization together with local control of heating system leads to significant energy saving results. A potential direction for the future research is to assess the effect of using double-glazed windows on indoor comfort indexes for the similar historic buildings having large windows to wall surface ratio and a centralized heating system.

Future Cities and Environment

Technical Articles

Evaluation of Thermal Comfort Inside an Office Equipped with a Fan Coil HVAC System: A CFD Approach

Abstract

1. Introduction

2. Numerical Method

3. Experimental measurements and the model validation

4. Results and discussion

5. Conclusions

Competing Interests

References