Application of Combined Micro- and Macro-Scale Models to Investigate Heat and Mass Transfer through Textile Structures with Additional Ventilation

: In this study, computational models of heat and mass exchange through textile structures with additional ventilation at the micro- and macro-scale were investigated. The finite element analysis of advanced textile materials provides a better understanding of their heat and mass transfer properties, which influence thermal comfort. The developed computational models can predict air permeability (AP), thermal resistance (R ct ), and heat transfer (h) coefficients at the micro-scale. Moreover, the mesh size was taken into consideration and validated with experimental data presented in the literature. In addition, computational models were extended to micro- and macro-scale forced ventilation models. Macro-scale finite element models require input parameters such as an effective heat transfer coefficient that are usually obtained experimentally. In this research, the heat transfer coefficients (h microlayer = 25.603 W/(K·m 2 ), h total = 8.9646 W/(K·m 2 )) were obtained numerically from the micro-scale model and were applied to a macro-scale model. The proposed methodology and developed models facilitate the determination of average temperature and temperature distributions through different through-thickness positions along the axis Oz. The simulations were carried out using Comsol Multiphysics and Matlab software.


Introduction
Clothing and textile products are an integral part of our everyday life. The main functions of textile materials are to provide physiological and psychological comfort and to protect the body from hazardous environments. In the human-clothing system, this protection has several functions such as maintaining the right thermal environment for the body and preventing the body from being injured by abrasion, radiation, wind, electricity, chemical, and microbiological toxic substances [1]. There are many different types of comfort, including thermal comfort, sensorial comfort, garment fit, and psychological comfort (e.g., esthetics) [1,2]. Thermal comfort is an essential component for maintaining thermal balance. Some of the most significant thermal comfort properties of textiles are water vapor resistance (breathability), thermal resistance, air permeability, water resistance (under hydrostatic pressure), and water repellency [3]. However, the theoretical investigations of thermal comfort properties of advanced textiles are complicated due to complex internal structures, coupled heat and moisture transport through fabrics, and other physical processes, which occur in different space and time scales [4]. systems (HPCS) incorporated with ventilation fans and phase change materials (PCMs) have been studied by Xu et al. [22]. Siddiqui and Sun [23,24] developed finite element models that allow the prediction of the effective thermal conductivity of plain weft knitted fabric, MicroPCM-coated cotton, wool, and Nomex composite fabrics. The numerical simulations were carried out using the ABAQUS software. An analytical clothing model that takes into account the effect of body movement on heat transfer was proposed by Joshi et al. [25,26]. Renard and Puszkarz [27] developed three-dimensional computational models of textile assemblies with complex morphology used in firefighters' multilayer protective clothing. Computer simulations of eleven textiles included in the four tested assemblies were created using SOLIDWORKS software. The heat transfer simulations were designed to determine the average temperature of the top surface [27]. Codau et al. [28] investigated heat transfer textile porous media at a micro-and macro-scale. In proposed models, the heat transfer due to convection and radiation was neglected. The models were built using Comsol Multiphysics software.
Taking into account the analysis presented above, various computational tools may be used to investigate textile structures. Computational modeling is a powerful and costeffective tool that enables numerical simulation of the physical behavior of textile structures during the early design stage, particularly when the structure is too complicated for analytical analysis. This saves money on raw materials, energy, labor, and the time it takes to produce them [27]. However, parameters such as heat and mass transfer coefficients, as well as evaporative resistances, are usually obtained experimentally and utilized as input for numerical simulations on heat and mass transport through textile structures [29,30].
The aim of this research is to develop and combine micro-and macro-scale models to investigate heat and mass transfer through textile structures with additional ventilation. This study consists of four parts. The first part is concerned with determining an optimal mesh. It is a crucial step, in order, to minimize simulation time while maintaining accuracy. A micro-scale computational model of the air permeability (AP) coefficient was created for validation purposes. The second part is dedicated to developing computational models of thermal resistance (Rct) and heat transfer (h) coefficient at the micro-scale in accordance with the LST EN 31902/ISO 11092:2014(en) standard [31]. The third part focuses on the forced ventilation model at a micro-scale. It should be noted that the proposed computational model requires a heat transfer coefficient obtained from the second part. Furthermore, additional analysis is required to identify boundary conditions on inlet and outlet surfaces due to air pressure and air flux density decrease by passing through numerous cells. The fourth part corresponds to a macro-scale finite element model that predicts temperature in the ventilation layer. In this macro-scale model, the heat transfer coefficient is required as input and it was determined from the forced ventilation model at a micro-scale. The numerical simulations were performed using the finite element method program Comsol Multiphysics and Matlab software.

The Representative Volume Element
In this study, the representative volume element (RVE) is used to determine the physical and geometrical properties of the 3D textile structure at a micro-scale. Usually, 3D textile structures are made up of two separate outer fabric layers (top and bottom layers), that are connected by the spacer yarns or a layer. Spacer fabrics such as 3D textile structures can be produced by using warp knitting, weft knitting, weaving, non-woven, or braiding technologies [32]. It should be noted that, in the case of a plain weave, the smallest periodic RVE consists of at least two warp and two weft yarns [33]. In this research, RVE was based on an experimental investigation conducted by Zupin et al. [34] in which construction parameters and the air permeability coefficient of one-layer woven structures were measured. The experiments on air permeability coefficient (AP) were done according to the ISO 9237:1995 (E) standard [35].
A representative volume element (RVE) of the 3D textile structures that were used in numerical simulations is depicted in Figure 1. Moreover, Model x_a, Model x_b illustrate close-to-real geometry of one-layer and two-layer textile structures. Model x_c and Model x_d depict simplified geometry models of one-layer and two-layer textile structures. In addition, Model x_e represents a simplified geometry 3D textile structure. The computational domain Ω consists of the fluid domain Ω1 and the textile domain Ω2. A more detailed explanation of RVE is based on the research carried out by the authors of this article [36].  Table 1 describes the RVE structural parameters that were used in this study. In this investigation, the main focus is on Model x_e which indicates a simplified geometry of 3D textile structure. Furthermore, Model 1_c was used to find an optimal number of mesh elements.

Mathematical Model and Assumptions
In this subsection, the forced ventilation model at a micro-scale is presented. The proposed model was based on an experimental cooler system (see Figure 2) [4]. It was assumed that the flow is laminar incompressible flow. In addition, the heat is transmitted due to conduction, convection, and the flow is non-isothermal flow. The mathematical model consists of a set of Navier-Stokes (see Equations (1) and (2)) [38] and energy (see Equation (3)) equations [38].
where u-velocity vector of the fluid in m/s, p-pressure in Pa, -dynamic viscosity of the fluid in Pa·s, ρ-density in kg/m 3 , I-identity matrix, -specific heat capacity in J/(kg·K), Q-overall heat transfer in W, and T-the temperature in K. Moreover, = −k∇T, where k indicates the thermal conductivity of the fluid-solid mixture in W/m·K [20].
It should be noted that Equation (1) consists of inertia term (u⋅∇)u, a divergence of the stress-pI+μ(∇u+(∇u) Τ ); −∇pI stands for the pressure difference forces [38,39]. However, the selection of boundary conditions is not straightforward. The air pressure and air flux density decrease by passing through numerous cells. In order, to find the approximate boundary conditions which must be applied on the representative cell, the models containing 1 to 7 cells connected in sequence were analyzed (see Figure 3). In each model, the input ventilation rate 0.8 dm 3 /min was used as the inlet boundary condition on the surface Ωinlet. Relative pressure p0 = 0 Pa was applied as outlet boundary conditions Ωoutlet_1 and Ωoutlet_2. The mathematical expressions of these boundary conditions are presented in Table 2. Moreover, the fluid velocity on the 3D textile surface Ωt was set to zero. The constant temperature T = 37 °C was used as an inflow condition to simulate the heat temperature released by the skin surface Ωinflow. On Ωoutflow surface heat flux condition was applied that boundary condition requires heat transfer (h) coefficient that was obtained from the computational model described in Section 2.2.1. Figure 4a illustrates the positions of the boundary conditions for heat transfer on solids and fluids (.ht) interface for 1 cell. Figure  4b illustrates the positions of the boundary conditions for the laminar flow (.spf) interface. A Nonisothermal flow (.nitf) interface was used to combine the set of equations.  Outflow. In a model with convective heat transfer, this condition states that the only heat transfer occurring across the boundary is by convection. The temperature gradient in the normal direction is zero [40].
on Ωoutflow_2 − • = q q = h T − T Heat flux. The h denotes a heat transfer coefficient. The h = 8.9646 W/(K·m 2 ) was selected from Table 3, which was numerically predicted.
is the external temperature. It was assumed that T = 20 °C.

− • = 0
Thermal insulation. This boundary condition is used as the default.
Numerical simulations were performed under these boundary conditions (see Table  2) which were adapted for different cell numbers (see Figure 3). The analysis demonstrated that at a given value of the inlet air flux rate the inlet pressure and the air flux rate are different depending on the number of cells. The larger number of cells increases the resistance to airflow; therefore, the inlet pressure also becomes larger at the same value of the inlet air flux rate. Simultaneously, the air flux rate in each cell also decreases because of the flux loss through the openings at the top of each cell. This enables to select the boundary conditions which could be suitable for a representative cell. The boundary conditions on the representative cell are based on the air flux difference between the rightand left-hand sides of a cell, see Figure 5. The values were taken from air fluxes between the boundaries of the 6th and 7th cells where the relationship of air flux against the number of cells is close to linear. Table 3 represents the updated boundary conditions that were used in further simulations. It should be noted that, on Ωinlet the average volume flow rate of 0.000014716 m 3 /s was applied. This value was determined from the simulation that the computational domain consists of 6 cells. Moreover, on Ωoutlet_1 the average volume flow rate was set to 0.000014552 m 3 /s and it was obtained from the simulation that the computational domain consists of 7 cells.

Determination of Thermal Resistance
The thermal resistance of textile fabrics experimentally might be evaluated using different standard test methods such as the tog test (BS 4745), clo test (ASTM D-1518), sweating-guarded hotplate (ISO 11092), and THL (ASTM F1868) [3]. In this subsection, a computational model that predicts the thermal resistance of textile structures is based on the sweating-guarded hotplate test method. The numerical prediction of the thermal resistance Rct was performed following standard LST EN 31902/ISO 11092:2014(en) [31]. ISO 11092 is usually referred to as the "skin model" and imitates the heat transfer through the skin [3]. The thermal resistance of textile structures is defined as the temperature difference between the two faces of a material (Tm and Ta), divided by the resultant heat flux H per unit area A in the direction of the gradient (see Equation (4)) [31]. The measurement of the thermal resistance of the fabric requires that the temperature of the measuring unit is 35 °C, the ambient air temperature is (20 ± 2) °C, and the relative humidity is (65 ± 5)% [31,41].
where Tm-mean temperature of the measuring unit K, Ta-air temperature in the test enclosure K, H-heating power W required to maintain the plate temperature at T = 35 °C, ΔHc-correction term for heating power is 0, and Rct0-the apparatus constant m 2 ·K/W for the measurement of thermal resistance without sample [31,41]. The thermal resistance Rct coefficient was modeled according to the energy equation (see Equation (3)) [38]. On the inflow boundary, the temperature was set to T = 35 °C, and on the outflow boundary T = 20 °C. The computational domain is depicted in Figure 6. In this simulation, the air domain height is approximately equal to the textile domain height.
The computational model was done under the assumption that the flow is stationary nonisothermal flow, and the heat exchange occurs due to conduction. During simulations, it was assumed that the 3D textile is made from polyester and the fluid is moist air. The moist air properties might be obtained from the literature [42] or interface Heat transfer in solids and fluids (.ht) using specific media that represent the region of a moist air domain. The thermal resistance Rct coefficient can be defined from Equation (5) [24] and the heat transfer h coefficient from Equation (6) [41]. In addition, the thermal resistance Rct coefficient in Comsol Multiphysics environment might be determined as (aveop1(T)-aveop2(T))/aveop1(ht.ntflux), where aveop1(), aveop2() are the average operators on surfaces Ωinflow and Ωoutflow, T-temperature, ht.ntflux denotes normal total heat flux [43].

Mesh Analysis
The finite element mesh represents a solution field. By increasing the number of mesh elements, a more accurate approximation and solution can be obtained. However, it takes more time and requires more memory to solve the problem [40].
In this study, the optimal mesh size was determined by evaluating and comparing the numerically predicted air permeability (AP) coefficient with experimental data from the literature [34]. The proposed computation model of the air permeability (AP) coefficient was built on the assumption that airflow is an incompressible Newtonian flow and has a low Re number. The computational model is based on the Navier-Stokes (see Equation (1)) [38] and Brinkman equations (see Equation (7)) [38], as well as the continuity Equation (see Equation (2)) [38], where ρ, μ, p, εp, k are the density, dynamic viscosity, pressure, porosity, and local permeability, respectively. The finite element model follows the ISO 9237:1995 (E) standard [35], where the pressure difference between the inlet and outlet boundaries is 200 Pa. A more detailed explanation of boundary conditions is presented in the research conducted by the authors of this article [36].
In order, to determine an optimal mesh five different meshes (extremely coarse, extra coarse, coarser, coarse, and normal) were taken into consideration. Then, the numerical solution of AP was compared with experimental data [34] taking into account computation time and relative error.

Development of a Micro-Scale Model
The computational model of the forced ventilation was built using COMSOL Multiphysics software. The interfaces of Laminar flow (.spf) and Heat transfer in solids and fluids (.ht) were applied. The approach of finite element analysis is depicted in Figure 7. Lagrange linear discretization was applied for the temperature, velocity, and pressure. Moreover, different mesh size was taken into account during finite element analysis. The simulations were performed using the steady-state nonlinear solver. A simplified case of equations and pseudo code is presented below. Newton's method is the most widely used technique for solving nonlinear equations. Given : ℝ → ℝ , = , … , and the task is to find the solution to nonlinear systems of Equation (8) [44].
It should be noted that F is well defined and has continuous partial derivatives on an open set of ℝ , J(x) represents the Jacobian matrix [44]. It was assumed that solution x can be written as x = x0 + s, where x0 is some known guess of x, and s is a correction term.
In this study, Equation (10) [44] was solved using PARDISO solver, which was well documented in the literature [45]. The pseudocode of Newton's algorithm is adapted from the literature [46] and presented below in Algorithm 1.  α can be used, with 0 < α < 1 where α is a damping parameter.

The Macro-Scale Model
This subsection describes the finite element model of forced ventilation through the textile layer on a macro-scale. The proposed model is based on the literature [4,30]. Figure  8 illustrates the numerical scheme of water vapor and airflow, heat fluxes through a multilayered textile structure. Figure 8 was adapted from the literature [30]. The one-dimensional model consists of three nodes, two of which are positioned in the gaps at the bottom and top of the 3D textile layer, and the third node represents the inside of the ventilation layer [30]. The finite element model might be described by the structural Equation (12) [30]. That consists of the capacity matrix of the structure , vector of nodal variables = m , m , T , … , m , m , T . The final equation system of the proposed model in matrix form with the left-hand side capacity matrix of dimension 9 × 9 can be described by Equation (13) [30]. More details might be found in the literature [4,30]. The model was developed using Matlab software.
The forced ventilation model in the macro-scale requires input parameters such as α -heat transfer coefficient of the three-dimensional textile layer; α -heat transfer coefficient of the ventilation layer. These characteristics were determined using the microscale models described in previous sections.

Mesh Analysis
One of the most difficult aspects of developing computational models is determining the appropriate mesh size. In order, to obtain an optimal number of elements, computation time and accuracy were taken into account. This mesh analysis is essential for reducing computation time and demonstrating developed model accuracy. In this study, the relative error between the finite element solution and the experimental air permeability (AP) coefficient was determined. A RVE of Model 1_c was used and a mesh was applied to it. The numerical simulations were conducted with different meshes: extremely coarse, extra coarse, coarser, coarse, and normal, which represent 3194, 6249, 18,270, 57,795, and 125,512 elements, respectively. It should be noted that the minimum element quality was greater than 0.01, in order, to avoid convergence problems. It was found that using a coarser mesh relative error between the numerical solution and experimental data is 4.66%. Moreover, when coarse (57,795 elements) and normal (125,512) mesh are used, the relative error decreases to 0.25% and 0.16%, respectively. However, computation time increases from 9 s to 111 s (see Figure 9). Figure 9b, demonstrates that computation time vs. the number of elements increases according to the polynomial function. It was obtained that the coefficient of determination R-square is 0.9985. For further simulations, the coarse mesh was employed in the simplified geometry domain as the optimal mesh for numerical simulations.

Thermal Resistance Evaluation at a Micro-Scale
Numerical predictions of thermal resistance Rct and heat transfer h coefficients were determined for 3D textile structures corresponding to Model 1_e, Model 2_e, and Model 3_e. It was assumed that 3D textile is made from polyester and the fluid is moist air. Table 4 summarizes the determined thermal resistance Rct and heat transfer h coefficients. It was found that the heat transfer coefficient varies from 8.53 to 8.96 W/(K·m 2 ). It might be noticed that the heat transfer coefficient slightly increases with the increased thickness. In the literature [30], the measured thermal resistance Rct and heat transfer h coefficients of 3D textile are presented and equal to 0.118 K·m 2 /W and 8.4746 W/(K·m 2 ), respectively. However, the main differences may be due to differences in the height of samples, different space yarn configurations, initial material properties, and other factors. Figure 10 depicts temperature distribution through the height (z-axis) of Model 1_e. It should be noted that the heat transfer h coefficient corresponding to Model 3_e was used in a forced ventilation micro-scale model.

Combined Forced Ventilation Micro-and Macro-Scale Models
In Section 2.2, a forced ventilation model at a micro-scale was presented. The computational model allows for predicting the average temperature of the surface. The stationary temperature distributions via 3D textile structure, considering the flow rate of 0 and 0.000014716 m 3 /s were investigated. Figure 11 depicts the surface average temperature T °C and the position of the cut plane z, mm, when the flow rate is 0.000014716 m 3 /s. It was found that in the middle of microclimate (z = −0.5 mm), when the flow rate is 0 m 3 /s, the temperature is T = 35.81 °C, and with the flow rate 0.000014716 m 3 /s, T = 30.8147 °C. It was obtained that on Ωoutflow (z = 3.996 mm) surface the temperature is T = 27.81 °C (no ventilation) and T = 20.93 °C (Qsv = 0.000014716 m 3 /s). Figure 12 depicts temperature distributions when the mass flow rate is zero. Moreover, a three-dimensional view of the object is presented. The main tendency is that temperature decreases due to increased mass flow rate.
Furthermore, the heat transfer h coefficients were determined. The effective heat transfer coefficients are summarized in Table 5. The hmicrolayer stands for heat transfer between the skin surface and the bottom of the textile layer. It was obtained by using average operators and follows Equation (6) which is equivalent to the expression aveop1(ht.ntflux)/(aveop1(T)-aveop5(T)). In addition, htotal is the total heat transfer coefficient between the skin and the outer layer, and it is determined from the expression: aveop1(ht.ntflux)/(aveop1(T)-aveop2(T)). It was obtained that an increase in the mass flow rate caused the rise of the heat transfer coefficient. On the other hand, it indicates a decrease in thermal resistance.   The computational model was built on the assumption that the flow is laminar and incompressible. As a result, the model is applicable in the range of small flow rates. Furthermore, 0.000014716 m 3 /s corresponds to a velocity of 7.1549 m/s, which can be considered a windy condition [47].
Finally, established effective heat transfer coefficients were applied to the finite element model on a macro-scale. In this investigation,α = hmicrolayer = 25.603 W/(K·m 2 ), α = htotal = 8.9646 W/(K·m 2 ) were used from Tables 4 and 5. Figure 13 depicts the obtained results. The temperature predicted numerically was compared to experimental data from the literature [4]. The determined results were quantitative and qualitatively similar to the experimental data. The mean absolute error between the finite element model and measurements in the 3D textile layer was approximately 1.29 °C. Figure 13. Temperature comparison between numerically predicted and experimental temperatures in a 3D textile layer.

Conclusions
The main purpose of this study was to develop combined micro-and macro-scale models to investigate heat and mass transfer through textile structures with additional ventilation. In this study, the wearing comfort properties such as air permeability, thermal resistance, and heat transfer coefficients might be predicted at a micro-scale. Moreover, computational models were extended to the forced ventilation models at a micro-and a macro-scale. The effective heat transfer coefficients of the ventilation layer from the microscale model were applied to a macro-scale model. Both finite element models simulate ventilation through the 3D textile in both steady-state and/or time-dependent modes. The numerical results demonstrated a good agreement with the experimental data presented in the literature. It was found that the mean absolute error of temperature between numerically predicted on a macro-scale and experimental data in the 3D textile layer was approximately 1.29 °C. However, the developed forced ventilation model at a micro-scale is valid for low mass flow rates.
The numerical simulations of heat and mass transfer through textile structures with additional ventilation were performed using COMSOL Multiphysics and Matlab software.
New technological developments in the field of textiles provide a wide spectrum of materials with different physical properties. For this reason, theoretical analysis of thermal comfort properties is essential for ensuring wearing comfort. The computational models proposed in this study can be used in the early stages of designing functional clothing such as protective clothing. Additionally, the proposed theoretical approach might be applied to other scientific and engineering purposes related to heat and mass transfer through three-dimensional structures.