Experimentally Veriﬁed Flow Distribution Model for a Composite Modelling System

: Requirements of modern process and power technologies for compact and highly efﬁcient equipment for transferring large heat ﬂuxes lead to designing these apparatuses as dense parallel ﬂow systems, ranging from conventional to minichannel dimensions according to the speciﬁc industrial application. To avoid operating issues in such complex equipment, it is vital to identify not only the local distribution of heat ﬂux in individual parts of the heat transfer surface but also the uniformity of ﬂuid ﬂow distribution inside individual parallel channels of the ﬂow system. A composite modelling system is currently being developed for accurate design of such complex heat transfer equipment. The modeling approach requires a ﬂow distribution model enabling to yield accurate-enough predictions in reasonable time frames. The paper presents the results of complex experimental and modeling investigation of ﬂuid ﬂow distribution in dividing headers of tubular-type equipment. Different modeling approaches were examined on a set of header geometries. Predictions obtained via analytical and numerical models were validated using data from the experiments conducted on additively manufactured header samples. Two case studies employing parallel ﬂow systems (mini-scale systems and a conventional-size heat exchanger) demonstrated the applicability of the distribution model and the accuracy of the composite modelling system. areas of 30 and 8 mm, respectively. Tubes that were 50 mm long were regularly spaced along the headers, with the relative length (L/D ratio) being kept at 10. The chosen L/D ratio in combination with the small tube pitch, on the one hand, had a reducing effect on friction in the distribution system; on the other hand, it negatively impacted the ﬂuid ﬂow maldistribution [44]. To conclude, although the geometries of the examined systems were relatively simple, their main dimensionless parameters (such as the L/D and cross-sectional ratios) were very close to those encountered in industrial practice.


Introduction
Various industrial applications demand higher compactness, improved energy efficiency, and prolonged service life of equipment for transferring large heat fluxes. This inevitably leads to designing these crucial components as dense parallel flow systems. Dimensions of such parallelized systems then vary from large scales encountered, e.g., in process and power industries (such as heat recovery boilers [1] or fin-and-tube heat exchanger [2]) to mini-and micro-scale apparatuses (e.g., minichannels [3], or waste heat recovery micro-units [4]) according to the needs of the respective application. Although the operating conditions can be radically different, a common feature is the reported nonuniformity of flow rates, also known as flow (mal)distribution, in individual parts of the complex systems. Accurate prediction of flow distribution within the heat exchanger tube and shell sides is one of the crucial aspects that directly influence the equipment's reliability and performance.
On the one hand, flow maldistribution is, in general, an undesirable effect of a nonuniform pressure field in a heat exchanger. Various geometry modifications are often employed to avoid such non-uniformity. One can encounter orifices in tubes [5], baffles in headers [6], non-standard header shapes (trapezoidal [7], concave and convex [8], or with triangular cross-section [9]), manifolds of different dimensions [10], or different positioning of the inlet and outlet zones in relation to the headers [11] Additionally, many researchers have reported on the influence of the system arrangement on the flow behavior (e.g., the review covering microchannel applications [12] by Ghani et al.). The most common are U- [13] and Z-arrangements [14], but steam generators often utilize other types of arrangements as well [15]. Moreover, the distribution systems may be very complex and employ a remarkable variety of channel shapes, as numerous studies reveal [16][17][18]. On the other hand, some authors intentionally developed uneven flow rates [19] in channels, or they took advantage of natural flow maldistribution as it may compensate non-uniform heating conditions in a mini-scale heat sink and mitigate hot spots [20].
Recently, new possibilities opened up by the development of minigap systems and by the usage of additive manufacturing technologies.
As for minigaps, such a single wide minichannel is a progressive option for dealing with high heat fluxes, mainly in electronic systems. Investigations into minigaps of different geometries (e.g., rectangular [21] or annular [22]) usually focus on the thermal performance (e.g., local heat transfer coefficients [23], Nusselt number [24], or temperature distribution [25]). Flow distribution in minigaps is discussed rarely, and if so, then usually as a side experimental observation [26] or in the form of numerical analyses using computational fluid dynamics (CFD) [27]. Several studies have also compared the performance of both mini-scale applications, i.e., minichannels and minigaps. Klugmann et al. [28] carried out experiments on minichannel and minigap evaporators of unusually large width and length. The standard dimensions combined with small thickness are supposed to widen the application possibilities of such heat exchangers. The authors concluded that pressure drop increased more significantly with a higher mass flow rate in the minichannel system, yet the obtained heat transfer coefficient was also higher in the minichannel system than in the minigap of the same size. Results of Dąbrowski's analyses [29] then showed similarities in flow distribution trends between minichannel and minigap systems, although minigap systems provided much worse flow uniformity.
Additive manufacturing enables adjustments of a wide range of geometrical parameters at once in order to yield better flow distribution in the system [30]. However, the benefits of additive technology lie not only in the relative freedom of geometrical complexity [31,32] but also in the ease of passive enhancement of heat transfer. For instance, the inner heat transfer surfaces, which could otherwise be modified only with considerable difficulty (e.g., [33][34][35]), can be easily enhanced in the process of additive fabrication, as shown, e.g., by Wei et al. [36] or Wu et al. [37].
Despite the fact that previous theoretical and experimental studies provided valuable insight into flow maldistribution in heat transfer equipment, modeling and predicting this phenomenon remains a challenge, especially given the necessary compromise between accuracy and computational cost. This has been proven absolutely necessary in the case of initial equipment design. Furthermore, the situation is also made more complicated by the lack of a standardized methodology for maldistribution evaluation. All these factors contribute to the common but erroneous assumptions of the flow being distributed evenly. Later on, such an approach inevitably leads to operating problems that may be the cause of an abrupt shutdown and damage to the apparatus itself or the associated devices. An ex-post solution, i.e., troubleshooting of the inappropriately designed equipment, is usually a costly and time-intensive procedure [38]. Still, the situation can be preempted by a proper design procedure that would respect flow maldistribution. This is also the main objective of the composite modelling system (CMS) that was first proposed in [39].
CMS is intended for initial design computations of heat transfer equipment; however, it is also suitable for quick equipment rating. The main principle is briefly recapitulated in Figure 1 (for a detailed description of CMS, the readers are kindly referred to [40]). It should be mentioned that CMS is still in development. Currently, the focus is on finding a suitable flow distribution model that would meet essential requirements: fast computation, sufficient accuracy, and easy implementation into the whole concept of CMS. Two differential models were selected for the final decision. Differential models can yield highly accurate results not only in terms of flow distribution but also the temperature and pressure distribution in the heat transfer equipment (see, e.g., [41]). However, as Turek [42] Energies 2021, 14, 1778 3 of 24 pointed out, with the increasing complexity of geometry also raises the computational time needed to obtain the solution of differential equations. highly accurate results not only in terms of flow distribution but also the temperature and pressure distribution in the heat transfer equipment (see, e.g., [41]). However, as Turek [42] pointed out, with the increasing complexity of geometry also raises the computational time needed to obtain the solution of differential equations.

Figure 1.
Basic principle of composite modelling system, which uses a set of interconnected submodels. The key objective of this work, i.e., to identify a suitable flow distribution model, is highlighted on the right side.
A closer investigation into the working of the distribution models revealed rather opposite results depending on geometry parameters of distribution systems [43]. This paper builds upon the previous work presented in [43] and provides a deeper view of the research into the modeling of flow distribution via analytical, numerical, and experimental means.
The main objective of the first part is the identification of the most suitable distribution model. At first, a thorough revision to both distribution models' performance is introduced. Then, predictions of flow distribution obtained using analytical models and results of CFD analyses are validated against data from physical experiments carried out on additively manufactured distributors. In the second part of this work, the selected (bestfitting) analytical model is utilized for modeling flow distribution in mini-scale systems as well as in equipment of conventional size. A preliminary study involving the possible application of the flow distribution model to minichannel and minigap systems is discussed. The suggested fast modeling approach may be employed to improve the prediction of heat transfer in mini-scale applications. Lastly, through an industrial case, the selected modeling approach is demonstrated in the broad context of CMS. The respective examples show the importance of flow distribution in heat transfer equipment and the risk of neglecting even minor maldistribution, which is especially common in engineering practice.

Materials and Methods
The 3 most common modeling approaches-analytical, numerical, and experimental-utilized for an investigation into flow distribution are described in the following subsections.
Analytical models, which are covered first, are based on partial differential equations that require several simplifications to yield the solution. Thus, to provide a reliable comparison of the modeling approaches, the same simplifying assumptions were introduced for all models: 1. steady-state conditions; Figure 1. Basic principle of composite modelling system, which uses a set of interconnected submodels. The key objective of this work, i.e., to identify a suitable flow distribution model, is highlighted on the right side.
A closer investigation into the working of the distribution models revealed rather opposite results depending on geometry parameters of distribution systems [43]. This paper builds upon the previous work presented in [43] and provides a deeper view of the research into the modeling of flow distribution via analytical, numerical, and experimental means.
The main objective of the first part is the identification of the most suitable distribution model. At first, a thorough revision to both distribution models' performance is introduced. Then, predictions of flow distribution obtained using analytical models and results of CFD analyses are validated against data from physical experiments carried out on additively manufactured distributors. In the second part of this work, the selected (best-fitting) analytical model is utilized for modeling flow distribution in mini-scale systems as well as in equipment of conventional size. A preliminary study involving the possible application of the flow distribution model to minichannel and minigap systems is discussed. The suggested fast modeling approach may be employed to improve the prediction of heat transfer in mini-scale applications. Lastly, through an industrial case, the selected modeling approach is demonstrated in the broad context of CMS. The respective examples show the importance of flow distribution in heat transfer equipment and the risk of neglecting even minor maldistribution, which is especially common in engineering practice.

Materials and Methods
The 3 most common modeling approaches-analytical, numerical, and experimentalutilized for an investigation into flow distribution are described in the following subsections.
Analytical models, which are covered first, are based on partial differential equations that require several simplifications to yield the solution. Thus, to provide a reliable comparison of the modeling approaches, the same simplifying assumptions were introduced for all models: 1. steady-state conditions; 2.
constant physical properties of the working fluid-water-in the measured temperature range.
As will be discussed further, 3 flow distribution systems were chosen. These systems differed primarily in the number of tubes, i.e., 7, 14, and 27. Figures 2-4 show schematics of the systems (further denoted as n07, n14, and n27, respectively) including the most important geometrical parameters. Headers and tubes were of uniform cross-sectional areas of 30 and 8 mm, respectively. Tubes that were 50 mm long were regularly spaced along the headers, with the relative length (L/D ratio) being kept at 10. The chosen L/D ratio in combination with the small tube pitch, on the one hand, had a reducing effect on friction in the distribution system; on the other hand, it negatively impacted the fluid flow maldistribution [44]. To conclude, although the geometries of the examined systems were relatively simple, their main dimensionless parameters (such as the L/D and cross-sectional ratios) were very close to those encountered in industrial practice.
2. constant physical properties of the working fluid-water-in the measured temperature range.
As will be discussed further, 3 flow distribution systems were chosen. These systems differed primarily in the number of tubes, i.e., 7, 14, and 27. Figures 2-4 show schematics of the systems (further denoted as n07, n14, and n27, respectively) including the most important geometrical parameters. Headers and tubes were of uniform cross-sectional areas of 30 and 8 mm, respectively. Tubes that were 50 mm long were regularly spaced along the headers, with the relative length (L/D ratio) being kept at 10. The chosen L/D ratio in combination with the small tube pitch, on the one hand, had a reducing effect on friction in the distribution system; on the other hand, it negatively impacted the fluid flow maldistribution [44]. To conclude, although the geometries of the examined systems were relatively simple, their main dimensionless parameters (such as the L/D and cross-sectional ratios) were very close to those encountered in industrial practice.   2. constant physical properties of the working fluid-water-in the measured temperature range.
As will be discussed further, 3 flow distribution systems were chosen. These systems differed primarily in the number of tubes, i.e., 7, 14, and 27. Figures 2-4 show schematics of the systems (further denoted as n07, n14, and n27, respectively) including the most important geometrical parameters. Headers and tubes were of uniform cross-sectional areas of 30 and 8 mm, respectively. Tubes that were 50 mm long were regularly spaced along the headers, with the relative length (L/D ratio) being kept at 10. The chosen L/D ratio in combination with the small tube pitch, on the one hand, had a reducing effect on friction in the distribution system; on the other hand, it negatively impacted the fluid flow maldistribution [44]. To conclude, although the geometries of the examined systems were relatively simple, their main dimensionless parameters (such as the L/D and cross-sectional ratios) were very close to those encountered in industrial practice.

Analytical Models
The 2 examined flow distribution models were introduced by Bajura [45] and Bajura and Jones [46], respectively. Normalized flow rate QN (i.e., the ratio of the header flow rate and the inlet flow rate) in the distribution system is expressed in these models as where x is the non-dimensional distance along the header, Φ1 and Φ2 the overall friction loss coefficients, and M1 and M2 describe the changes in fluid momentum. If the distribution system in question consists of a dividing header and a bundle of tubes, the coefficients in Equation (1) can be substituted by expressions listed in Table 1. For the sake of simplicity, Bajura's model is further denoted as Model B, the model by Bajura and Jones as Model BJ.
Ar stands for the cross-sectional ratio, DD the header diameter, fD friction factor, H the lateral flow resistance, LD header length, γD the static pressure regain coefficient, and θD is the momentum coefficient for the header flow. For a detailed description of deriving the respective coefficient, please refer to either of the two papers.

Analytical Models
The 2 examined flow distribution models were introduced by Bajura [45] and Bajura and Jones [46], respectively. Normalized flow rate Q N (i.e., the ratio of the header flow rate and the inlet flow rate) in the distribution system is expressed in these models as d 2 where x is the non-dimensional distance along the header, Φ 1 and Φ 2 the overall friction loss coefficients, and M 1 and M 2 describe the changes in fluid momentum. If the distribution system in question consists of a dividing header and a bundle of tubes, the coefficients in Equation (1) can be substituted by expressions listed in Table 1. For the sake of simplicity, Bajura's model is further denoted as Model B, the model by Bajura and Jones as Model BJ.
A r stands for the cross-sectional ratio, D D the header diameter, f D friction factor, H the lateral flow resistance, L D header length, γ D the static pressure regain coefficient, and θ D is the momentum coefficient for the header flow. For a detailed description of deriving the respective coefficient, please refer to either of the two papers.
Equation (1) can, therefore, be rewritten as For Models B and BJ, respectively. The boundary conditions of both analytical models are Despite the fast performance and simplicity of the models, which enables their easy implementation into a computational procedure, the models suffer from an obvious drawback. Namely, the flow distribution prediction depends on empirical coefficients, whose recommended value ranges are defined in [45,46]. Furthermore, as can also be seen in Equations (6) and (7) above, the dependence on cross-sectional ratio A r (9) that compares the outlet areas A T of n tubes and the inlet area A D of the distribution system is entirely different in each model.
Before proceeding to the experimental work (whose purpose was to identify the more accurate model), we carried out a sensitivity analysis to estimate the influence of A r and the lateral flow resistance H. Additionally, the response of the models to the changes of the static pressure regain coefficient γ D and the momentum coefficient for the header flow θ D observed in a matter of values recommended by authors [45,46], respectively. Table 2 summarizes the parameters of the tested distribution system. The effect of a single parameter on flow distribution was determined in each case shown in Figure 5; others were kept constant as listed for the basic setup. Results of the conducted sensitivity study clearly demonstrate that A r is the most influential parameter (especially in the case of Model BJ), which agrees with the flow distribution-related findings presented, e.g., in [14] or [44]. However, both coefficients have a relatively small effect on flow distribution. Predicted flow maldistribution grows as the coefficient γ D decreases in the case of Model B. The same effect on non-uniformity of the fluid flow can be observed with the increasing coefficient θ D if Model BJ is utilized. Conversely, the effect of changes in lateral flow resistance H is almost identical in both models, i.e., with a higher H (due to extensive length of the tubes or relatively small tube diameters), distribution of fluid becomes more uniform. In consequence of the data obtained in the sensitivity study, a set of 3 representative geometries (with cross-sectional ratios approximately 0.5, 1.0, and 2.0) was selected for a more thorough investigation of the behavior and accuracy of Models B and BJ.

CFD Models
Numerical simulations were chosen as the second approach to investigate flow behavior in distribution systems. 3D CFD models can offer a detailed view into flow geometry, especially in the critical parts of the distribution systems where fluid properties and flow parameters would be measurable with much difficulty or not at all [43]. CFD also provides a convenient comparison to the performance of simpler mathematical models described in the previous part.  CFD simulations (including the computational grid preparation) were carried out in commercial software Ansys Fluent [47]. Regarding the model assumptions specified above, the following CFD setup was used: • Pressure-based solver with absolute velocity formulation and double-precision; • Enhanced wall treatment with realizable k-ε model; • SIMPLE algorithm for pressure-velocity coupling; • Green-Gauss node based gradient calculation; • Spatial discretization: second order for pressure, second order upwind for density, momentum, turbulent kinetic energy, and turbulent dissipation rate.
As for boundary conditions, gauge pressure at the outlet was set to 0 Pa, while a mass flow rate boundary condition was used at the inlet. Given the set of test cases, the target water velocity in the header was 1 m/s. This corresponded to the inlet mass flow rate of approximately 0.7056 kg/s. Otherwise, the inlet mass flow rate was set with respect to the flow rates observed during laboratory experiments. Only the no-slip condition was applied at the wall boundary since wall roughness cannot be specified for the utilized near-wall modeling method.
The most significant feature of the presented CFD models was the usage of the porous jump boundary condition (PJBC) [47] (p. 7.3.21) in the tubes. As a result, the modeled tubes could be shortened in comparison with the experimental apparatus. It should be noted that each flow system had a specific setup of PJBC according to the length of the substituted tube portions and average fluid velocities therein. This modeling simplification led to considerably smaller cell counts in the areas where gradients of flow variables were likely to be moderate. On the other hand, more validation steps are still necessary to obtain reliable and high-quality results.
Another technique applied to reduce the cell counts even further with the same (or even better) cell quality in the fluid domains was the variable mesh density. Apart from near-wall regions, the most refined mesh was created around the tube inlets where the highest gradients of flow variables were expected. Figure 6 shows an example of such a grid.

CFD Models
Numerical simulations were chosen as the second approach to investigate flow behavior in distribution systems. 3D CFD models can offer a detailed view into flow geometry, especially in the critical parts of the distribution systems where fluid properties and flow parameters would be measurable with much difficulty or not at all [43]. CFD also provides a convenient comparison to the performance of simpler mathematical models described in the previous part.
CFD simulations (including the computational grid preparation) were carried out in commercial software Ansys Fluent [47]. Regarding the model assumptions specified above, the following CFD setup was used: • Pressure-based solver with absolute velocity formulation and double-precision; • Enhanced wall treatment with realizable k-ε model; • SIMPLE algorithm for pressure-velocity coupling; • Green-Gauss node based gradient calculation; • Spatial discretization: second order for pressure, second order upwind for density, momentum, turbulent kinetic energy, and turbulent dissipation rate.
As for boundary conditions, gauge pressure at the outlet was set to 0 Pa, while a mass flow rate boundary condition was used at the inlet. Given the set of test cases, the target water velocity in the header was 1 m/s. This corresponded to the inlet mass flow rate of approximately 0.7056 kg/s. Otherwise, the inlet mass flow rate was set with respect to the flow rates observed during laboratory experiments. Only the no-slip condition was applied at the wall boundary since wall roughness cannot be specified for the utilized nearwall modeling method.
The most significant feature of the presented CFD models was the usage of the porous jump boundary condition (PJBC) [47] (p. 7.3.21) in the tubes. As a result, the modeled tubes could be shortened in comparison with the experimental apparatus. It should be noted that each flow system had a specific setup of PJBC according to the length of the substituted tube portions and average fluid velocities therein. This modeling simplification led to considerably smaller cell counts in the areas where gradients of flow variables were likely to be moderate. On the other hand, more validation steps are still necessary to obtain reliable and high-quality results.
Another technique applied to reduce the cell counts even further with the same (or even better) cell quality in the fluid domains was the variable mesh density. Apart from near-wall regions, the most refined mesh was created around the tube inlets where the highest gradients of flow variables were expected. Figure 6 shows an example of such a grid. Figure 6. Example of a mesh with variable density (mesh arrangement "03") at the entrance of a distributor. Extra inlet and outlet volumes were added to the dividing systems to avoid computational issues due to reversed flow (i.e., a light blue mesh visible on this figure's far-right belongs to an inlet zone). Example of a mesh with variable density (mesh arrangement "03") at the entrance of a distributor. Extra inlet and outlet volumes were added to the dividing systems to avoid computational issues due to reversed flow (i.e., a light blue mesh visible on this figure's far-right belongs to an inlet zone).
The first step of the validation process was adjusting the mesh and porous jump setting in a single tube to yield the same pressure loss as in a 1.55 m long tube (as per the setup of the physical experiments). The necessary data were obtained using CFD analysis of flow in a full-scale model of tube employing a fine-tuned mesh (5.5 M cells) and Shear-Stress Transport k-ω viscous model.
In the second step, similarly to the previous subsection, a geometry with 14 tubes was utilized for the mesh independence test. The examined polyhedral grids were based Energies 2021, 14, 1778 9 of 24 on 6 promising tube mesh arrangements differing mostly in the level of refinement near the walls. Details including the mesh quality metrics are listed in Table 3. Since mass flow rate distribution proved ineffective for the evaluation of the effect of mesh size (see, e.g., [9] or [14]), the overall pressure drop was taken as the metric instead. As Figure 7 clearly shows, the overall pressure drop in the distribution system converged to the value of ≈4680 Pa; thus, the best compromise between grid size and accuracy was achieved using mesh "04". The first step of the validation process was adjusting the mesh and porous jump setting in a single tube to yield the same pressure loss as in a 1.55 m long tube (as per the setup of the physical experiments). The necessary data were obtained using CFD analysis of flow in a full-scale model of tube employing a fine-tuned mesh (5.5 M cells) and Shear-Stress Transport k-ω viscous model.
In the second step, similarly to the previous subsection, a geometry with 14 tubes was utilized for the mesh independence test. The examined polyhedral grids were based on 6 promising tube mesh arrangements differing mostly in the level of refinement near the walls. Details including the mesh quality metrics are listed in Table 3. Since mass flow rate distribution proved ineffective for the evaluation of the effect of mesh size (see, e.g., [9] or [14]), the overall pressure drop was taken as the metric instead. As Figure 7 clearly shows, the overall pressure drop in the distribution system converged to the value of ≈4680 Pa; thus, the best compromise between grid size and accuracy was achieved using mesh "04".  Considering the convergence of the solution process, it was monitored via scaled residuals of the usual variables (all set to 10 -3 ) and static pressure, mass flow rate, and magnitude of the average facet velocity. The former 2 quantities were monitored on the tube cross-sectional faces, while the facet velocity monitors were located at 4 points shown in Figure 8. A large number of iterations was needed to reach a converged average facet velocity; therefore, unless a convergent solution was reached, a threshold of the scaled residual of continuity was lowered to 10 -6 . Considering the convergence of the solution process, it was monitored via scaled residuals of the usual variables (all set to 10 −3 ) and static pressure, mass flow rate, and magnitude of the average facet velocity. The former 2 quantities were monitored on the tube cross-sectional faces, while the facet velocity monitors were located at 4 points shown in Figure 8. A large number of iterations was needed to reach a converged average facet velocity; therefore, unless a convergent solution was reached, a threshold of the scaled residual of continuity was lowered to 10 −6 .

Experiments
Three experimental distributors were additively manufactured via Fused Deposition Modelling technology from the ABS (acrylonitrile butadiene styrene) material using a Trilab DeltiQ XXL machine (TriLAB Group s.r.o., Brno, Czech Republic). This printer features an 800 mm tall print area with a diameter of 250 mm. This allowed the distributors to be manufactured in a vertical position, which positively impacted the quality of printed parts. To minimize the inner surface imperfections, we tilted the distributors at the angle of 10 • between the header's longitudinal axis and the printing direction (see Figure 9). Such orientation also significantly reduced residual stresses in the material; therefore, distortions of the final parts were minimal. The distributors were printed with a layer height of 0.12 mm and 100% infill, which ensured sufficient stiffness. During the leakage tests, the printed distributors were filled with water for over 3 h and no leakage was observed; i.e., they were found waterproof. The build time differed from 29 h and 28 min for the n07 configuration to 53 h and 25 min for the n27 configuration. After additive manufacturing, the supporting structures were removed, and the inner diameters of the tubes were machined.

Experiments
Three experimental distributors were additively manufactured via Fused Deposition Modelling technology from the ABS (acrylonitrile butadiene styrene) material using a Trilab DeltiQ XXL machine (TriLAB Group s.r.o., Brno, Czech Republic). This printer features an 800 mm tall print area with a diameter of 250 mm. This allowed the distributors to be manufactured in a vertical position, which positively impacted the quality of printed parts. To minimize the inner surface imperfections, we tilted the distributors at the angle of 10° between the header's longitudinal axis and the printing direction (see Figure 9). Such orientation also significantly reduced residual stresses in the material; therefore, distortions of the final parts were minimal. The distributors were printed with a layer height of 0.12 mm and 100% infill, which ensured sufficient stiffness. During the leakage tests, the printed distributors were filled with water for over 3 hours and no leakage was observed; i.e., they were found waterproof. The build time differed from 29 h and 28 min for the n07 configuration to 53 h and 25 min for the n27 configuration. After additive manufacturing, the supporting structures were removed, and the inner diameters of the tubes were machined.

Experiments
Three experimental distributors were additively manufactured via Fused Deposition Modelling technology from the ABS (acrylonitrile butadiene styrene) material using a Trilab DeltiQ XXL machine (TriLAB Group s.r.o., Brno, Czech Republic). This printer features an 800 mm tall print area with a diameter of 250 mm. This allowed the distributors to be manufactured in a vertical position, which positively impacted the quality of printed parts. To minimize the inner surface imperfections, we tilted the distributors at the angle of 10° between the header's longitudinal axis and the printing direction (see Figure 9). Such orientation also significantly reduced residual stresses in the material; therefore, distortions of the final parts were minimal. The distributors were printed with a layer height of 0.12 mm and 100% infill, which ensured sufficient stiffness. During the leakage tests, the printed distributors were filled with water for over 3 hours and no leakage was observed; i.e., they were found waterproof. The build time differed from 29 h and 28 min for the n07 configuration to 53 h and 25 min for the n27 configuration. After additive manufacturing, the supporting structures were removed, and the inner diameters of the tubes were machined. Figure 9. Manufacturing of the n07 distributor. Figure 9. Manufacturing of the n07 distributor.
As mentioned above, the surface quality of additively manufactured parts highly depends on surface orientation during fabrication. Therefore, roughnesses of the header and tube inner surfaces were investigated separately. Measurements were taken using a digital microscope Keyence VHX-6000 (Keyence International NV/SA, Mechelen, Belgium). Three profile measurements were made for both samples-header and tube-and the maximum height values Rz (i.e., maximum peak to valley height of the profile, within a single sampling length) were obtained. The profile line was always oriented perpendicular to the layers, as shown in the example in Figure 10. Figure 11 then shows an example of the surface topography in the header of the n14 distributor. The roughness of the surface was then obtained as the average value, i.e., Rz was 53.83 and 33.98 µm in the header and the tubes, respectively. Although the tubes needed support structures due to their orientation and worse surface quality was expected, the roughness value was lower because of the subsequent machining. single sampling length) were obtained. The profile line was always oriented perpendicular to the layers, as shown in the example in Figure 10. Figure 11 then shows an example of the surface topography in the header of the n14 distributor. The roughness of the surface was then obtained as the average value, i.e., Rz was 53.83 and 33.98 µm in the header and the tubes, respectively. Although the tubes needed support structures due to their orientation and worse surface quality was expected, the roughness value was lower because of the subsequent machining.   maximum height values Rz (i.e., maximum peak to valley height of the profile, within a single sampling length) were obtained. The profile line was always oriented perpendicular to the layers, as shown in the example in Figure 10. Figure 11 then shows an example of the surface topography in the header of the n14 distributor. The roughness of the surface was then obtained as the average value, i.e., Rz was 53.83 and 33.98 µm in the header and the tubes, respectively. Although the tubes needed support structures due to their orientation and worse surface quality was expected, the roughness value was lower because of the subsequent machining.     Figure 12b then shows a photograph taken during experiments performed on the n27 apparatus. The installed measuring devices and their maximum uncertainties (according to manufacturers' specifications) are listed in Table 4. Table 4 also contains the operating condition ranges for all the experiments. It should be mentioned that the source of water could supply at most 5.5 m 3 /h at a maximum pressure of 3.8 bar (g). m 3 /h at a maximum pressure of 3.8 bar (g).
The tested distribution systems contained a relatively large number of tubes; hence, the flow rate measurement had to be as simple as possible, albeit accurate enough, in order to reduce the cost while still sustaining the results' reliability. The fluid flow distribution into individual tubes was measured via cumulative flow rate within a 90-s time span. The measurements were repeated 3 times for each experimental setup. Consequently, the respective averaged value of flow rate was specified for each tube.

Results and Discussion
Predictions of flow rate distribution yielded by the analytical models B (6) and BJ (7) and the data obtained via CFD simulations are analyzed in this section. The mathematical models (analytical as well as the detailed numerical ones) were validated using flow rates from the physical experiments.
Considering the evaluation of flow distribution non-uniformity, the methodology is not standardized. Several different approaches are, therefore, utilized by researchers. The evaluation methods can be divided into two groups.
1. Continuous criteria, e.g., dimensionless (also called normalized) flow rate [13] or relative extent of non-uniformity [48], describe the distribution of flow rate by a set of values for each branch of the distribution system. 2. One-value criteria provide a single value determining flow maldistribution for the entire system, usually based on flow velocities, mass flow rates, or volumetric flow  The tested distribution systems contained a relatively large number of tubes; hence, the flow rate measurement had to be as simple as possible, albeit accurate enough, in order to reduce the cost while still sustaining the results' reliability. The fluid flow distribution into individual tubes was measured via cumulative flow rate within a 90-s time span. The measurements were repeated 3 times for each experimental setup. Consequently, the respective averaged value of flow rate was specified for each tube.

Results and Discussion
Predictions of flow rate distribution yielded by the analytical models B (6) and BJ (7) and the data obtained via CFD simulations are analyzed in this section. The mathematical models (analytical as well as the detailed numerical ones) were validated using flow rates from the physical experiments.
Considering the evaluation of flow distribution non-uniformity, the methodology is not standardized. Several different approaches are, therefore, utilized by researchers. The evaluation methods can be divided into two groups.

1.
Continuous criteria, e.g., dimensionless (also called normalized) flow rate [13] or relative extent of non-uniformity [48], describe the distribution of flow rate by a set of values for each branch of the distribution system. 2.
One-value criteria provide a single value determining flow maldistribution for the entire system, usually based on flow velocities, mass flow rates, or volumetric flow rates. Examples of such quantification can be the minimum to maximum ratio, ratio of the minimal (maximal) to average value [49], maximum deviation from uniform distribution [50], or relative standard deviation from uniform distribution, to name just a few.
Despite the useful outline of flow distribution along the header that can be given by continuous criteria, one-value criteria are more favorable for fast decisions as to whether the flow maldistribution is or is not at an acceptable level. One of the well-known criteria is non-uniformity (NU) percentage (10); however, probably the most commonly used flow distribution criterion is the relative standard deviation (RSD) of the tube flow rate (used, e.g., in [1]). A significant benefit of RSD (11) is its suitability even if backflow occurs in some of the distribution system branches [51]. Moreover, because this method determines the extent of variability in relation to the mean value, it allows one to compare systems of different sizes (numbers of tubes).
To provide a complex description of the examined flow distribution systems, we determined the non-uniformity of volumetric flow rates in tubes Q Ti by both methods (Equations (10) and (11)). Please note that Q T,id denotes the ideal (mean) flow rate in a single tube, and n stands for the number of tubes. In the case of both criteria, a value closer to zero represents a lower level of flow maldistribution, and hence a more suitable flow system configuration. Table 5 lists the geometry and flow parameters required particularly by Models B and BJ. The comparison of flow rate distribution predicted by mathematical models is shown together with experimental data in Figures 13-15 for each distributor configuration separately.  Figure 13 visually confirms that the most uniform distribution was reached in the n07 system. This also agrees with the assumptions mentioned in Section 2, i.e., that small crosssectional ratio and lateral flow resistance are mitigating factors. Accordingly, considerably higher maldistribution was observed in the n14 and n27 systems, whose volumetric flow rates are presented in Figures 14 and 15, respectively. Furthermore, remarkably low flow rates were measured in the third and seventh tubes of the n14 configuration, similarly to the first and second tubes of the n27 configuration. Nevertheless, no flow obstacles or surface flaws were observed in either distributor; thus, it can be concluded that a larger data set would be necessary to discover the root cause of the outlying flow rates. Table 6 summarizes the observed levels of flow non-uniformity. Apart from the RSD and NU criteria, the maximum relative differences between the measured and predicted tube flow rate are listed as well. However, especially in the case of the n27 distributor, these discrepancies mainly reflect on the outlying flow rates mentioned before.         Table 6 summarizes the observed levels of flow non-uniformity. Apart from the RSD and NU criteria, the maximum relative differences between the measured and predicted tube flow rate are listed as well. However, especially in the case of the n27 distributor, these discrepancies mainly reflect on the outlying flow rates mentioned before.   Table 5); BJ = Analytical Model BJ (best fitting, see Table 5); CFD = Steady CFD model; E = experimental data.

Discussion
Results presented in Section 3.1 show that, if appropriately fitted, analytical models can reach a comparable or even better agreement with experimental data than detailed numerical simulations. In general, the best-matching predictions of flow distribution were obtained using the fitted Model BJ. If we exclude the outliers from the flow rate data sets, the relative difference between the experimental and predicted RSDs does not exceed 8%. If the recommended setup of Model BJ is used, the discrepancy is 3% to 13% (Model BJ underestimates the actual flow non-uniformity). Conversely, Model B provides a more accurate result in the basic configuration n14 without any modifications. This model's predictive ability, however, is unsatisfactory in a wider range of geometries, as proven by the described experiments.
A tremendous advantage of the analytical models-their extremely short evaluation times-stands out, particularly when compared to the times required by CFD simulations to yield steady-state solutions (units of hours plus the extra time needed for mesh preparation and post-processing). Moreover, the observed differences between experimental data and the results of the CFD analyses suggest that steady simulations cannot offer a sufficient degree of accuracy. The authors of a previous work that focused on dense tube bundles [52] argued that this traditional modeling approach is insufficient and recommended transient simulations in the case of more complex systems. The findings of this study lead to a similar conclusion.
As a result, analytical Model BJ appears to be a reasonable trade-off between accuracy of the flow distribution prediction and computational cost, especially given the intended purpose, that is, utilizing the distribution model in a composite modeling system for initial design of heat transfer equipment.

Case Studies
The following text presents two examples in which the analytical Model BJ was employed to predict flow distribution. The first study concerns a possible utilization of the distribution model in mini-scale systems. Performance of Model BJ is evaluated using data from open literature. The second, rather conventional application showcases the distribution model in the context of the composite modelling system.

Minichannel/Minigap Systems
Knowledge of the flow distribution is essential in the mitigation of hot spots in miniscale heat transfer equipment because any larger non-uniformity may fatally influence the associated devices [53]. However, minigap and minichannel systems often work in a laminar flow regime, which is not the originally intended application of the model proposed by Bajura and Jones [46]. The applicability of Model BJ in the minisystem field is hereby examined.
Regarding minichannels, two distribution systems with the U-and Z-arrangements were considered. Validation data were derived from CFD simulations performed by Kumar and Singh [20]. The modeled systems with 28 parallel rectangular channels were of a total area of 55 × 47 mm, and flow behavior was investigated under isothermal conditions with the inlet volumetric rate being 0.5 LPM, which corresponded in the channels to the Reynolds number of 223. It should be noted that Kumar and Singh named the reverse flow system as a C-type arrangement; however, a more common designation is "U-type".
For each system arrangement, two analytical models were prepared. The first model (BJ, second column in Table 7) employed the momentum coefficients for both headers from the recommended intervals [46], namely, 1.00 in the distributor and 2.65 in the collector. The second model (BJ-Mod) was adjusted to fit the CFD data. The ranges of the particular momentum coefficients had to be widened considering both types of arrangements. The best fit was obtained with the collector momentum coefficient of 3.1 and the distributor momentum coefficients of 0.85 and 0.65 in the U-and Z-systems, respectively. The predicted flow non-uniformity (specified in Table 7) suggests that, in the case of the Z-arranged system, the recommended setup of Model BJ leads to a good agreement with numerical data. Yet, flow rates in peripheral zones differ considerably, as shown in Figure 16a. Model BJ-Mod addresses this issue at the cost of a slightly worse prediction of RSD and NU. More significant is the effect of the modification of Model BJ in the case of the reverse flow system (Figure 16b). The observed maximal relative difference of 4.91% between analytically and numerically predicted flow rates indicates that Model BJ-Mod might be suitable for minichannel applications. A pilot study of possible utilization of analytical models for minigap geometries was conducted using the numerical data reported on the Z-system by Dąbrowski [29]. Dąbrowski's paper covers numerous cases of minigap geometries and inlet flow rates, of which three are discussed here.
Due to the employed analytical model's nature, the gap area (with dimensions 100×99×1 mm) was discretized into 50 branches, and dividing and combining headers were assumed rectangular with the depths of 1, 2, or 8 mm. In all three cases, the inlet mass flow rate was 0.05 kg/s, which corresponded in the discretized branches to the Reynolds number of 1950. For the sake of consistency, the following figures show the obtained A pilot study of possible utilization of analytical models for minigap geometries was conducted using the numerical data reported on the Z-system by Dąbrowski [29]. Dąbrowski's paper covers numerous cases of minigap geometries and inlet flow rates, of which three are discussed here.
Due to the employed analytical model's nature, the gap area (with dimensions 100 × 99 × 1 mm) was discretized into 50 branches, and dividing and combining headers were assumed rectangular with the depths of 1, 2, or 8 mm. In all three cases, the inlet mass flow rate was 0.05 kg/s, which corresponded in the discretized branches to the Reynolds number of 1950. For the sake of consistency, the following figures show the obtained normalized velocities (ratio of fluid velocity in the i-th branch and the average value) and the local maldistribution coefficient (12) together with those provided by the respective study [29].
It is clear that the results yielded by the analytical models can only provide a general distribution trend in the minigap. Figures 17 and 18 display fluctuating CFD data with the associated predictions obtained using the modified analytical models. In spite of the insufficient accuracy of flow prediction for the first system configuration (header depth of 1 mm, the same as the gap itself), results obtained for the third configuration (header depth of 8 mm) fit the CFD data very well as the relative differences of normalized velocity were within the ±5% range. Apart from the rise in flow uniformity, it can be stated that increasing the header depth also improves the accuracy of the modified Model BJ. Despite this improving tendency, no relation between the header depth and the momentum coefficient for distributor flow was observed. The momentum coefficient for collector flow decreased in the case of deeper headers.

Industrial Steam Superheater
As for conventional heat transfer equipment, such as sectional tubular heat exchangers in boilers or steam generators, some maldistribution is always present but usually it is too low to concern the operators. However, this is only true if the equipment is designed properly and operates under stable conditions, as any process fluid flow instability may induce serious problems. It will be demonstrated here that even a minor flow non-uniformity can cause operating problems in an industrial steam superheater. The respective steam superheater was introduced by previous studies [39,40,54]; therefore, only brief information follows:

•
The steam superheater is the first convection heat transfer area in the second stage of a boiler in a chemical plant. The apparatus was in a rather unsatisfactory state due to significant discrepancies between design parameters and the actual operating conditions and also because of unsuitable tube material being used. This led to the ruptures of multiple tubes, intensive fouling, and shutdowns of the Y01 and Y03 channels (named after the "Y"-shaped elements dividing each channel into two tubes) [39].

•
Only a limited amount of information from the operator of the apparatus was available; therefore, a previous work [54] focused on a CFD simulation of the flue gas flow in the boiler. It should be mentioned that it was not possible to carry out physical experiments that would provide data for verification of the computational model.

•
Authors of the present work described in [40] a new methodology of composite modeling that is intended for usage in the initial design of heat transfer equipment. The study [40] also reported drawbacks of modeling simplifications used by commercial software. The distinct advantage of CMS over commercial tools lies in the possibility to easily include data on tube-side fluid flow and, consequently, on more accurate temperature distribution. This can be of vital importance if the equipment is expected to suffer from, e.g., excessive thermal loading due to extreme operating conditions.

Industrial Steam Superheater
As for conventional heat transfer equipment, such as sectional tubular heat exchangers in boilers or steam generators, some maldistribution is always present but usually it is too low to concern the operators. However, this is only true if the equipment is designed properly and operates under stable conditions, as any process fluid flow instability may induce serious problems. It will be demonstrated here that even a minor flow non-uniformity can cause operating problems in an industrial steam superheater. The respective steam superheater was introduced by previous studies [39,40,54]; therefore, only brief information follows: The present study provides the results of thermal analysis of selected tube sections in the superheater that is also shown in Figure 1. Data from general CFD simulations and information on flow distribution in the tubular system obtained using Model BJ provided boundary conditions for a 2D thermal analysis of a single section in the heat exchanger. In addition, CMS results were compared with data obtained using the commercial software HTRI Xchanger Suite [55] to present the current capabilities of the model.
A single tube in the problematic Y03 channel and its equivalent in the Y31 channel were selected for the thermal analysis. The superheater arrangement was unknown; therefore, both U-and Z-arrangements were evaluated. In total, five distribution scenarios were analyzed using CMS. The simplifying assumptions of CMS were as follows:

•
Flue gas flow (derived from velocity and temperature fields obtained using CFD simulations) was discretized into five streams ( Figure 19). • Superheated steam was divided uniformly in the "Y"-shaped element, i.e., every pair of tubes featured the same flow rate (resulting in a staircase-like shape of the distribution graphs in Figure 20).
information on flow distribution in the tubular system obtained using Model BJ provided boundary conditions for a 2D thermal analysis of a single section in the heat exchanger. In addition, CMS results were compared with data obtained using the commercial software HTRI Xchanger Suite [55] to present the current capabilities of the model. A single tube in the problematic Y03 channel and its equivalent in the Y31 channel were selected for the thermal analysis. The superheater arrangement was unknown; therefore, both U-and Z-arrangements were evaluated. In total, five distribution scenarios were analyzed using CMS. The simplifying assumptions of CMS were as follows: • Flue gas flow (derived from velocity and temperature fields obtained using CFD simulations) was discretized into five streams ( Figure 19). • Superheated steam was divided uniformly in the "Y"-shaped element, i.e., every pair of tubes featured the same flow rate (resulting in a staircase-like shape of the distribution graphs in Figure 20). Figure 19. Flue gas flow distribution. Left: velocity field just above the tube bundle obtained via earlier CFD simulations [54]; right: discretized flow rates together with the average flow rate in the respective parts of the shell. Figure 19. Flue gas flow distribution. Left: velocity field just above the tube bundle obtained via earlier CFD simulations [54]; right: discretized flow rates together with the average flow rate in the respective parts of the shell. The observed differences in flue gas flow fields at the inlet were negligible, but a major discrepancy was spotted in the obtained temperature fields in the respective portions of the Y03 and Y31 channels. Obviously, the assessments of tube-side flow non-uniformity produced opposite distribution curves (Figure 20) in the case of U-and Z-arrangements. Hence, five possible situations were investigated as to how the distribution of both fluids influenced the heat transfer prediction. The resulting outlet temperatures of flue gas and superheated steam, heat duty, and estimated thermal effectiveness are summarized in Table 8. The observed differences in flue gas flow fields at the inlet were negligible, but a major discrepancy was spotted in the obtained temperature fields in the respective portions of the Y03 and Y31 channels. Obviously, the assessments of tube-side flow non-uniformity produced opposite distribution curves (Figure 20) in the case of U-and Z-arrangements.
Hence, five possible situations were investigated as to how the distribution of both fluids influenced the heat transfer prediction. The resulting outlet temperatures of flue gas and superheated steam, heat duty, and estimated thermal effectiveness are summarized in Table 8. Although the steam maldistribution was low-the predicted RSD was 0.14% and 1.22% in the U-system and Z-system, respectively-it negatively influenced the outlet temperatures of both fluids. Contrary to expectations, data yielded by CMS shows that the high shell-side flow maldistribution along the tube length (RSD of 17.10%) did not lead to a significant change in the output parameters. Additionally, a combination of higher flow rates of both fluids (U-system Y03) had a calming effect on the yielded outlet temperature maxima. Results of the thermal analysis of the Y03 Z-configuration, which was characterized by the lowest tube flow rate and the highest average temperature of flue gas, were in line with the actual (least suitable) temperature conditions and the reported problems of the Y03 channel. Certainly, the outlet steam temperature in the most exposed tube row was considerably higher, and Table 8 lists only the mean values of the respective parameters. Data obtained for the non-uniform scenarios further indicated that the apparatus's heat duty and effectiveness were affected only slightly by flow maldistribution because the variation of steam flow rates was still at a reasonable level.
Although it was evident that the flow non-uniformity was not the principal cause of the poor condition of the superheater, flow maldistribution might have aggravated other operating problems, especially in the Y03 channel. Unfortunately, it often is impossible to monitor temperatures in every single tube. Therefore, the deterioration of equipment must be slowed first through a proper design that would take flow distribution into account, and then by the utilization of flow improving modifications and regular inspection of the most exposed channels.

Conclusions and Future Work
The study discussed a complex investigation into flow distribution modeling approaches with the main aim of finding a suitable model for further utilization in the composite modelling system that the authors are developing. Two analytical distribution models were examined by using three distinctive geometrical configurations. The conducted experiments (physical and numerical) provided verification data. According to the results, the following conclusions were drawn: 1.
The modified Model BJ achieved the best fit with experimental data in two modeled configurations with a cross-sectional ratio of ≈0.5 and 2.0. In the case of the basic configuration (A r of ≈1.0), the model yielded the second closest match. The maximal relative error in the RSD was 8%, which is an acceptable level of accuracy given the extremely low computational time.

2.
Modifications of Model BJ lay in the adjustment of the momentum coefficient. The best-performing setup was with the momentum coefficient of 1.00 for the distribution system with A r < 1.0, and with the momentum coefficient of 1.10 for the distribution system with A r > 1.0. Flow distribution in the basic system was evaluated most accurately when the momentum coefficient attained the value of 1.078, which is higher than the general recommendation (1.05) of the authors of the original model [46].

3.
As a side effect of the investigation, a notable lack of accuracy was revealed in the maldistribution predictions of steady CFD analyses compared to experimental data. This questions suitability of steady-state data for modeling distribution systems. These inaccuracies, as well as the influence of dynamic flow conditions that were not addressed in this work, will be a part of further research.
Considering the widening of the application area of CMS, the best-performing Model BJ was also examined using mini-scale systems. Data from the literature demonstrated that such equipment suffers from more excessive maldistribution than is usual in conventional applications. The pilot case studies suggested the following:

4.
In minichannel systems, the utilization of modified analytical models is promising, even though the momentum coefficients had to be significantly modified. As for the distributor flow, the best performing momentum coefficient was 0.65 in the most extreme case (Z-arrangement), i.e., almost 40% smaller compared to the mean value recommended for conventional applications. Flow in collectors of both system arrangements was most precisely predicted with the momentum coefficient of 3.10 (approximately 20% larger than recommended). The maximal relative error in the RSD (5%) was observed in the U-system.

5.
Usage of analytical models for modeling of rectangular minigaps is limited. However, Model BJ-Mod was at least able to predict the flow non-uniformity trend (up to 5% discrepancy) if the minigap contained deeper headers. On the other hand, no apparent dependence of the momentum coefficients on geometrical parameters was observed.
The industrial example showed the importance of flow maldistribution and its implementation into design procedures. The analytical method (Model BJ) was utilized to evaluate flow rates in the steam superheater's tubular system. The thermal analysis, carried out using CMS, followed.

6.
It can be concluded that low maldistribution affects the performance of equipment only slightly. Moreover, the obtained data showed a significant discrepancy between temperature fields if uniform and non-uniform flow of both process fluids was assumed. The results of the thermal analysis also agreed with the reported operating issues of the respective apparatus.
Future work on CMS, which utilizes the distribution model discussed in this study, will focus on two main aspects. Firstly, the usage of analytical models in mini-scale applications requires thorough validation. Secondly, major thermal computation modifications are necessary to get a versatile tool for predicting and evaluating non-uniform heat transfer in conventional and mini-scale equipment.