Model Validation for the Heat Transfer in Gasket Plate Heat Exchangers Working with Vegetable Oils

: Many models for accurately predicting the performance of gasket plate heat exchangers were developed in the last decades, grouped in three categories: empirical, semi-analytical or theoretical/numerical, with the view to saving materials and energy through correct design of industrial equipment. This work addresses one such model, namely L é v ê que correlation modiﬁed by Martin and by Dovi´c, which is promising due to the correct assumption of the ﬂow in sine duct channels and the consideration of energy losses caused by ﬂow reversal at plate edges and the ﬂow path changing when entering the chevron angle. This model was validated by our own experimental data under industrial conditions for vegetable oils processing, both in laminar ﬂow ( Re = 8–42) and fully developed turbulent ﬂow ( Re = 446–1137). Moreover, in this study, particular values for constants/parameters of the model were determined for the corrugation inclination angle relative to vertical direction equal to 30 ◦ . Through statistical analysis, this study demonstrates that this particularized form of the generalized L é v ê que correlation can be used with conﬁdence.


Introduction
In technological installations, heat exchangers are parts of the process flow; they ensure the heating, cooling or heat recovery from hot fluxes, as the technology requires. The use of plate heat exchangers (PHE) is desirable since they have a simple compact construction, are reliable, and very efficient at small temperature difference between fluids. Additionally, they have low fabrication and maintenance costs. This type of heat exchanger is a good choice in vegetable oils refining where the high viscosity of the fluids represents a high resistance to flowing and heat transfer.
Gasket plate heat exchangers were introduced in vegetable oil processing due to their specific geometry [1], with small hydraulic diameter and corrugations, therefore inducing enhanced Reynolds numbers and enlarging the heat transfer area, resulting in the intensification of heat transfer [2].
Numerous studies have been carried out to explore the influence of geometric characteristics on the heat transfer [3][4][5][6][7][8]. The results show that geometric parameters such as chevron/herringbone angle, surface enlargement factor, and the corrugation profile are important factors influencing the thermal-hydraulic behavior. Another influencing factor is the flow maldistribution [9][10][11]. S. Al-Zahrani, et al. [12] showed that the maldistribution of the fluids in ports contributes to the performance of the apparatus and their position should be rethought by introducing a longitudinal baffle, so that the cross sectional area changes radically. The simulations proved that friction factors for the modified PHE are 4.5-7 fold higher than those of conventional corrugate plate apparatus with effect on the overall thermal performance.
The experiments on fluid hydrodynamics and thermal transfer in this type of heat exchanger allowed for the development of correlations in accordance with the similitude

Equipment and Data Collection Procedure
For the development of a mathematical model describing the heat transfer in chevron plate heat exchangers, four gasket chevron plate exchangers (HE) were tested, with characteristics described in Table 1 and Figure 1. The flow of both fluids in the pre-heater HE#1 is in counter-current with two passes, and in the other three HE is one-pass counter-current.  The advantage to use large scale equipment for modeling is that models proceeding from laboratory scale experiments can be validated on it, on one hand, and newly devel- The advantage to use large scale equipment for modeling is that models proceeding from laboratory scale experiments can be validated on it, on one hand, and newly developed models are more accurate, on the other hand.
Data collection was performed in a vegetable oil processing installation. It is a fully automated process. Temperatures are controlled at constant values by the cooling water flowrate in case of HE#2, HE#3, HE#4 and by the bleached oil flowrate in case of HE#1. Measured temperatures are shown in the third column of Table 2. The average temperature at the wall is calculated by considering no temperature gradient in it, so the temperature should be the same on both sides of the wall. In each campaign, the measurements were performed during the steady-state regime which lasts for hours and days. In the case of first two campaigns, there were four such regimes when flowrates and temperatures remained constant. In the third campaign, the installation worked in only one steady-state regime.

Fluids
It is of great importance to determine the values of the physical properties of fluids with precision, either experimentally or by calculating with reliable mathematical correlations. Therefore, we determined, experimentally, the variation with temperature of density and dynamic viscosity of all three types of oils (raw, bleached, and winterised), covering the whole range in functioning (20-110 • C). The analyses were made with the apparatus Anton Paar, model SVM 3000 (Anton Paar OptoTech GmbH, Seelze-Letter, Germany), using the standard method ASTM D445/ISO 121852.
There were three processing campaigns in the factory, working with three crudes as follows: sunflower oil 1 (SO1), sunflower oil 2 (SO2) and rapeseed oil (RO).
By processing these data with the least squares' method, the correlations allowing the calculation of oil density and viscosity at the working temperature were obtained. These correlations are given in Tables 3 and 4. Even though the equations parameters seem very close for different type of oils, small differences can lead to bigger differences in similitude criteria values, especially for Prandtl numbers, with an effect on the accuracy of the mathematical model predicting the global heat transfer coefficients.
The specific heat capacity (c p ) and the thermal conductivity (λ) were estimated by the Equations (1) and (2) of Choi and Okos [24], developed for the vegetable oils: (1) For cooling water, the physical properties at different temperatures (density, viscosity, specific heat capacity and thermal conductivity) are found in reference [25].

Modeling
The predominant heat transfer mechanism involved in this type of heat exchanger is the forced convection in stationary regime. Most correlations that are used for modeling the heat transfer describe the phenomenon in terms of similitude criteria, following the similitude theory, and the most frequent correlation is exponential (Equation (3)): where: Nu, Re and Pr are three similitude criteria: Nusselt, Reynolds and Prandtl and coefficients a, b and c are to be determined for different flow regimes, for fluids with and without change in phase, for different geometries. Reynolds number is calculated in relation to mass velocity in a channel G ch [kg m −2 s], hydraulic diameter d h and dynamic viscosity µ (Equation (4)).
Prandtl number is calculated from physical properties of the fluids with Equation (5).
Nu number contains the heat transfer coefficient α [W m −2 • C −1 ] which can be determined from Equation (6), after calculating Nu with Equation (3): For the forced convection in stationary regime, the coefficients a, b, c in Equation (3) were determined in many experimental studies but, the Kumar correlation (Equation (7)) is mostly recommended [18,26] to calculate Nu number when starting the development of a new and more accurate model.
where µ w is the dynamic viscosity of the fluid at the wall temperature. The Equation (7) is valid when referring to the cross-sectional free area. Lévêque gave another solution for the heat transfer during the laminar flow in the gap between two plates of length L and constant spacing H (Equation (8)), the so-called "1/3 power law", confirmed by the numerical results of Holzbecher's simulations [27]: The original Equation (8) was modified for micro-channels with sine ducts, after Martin [21] showed that the heat coefficient depends on pressure drop which is directly proportional to f app × Re 2 , where f app is the apparent friction factor (Equation (9)): Later, Dović and co-authors [22] adjusted this equation, taking into account the influence of factors ignored in the Equation (9) (swirling, real effective flow length), and in accordance with their own experimental data corroborated by other authors' data [28,29] (Equation (10)): In plate heat exchangers with chevron corrugations, the high heat transfer is due to the flow in small hydraulic diameter sections and interractions between streams flowing and generating swirl secondary flows. Dović and co-authors [22] describe the flow pattern, following an experiment with aqueous glycerol (62% wt.), with Re in range 0.1 to 250.
This experiment allowed us to visualize two flow components: the furrow component (Equation (11)) and the longitudinal component (Equation (12)), and their partial mixing.
where l is the corrugation wavelength.
In the classical correlation (Equation (3)), Nu, Re and average velocity u [m/s] are related to the total chevron channel cross sectional area transverse to the main flow direction. Taking into account the sinusoidal flow, by means of analytical approach, the average velocity in the cell's sine duct in furrow direction u sine is recalculated with Equation (13): where . m ch is the mass flowrate in the channel [kg/s], and A ch, sine is the channel crosssection transverse to the furrow (Equation (14)): L w and β are geometrical characteristics from Table 1. Then, Reynolds sine is calculated with Equation (15).
where d h,sine is the the hydraulic diametre of a sine duct; it is related to the independent variable x-the ratio corrugation depth: corrugation wavelength (x = b/l), in Equation (16) [30]: Nu sine referred to the cell's sine duct is correlated with Nu related to the the whole channel through Equation (17): Nu Nu sine = dh dh sine (17) L furr in Equation (11) is the effective cell length L cell , in case of β < 60 • . In case of β > 60 • , L cell = L long , following the prevailing pattern flow; f app is the apparent friction coefficient which, for a given geometry, has the generalized form (Equation (18)): Constants B (Equation (19)) and C (Equation (20)) are functions of channel geometry, and can be calculated as a function of the aspect ratio (x = b/l), by polynomial functions Equations (19)-(23), developed by Dović and co-authors [22] based on the work of Shah [30]: C = 2.6624x 4 − 10.586x 3 + 11.262x 2 − 0.1036x + 9.6 (20) K e (∞) = −5.888x 4 + 9.46113x 3 − 4.248x 2 − 0.1333x + 2.648 (21) where K(∞) represents the incremental pressure drop number, K e -kinetic energy correction factor, and K d is the momentum flux correction factor. This model (Equations (10)- (23)) is intended to be validated by our own experimental data. The challenge is to find generalized correlations fitting to all flow regimes (laminar, intermediate and turbulent).

Results and Discussion
As mentioned before, there were three processing campaigns of vegetable oils. In two of them, the sunflower oils I and II were processed at four different flowrates, and in the third, the rapeseed oil was processed at one flowrate. In total, 72 sets of experimental data resulted. In Table 5, the total mass flowrate refers to the whole flowrate entering the circular port of the heat exchangers, since Re is calculated taking into account the hydraulic diameter of the channel, d h (see also Table 1). The Nu number is calculated by applying the Kumar correlation (Equation (7)) from [18]. As seen in Table 5, the experimental data covers a range of laminar flow regime (Re = 18-90, Re < 100) and another range of fully developed turbulent regime (Re = 921-2391, Re > 250). Re and Nu numbers, related to the chevron channel cross sectional area transverse to the main flow direction, served later to calculate the experimental values of Re sine and Nu sine , related to cell's sine duct in furrow direction.
Both Equations (9) and (10) consider that heat transfer depends on pressure drop, which is directly proportional to product fRe 2 . Additionally, u sine , dh sine , f app , L furr , and constants B and C were determined in each experimental point, for the Nu sine,calc estimation, according to the model of Lévêque adjusted by Dović et al. [22].
The constants B and C were calculated with Equations (18)-(23), for each heat exchanger and the average values were B = 0.19952 ± 0.00006 and C = 12.4239 ± 0.0121. These values are for the first time calculated in case of gasket plate heat exchangers with corrugation angles 30 • and the aspect ratio b/l = 0.8, since the model to be validated recommends application for b/l < 0.5.
The friction factor, f app , determined in each experimental point with Equation (17), was plotted vs. Re sine (Figure 2) and compared with the model, on the whole range of Re sine = 8-1137. The model fits almost perfectly at Re sine > 25 with experimental data of Muley et al. [13] and Heavner [31] but has errors of 5-15% at Re sine = 9-25. These errors can be explained by the different geometry (the aspect ratio b/l > 0.5 in contrast with b/l < 0.5 for all previous experiments), which can affect the friction coefficients at Re sine < 25.  Table 6. This table also includes the relative errors between the experimental values of Nu sin,exp (Equation (14)) and those computed with adjusted correlation Lévêque (Equation (10)).
For the whole set of data, the relative errors are between −18.8% and +28.5%, both positive and negative values, demonstrating that they are not systematical errors. The medium error is |ε| = 9.56%.
The plot Nu sine,calc vs. Nu sine,exp (Figure 3), showing the dispersal of data around the bisector, uncovers that points spreading is minimal at lower values of Nu sine , corresponding to the laminar flow, since at values corresponding to the turbulent flow, Nu sine are spread at a bigger distance; however, dispersal on both sides of the bisector for these values too, indicates that the model can be used up to Re sine = 1137, with confidence.  For the forced convection in stationary regime, it is customary to look for corelations between Re and Nu/Pr 1/3 (μ/μw) 0.14 . We used the facility Data analysis in Excel. The computed results for the similitude criteria taking into account the sine duct flow are shown in Table 7. Table 7. Statistical data analysis for Resine vs. Nusine /Pr 1/3 (μ/μw) 0.14 correlation.

SUMMARY OUTPUT
Regression Statistics Multiple R 0.9924 For the forced convection in stationary regime, it is customary to look for corelations between Re and Nu/Pr 1/3 (µ/µ w ) 0.14 . We used the facility Data analysis in Excel. The computed results for the similitude criteria taking into account the sine duct flow are shown in Table 7. Table 7. Statistical data analysis for Re sine vs. Nu sine /Pr 1/3 (µ/µ w ) 0.14 correlation. As seen in Table 7, both regression statistics and ANOVA analysis show that Re sine and Nu sine /Pr 1/3 (µ/µ w ) 0.14 are strongly correlated through linearity. The coefficient of determination is close to 1 (R 2 = 0.9849). The ANOVA variance test demonstrated that null hypothesis is false, so the mean squares for treatments MST is much larger than the mean squares for errors MSE (3770.1214 0.8269); the MST/MSE ratio is F = 4559.34, a very high value, far from F = 1, the value for non-correlated series of data. In conclusion, the Fisher test confirms the strong correlation between the series of data. The linear correlation (Equation (24)):

Regression
Nu sine /Pr × (µ/µ w ) 0.14 = 1.3305 + 0.0195 Re sine The p-values < 0.05 for the intercept and the slope of the line demonstrate that the model works with the confidence in the interval ±95%. Nu sine includes the partial heat transfer coefficient α in the fluid on one side of the plate; the coefficient on the other side of the plate can be calculated with the corresponding value of Nu sine . Taking also into account the conductivity of the plate material, the overall heat transfer coefficient can be found. It is very important to have an accurate model for the estimation of Nu sine , because the errors spread in the values of the overall coefficient.

Conclusions
In this study, the semi-analytical model of Lévêque for heat transfer in gasket plate heat exchanger channels was validated in our experimental data obtained from industrial equipment for vegetable oil processing. The vegetable oils had different physical properties proceeding from three different raw materials in different stages of the process (raw, bleached, winterized). There were four exchangers with similar geometry but different in size. The flow regime was laminar (Re < 100) for oils and turbulent (Re > 250) for water, allowing us to check the model in a large range of Re numbers, related to the cross-sectional free area of the channel.
The model of Lévêque takes into consideration the flow in cell's sine duct in furrow direction. Consequently, the similitude criteria Re sine and Nu sine are calculated in new considered conditions. They are correlated with a more complex relationship, considering the construction details of the channels. In this work, the constants C and D of the model were determined, for corrugation inclination angle relative to vertical direction equal to 30 • , for the first time. Thus, B = 0.19952 ± 0.00006 and C = 12.4239 ± 0.0121. As regards the constant C 1 of the model, that is the same as determined by Dović et al. [22] when adapting Lévêque correlation to the corrugation inclination angle 60 • (C 1 = 0.1534).
Both the analysis of relative errors and the statistical analysis concluded that Lévêque correlation adapted by Dović et al., with model parameters, determined here for corrugation angle 30 • , can be used with confidence for predicting Nu criterium as a basis for heat transfer calculations in gasket plate exchangers.
Author Contributions: Conceptualization, methodology, software, validation, investigation, A.-A.N. and C.I.K.; project administration, data curation and writing-original draft preparation, A.-A.N.; supervision, writing-review and editing, C.I.K. All authors have read and agreed to the published version of the manuscript.
Funding: No external fundings were received for this study.
Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.
Data Availability Statement: Primary data were obtained from our original experiment. All processed data are included in this article and are available for further processing and interpretation by other authors.

Conflicts of Interest:
The authors declare no conflict of interest.