Optimizing heat exchanger efficiency: Predictive modeling to minimize fouling in crude oil refining

Abstract

Fouling in heat exchangers significantly compromises energy efficiency in crude oil refining, leading to increased operational costs and environmental impacts. This study presents a predictive model aimed at enhancing heat exchanger performance by minimizing fouling resistance. Model fitting was conducted using approaches of varying complexity, with measures taken to avoid overfitting. The models were subsequently refined to incorporate key variables, such as inlet and outlet temperatures and mass flow rates, ensuring robustness and generalizability. The final simplified model comprises only 19 terms, yet achieved high predictive performance (R² = 0.961; predicted R² = 0.956) and effectively addressed multicollinearity. The selected model identified significant linear, quadratic, and interaction effects among thermal and flow parameters, with the mass flow rates of the hot fluid (MFH) and cold fluid (MFC) emerging as particularly influential. Notably, the model demonstrated that fouling resistance decreases substantially with increasing hot fluid flow rate. Optimization using a desirability function identified 37 parameter combinations that achieved a fouling resistance (RFC) of 0.001 m²°C/W with a maximum desirability score of 1.00, consistently favoring high MFH values (~93.00 kg/s) and hot fluid outlet temperatures (THO) near 43.00°C. These findings confirm the model’s robustness and practical applicability, providing actionable insights for operational strategies aimed at minimizing fouling while maintaining thermal efficiency.

Keywords

heat exchanger efficiency model optimization energy efficiency predictive modeling flow rate optimization multicollinearity fouling resistance minimization

Introduction

Amid escalating energy prices and the gradual depletion of conventional energy resources, enhancing energy efficiency has become a strategic imperative across industrial sectors. This urgency has driven significant advancements in the design and optimization of energy systems, particularly those involving thermal management. The development and application of heat exchangers—key components in thermal management systems—span a broad range of industrial contexts, underscoring their essential role.¹ Heat exchangers facilitate the transfer of thermal energy between fluids without mixing, thereby directly influencing the overall efficiency of numerous processes.² To improve heat transfer efficiency, optimize energy recovery, and enhance cost-effectiveness, researchers and manufacturers have explored various strategies, including optimal design, precise control of operating conditions, and routine maintenance.³ Hameed et al.⁴ introduced a novel oscillation technique for multi-tube heat exchangers, demonstrating a fivefold improvement in heat transfer and enhanced thermal performance at specific Reynolds numbers, thus emphasizing the potential of passive technologies for efficiency enhancement.

Despite their advantages, heat exchangers are susceptible to fouling—the accumulation of solid matter on their surfaces—which significantly reduces their operational efficiency. Fouling impairs heat transfer, restricts fluid flow, promotes corrosion, and contaminates working fluids.^5,6 These effects result in increased equipment requirements, production losses due to fouling-related downtime, elevated costs for deposit removal, and greater consumption of fuel, water, and electricity across various industrial sectors.⁷ The petrochemical industry is among the most severely affected. Deposits can be identified through filtration and chemical analysis, and subsequently addressed via chemical treatments.⁸ Following fouling detection, systematic cleaning is required. Short cleaning intervals lead to increased production downtime, whereas extended intervals exacerbate energy consumption and environmental impacts.⁹ These challenges highlight the growing importance of predictive maintenance over routine, systematic maintenance.

In response, researchers are increasingly prioritizing the prediction of fouling, prompting the development of methods to enhance predictive accuracy.^10,11 Traditional approaches—such as experimental investigations and computational fluid dynamics (CFD) simulations—often face limitations in terms of time efficiency and predictive accuracy when applied to fouling analysis.¹² Consequently, research efforts have increasingly focused on the use of statistical modeling algorithms. For instance. Elghool et al.¹³ employed multi-objective optimization for heat pipe–heat sink (HP–HS) systems in thermo-electric generators (TEGs), achieving a 36.7% improvement in thermal efficiency and a 17.9% reduction in cost, thereby demonstrating the potential of integrated design approaches in predictive modeling.

Design of Experiments (DOE) methods are statistical techniques employed to investigate the nonlinear behavior of complex heat transfer systems. Sahin¹⁴ applied the Taguchi method to evaluate the influence of design parameters—such as Reynolds number, fin height, and fin pitch—on the performance of heat exchangers equipped with circular fins. Performance was assessed using the Nusselt number and friction expanded this line of research by analyzing perforated fin configurations under forced convection, demonstrating that square perforations reduced thermal resistance by 16°C and enhanced turbulence, while circular perforations increased heat transfer by up to 51.29% compared to solid fins. Wang et al.¹⁵ addressed triangular leakage zones between tubular heat exchanger baffles by employing response surface methodology (RSM) with a central composite design to study the effects of skew angle, overlap, and inlet velocity on the Nusselt number and pressure drop. Strong correlations were observed, with coefficients of determination of 0.943 for the Nusselt number and 0.999 for the friction factor. Maddah et al.¹⁶ investigated the effects of nanofluid concentration, Reynolds number, and twist ratio on the thermal performance of a double-tube heat exchanger using an Al₂O₃–TiO₂ hybrid nanofluid. Al Kumait et al.¹⁷ further explored nanofluid applications, demonstrating that TiO₂/water nanofluids in helical-ribbed tubes enhanced turbulent kinetic energy and Nusselt numbers by approximately 2%, with minimal deviation in friction factor models.

The most critical efficiency characteristics were identified using a full factorial design combined with statistical analyses, including t-tests, ANOVA, and F-tests. Parameters such as nanoparticle concentration, Reynolds number, and swirl ratio were shown to enhance heat exchanger efficiency. Wang et al.¹⁸ investigated heat transfer enhancement in drag-reducing surfactant fluids using photo-rheological counterions, identifying optimal surfactant/counterion concentrations (4 mM/5 mM) that improved thermal performance. Bisognin et al.¹⁹ employed a fractional factorial 2⁵⁻¹ design to evaluate both the individual and interactive effects of particle size, gas velocity, particle thermal conductivity, tube diameter, and tube spacing on the heat transfer coefficient in a fluidized bed heat exchanger with a horizontal tube arrangement. Predictive models were developed based on 16 experimental trials to assess the statistical significance of each coefficient. Chowdhury and Borah²⁰ analyzed how the heat transfer rate, number of transfer units (NTUs), and heat exchanger efficiency are influenced by the inlet temperature of the hot fluid and the mass flow rates of both hot and cold fluids. These performance characteristics were modeled, analyzed, and optimized using a Box–Behnken design and response surface methodology. The resulting regression models for heat transfer rate, NTU, and efficiency showed strong agreement with experimental data, with R² and adjusted R² values exceeding 90%. The study concluded that optimal performance was achieved with hot and cold fluid mass flow rates of 50 and 200 L/h, respectively, and a hot fluid inlet temperature of 80°C. Jradi et al.²¹ utilized a central composite design to assess the thermal performance of a cross-flow heat exchanger. Input variables included acid inlet and outlet temperatures, vapor temperature, acid density, volumetric flow rate, and time, while fouling resistance served as the response variable. Both the main effects and interactions of operational parameters significantly influenced fouling resistance. The statistical indicators R², adjusted R², and predicted R² approached 1.0, indicating a high degree of model accuracy, and the regression curves demonstrated strong fit. Finally, Hallaji et al.²² conducted a comprehensive Taguchi experimental design study to evaluate the reduction in heat transfer and fouling rate of aqueous CaSO₄ solutions. They investigated the effects of flow rate, CaSO₄ concentration, heat flux, and bulk temperature under conditions of forced convection and supercooled flow boiling. ANOVA results confirmed strong agreement between experimental outcomes and model predictions.

The primary objective of this study is to develop and validate an optimized predictive model for heat exchanger fouling in crude oil refining, utilizing factorial design and multistep regression analysis. This model aims to enable refineries to minimize energy losses, reduce environmental impacts, and optimize maintenance schedules while maintaining operational efficiency. The novelty of this work lies in a comprehensive approach that balances model complexity and predictive accuracy through the systematic reduction of higher-order polynomial terms, while preserving essential interaction effects. In contrast to earlier studies that prioritized maximizing R² through complex polynomial expansions—often at the expense of severe multicollinearity—this study explicitly addresses overfitting and multicollinearity by employing term elimination and variance inflation factor (VIF) analysis, thereby resolving a critical gap in fouling prediction methodologies.

Methodology

Box–Behnken method

The Box–Behnken design (BBD) is a robust and efficient technique used in response surface methodology (RSM) to optimize experimental processes. Developed by George E. P. Box and Donald W. Behnken in 1960, this design is particularly advantageous because it reduces the number of experimental runs while still providing high-quality predictions for the response variable within the experimental region.²³ BBD is widely applied in various fields, including engineering, chemistry, and biotechnology, due to its ability to efficiently investigate quadratic response surfaces without requiring a full three-level factorial experiment.²⁴

The relationship between the model response ( $Y$ ) and the process parameters (X_i) is defined by a general-order polynomial regression equation:

Y = β_{0} + \sum_{i = 1}^{k} β_{i} X_{i} + \sum_{i = 1}^{k - 1} \sum_{j = 2}^{k} β_{ij} X_{i} X_{j} + \sum_{i = 1}^{k} β_{ii} X_{i}^{n} + ε

(1)

where $β_{0}$ is the intercept term, $β_{i}$ are the linear coefficients, $β_{ii}$ represent the quadratic, cubic, or quartic coefficients, $β_{ij}$ are the interaction coefficients, $X_{i}$ denotes the coded value of the i-th independent variable, k is the number of independent variables, and $ε$ is the error term.

The Box–Behnken design was employed to evaluate the experimental variables, and Design-Expert software was utilized for regression and graphical analysis. A multicriteria optimization approach based on the desirability function was also applied.²⁵ The objective of the optimization procedure was to minimize fouling resistance to 0.001 m²°C/W while maximizing heat exchanger efficiency, without compromising other performance criteria.

Evaluation of model effectiveness

To effectively validate a model’s performance, it is essential to employ a range of statistical metrics that assess various aspects, including explanatory power, predictive accuracy, and robustness.^26–30 Metrics such as the coefficient of determination (R²), adjusted R², and predicted R² evaluate how well the model accounts for variance in the data, thereby offering insights into its explanatory capabilities. Additionally, advanced metrics such as adequate precision, F-statistics, and p-values help assess the model’s reliability and statistical significance. To further examine the model’s fit and complexity, criteria such as the prediction error sum of squares (PRESS), −2 log-likelihood, Bayesian information criterion (BIC), and corrected Akaike information criterion (AICc) are critical. Collectively, these metrics provide a comprehensive evaluation of the model’s performance and its overall adequacy in representing the data.^31–35

Experimental procedure and data collection

The atmospheric distillation unit U100 serves as the primary processing unit at a refinery in Algiers. It fractionates crude oil into various end products, including kerosene, diesel, fuel, liquefied gas, and light and heavy solvents. These products may be marketed directly or subjected to further treatment processes. One of the three centrifugal pumps, P101, transfers crude oil from storage tanks at ambient temperature to the atmospheric distillation unit. The oil then flows through two circuits of the E101 heat exchanger battery (CBA and FED).³⁶ On the tubular side of the battery, the oil is preheated via overhead reflux (RT), consisting of a mixture of light products collected from the top of the C101 distillation column at tray No. 46. Subsequently, the oil passes through an electrostatic desalter, where treated water and caustic soda are added. This mixture is injected both at the inlet of the heat exchanger E101 and at the inlet of the desalter to wash the crude oil and remove the salts present (Figure 1). The inlet and outlet temperatures of the two fluids are measured at the ends of the heat exchanger using four thermocouples. Simultaneously, the flow rates of the crude oil and the return flow at the heat exchanger inlet are monitored using flow meters. The physical properties of the two fluids are provided by the refinery’s control room.

Figure 1.

Simplified representation of the crude oil preheating circuit.

The present study utilized data collected over 290 days, from March 14, 2019, to December 17, 2019, from the E101 CBA heat exchanger cell, located within the preheating circuit of the Algiers refinery. This cell comprises three counterflow shell-and-tube heat exchangers connected in series. The characteristics of these heat exchangers are presented in Table 1. The operating variable ranges corresponding to the heat exchanger used in this study are provided in Table 2.

Table 1.

Characteristics of the E101 CBA heat exchanger.

Parameter	Value
Construction material	Carbon steel
Shell diameter (m)	1.067
Baffle spacing (m)	0.465
Number of shells	3.000
Tube outer diameter (m)	0.020
Tube specification	BGW14
Tube length (m)	5.740
Total number of tubes	6600
Pitch: staggered (m)	0.025
Total heat exchange surface area (m²)	2322.770
Overall heat transfer coefficient (kW/m²°C)	36.680
Operation conditions	Side shell	Tube side
Circulating fluid	Head reflux	Crude oil
Mass flow rates (kg/s)	126.000	90.120
Viscosity (m²/s) inlet/outlet	—	2.4 × 10⁻⁶–9.6 × 10⁻⁷
Inlet temperature (°C)	115.556	26.667
Outlet temperature (°C)	65.56	104.44
Number of passes	1.000	4.000
Fouling factor	0.001	0.002

Table 2.

Parametric ranges of the operating variables.

Variable	Symbol	Abbreviation	Range
Inlet temperature of crude oil (°C)	t _e	TCI	31–17
Outlet temperature of crude oil (°C)	t _s	TCO	92–110
Inlet temperature of head reflux (°C)	T _e	THI	111–130
Outlet temperature of head reflux (°C)	T _s	THO	44–64
Mass flow rate of crude oil (kg/s)	${\overset{\cdot}{m}}_{t}$	MFC	23.50–46.09
Mass flow rate of head reflux (kg/s)	${\overset{\cdot}{m}}_{c}$	MFH	38.98–80.10

The experimental calculations of the overall heat transfer coefficient (U_t) and fouling resistance ( $R_{f}$ ) are based on three simplifying assumptions: (i) both fluid flows (crude oil and head reflux) are arranged in a countercurrent configuration; (ii) thermal losses are neglected; and (iii) fouling occurs only on the crude oil side. Although the specific heat of a fluid typically varies with temperature, it can be considered constant within a defined temperature range by assuming an average value, with minimal loss of accuracy. Additionally, the changes in velocity and elevation of the working fluid streams are negligible, thus minimizing the impact on kinetic and potential energy terms. The outer surface of the tested heat exchanger is fully insulated, preventing heat loss to the surroundings and ensuring that heat transfer occurs exclusively between the two working fluids. For the control volume enclosing the heat exchanger, as depicted in Figure 1, the only work performed is flow work at the inlet and outlet boundaries. Therefore, the term W_cv is eliminated, reducing the energy rate balance to:

\begin{matrix} \underset{0}{\underset{︸}{\frac{d E_{cv}}{dt}}} = {\overset{\cdot}{Q}}_{cv} - \underset{0}{\underset{︸}{{\overset{\cdot}{W}}_{cv}}} + \sum_{i} {[{\overset{\cdot}{m}}_{i} (h_{i} + \underset{0}{\underset{︸}{\frac{{V_{i}}^{2}}{2} + g z_{i})}}]}_{i} \\ - \sum_{e} {[{\overset{\cdot}{m}}_{e} (h_{e} + \underset{0}{\underset{︸}{\frac{{V_{e}}^{2}}{2} + g z_{e})}}]}_{e} = 0 \end{matrix}

(2)

For a quasi-isobaric process, the heat transfer rate from the head reflux to the crude oil can be expressed via an energy balance as:

{\overset{\cdot}{Q}}_{r} = {\overset{\cdot}{m}}_{C} C_{pr} (T_{s} - T_{e})

(3)

The overall heat transfer coefficient in the presence of fouling is defined by:

U_{s} = \frac{Q}{A . F . Δ T_{TM}}

(4)

The logarithmic mean temperature difference (ΔT_LM) for a system of three counterflow shell-and-tube heat exchangers is given by:

Δ T_{LM} = \frac{(T_{s} - t_{e}) - (T_{e} - t_{s})}{\ln \frac{(T_{s} - t_{e})}{(T_{e} - t_{s})}}

(5)

The heat exchanger under investigation undergoes mechanical cleaning during transitions between operational cycles. Consequently, each new process begins without any fouling. The initial heat transfer coefficient at the start of each process is regarded as the clean design value (U_p).

The fouling resistance over time (RFC) is calculating using the following equation:

RFC = \frac{1}{U_{s}} - \frac{1}{U_{p}}

(6)

To select the optimal operating parameters, the thermal efficiency (ε) is determined by the ratio of the actual heat transfer rate to the maximum possible heat transfer rate, as given by:

ε = \frac{\overset{\cdot}{Q}}{{\overset{\cdot}{Q}}_{\max}}

(7)

The actual heat transfer rate is expressed as:

\overset{\cdot}{Q} = {\overset{\cdot}{m}}_{c} C_{p, c} (T_{e} - T_{s}) = {\overset{\cdot}{m}}_{t} c_{p, t} (t_{s} - t_{e})

(8)

The maximum possible heat transfer rate is calculated by:

{\overset{\cdot}{Q}}_{\max} = C_{p, \min} (T_{e} - t_{e})

(9)

To determine the maximum possible heat transfer rate $({\overset{\cdot}{Q}}_{\max})$ in the tested heat exchanger, the initial step is to verify whether the maximum temperature difference occurs between the inlet temperatures of the hot and cold fluids.

Table 3 presents the thermal and flow statistics of the heat exchanger used for crude oil processing. The tubing section includes four temperature variables (A–D). The crude oil inlet temperature (TCI) ranged from 13.00°C to 31.00°C, with an average of 24.25°C and a standard deviation of 4.08°C. This range and variability reflect fluctuations in the incoming crude oil temperature. The crude oil outlet temperature (TCO) varied between 92.00°C and 119.00°C, with an average of 104.16°C, indicating an average thermal increase of approximately 80°C due to heat exchange. The average inlet temperature of the head reflux (THI) was 121.08°C, ranging from 108.42°C to 136.00°C, while the average outlet temperature (THO) was 63.39°C, with a range of 42.00°C–81.00°C. The average temperature drop of 57.69°C in the head reflux highlights its critical role in thermal energy transfer. The shell-side mass flow rates of head reflux and crude oil are represented by MFH and MFC, respectively. The crude oil mass flow rate (MFC) ranged from 23.50 to 92.72 kg/s, with an average of 40.60 kg/s and a standard deviation of 10.46 kg/s. The head reflux mass flow rate (MFH) averaged 54.20 kg/s, with a range of 23.22–91.42 kg/s. Figure 2 illustrates the distribution and behavior of fouling resistance (RFC) in the crude oil refining system. The distribution is right-skewed, with a higher concentration of observations at lower RFC values. The highest frequency was observed in the 0.25 m²°C/W bin, accounting for 46.21% of the data. This indicates that fouling resistance is low nearly half the time, suggesting efficient operational performance.

Table 3.

Statistical analysis of heat exchanger parameters in crude oil processing.

Side	Factor	Name	Minimum	Maximum	Mean	Std. dev.
Tube	A	TCI	13.00	31.00	24.25	4.08
	B	TCO	92.00	119.00	104.16	3.86
	C	THI	108.42	136.00	121.08	3.07
	D	THO	42.00	81.00	63.39	7.40
Shell	E	MFC	23.50	92.72	40.60	10.46
	F	MFH	23.22	91.42	54.20	10.10

Figure 2.

Relative frequency distributions for various RFC categories.

As RFC values increase, their frequency decreases. The 0.75 and 1.25 m²°C/W bins account for 16.77% and 15.16% of the observations, respectively. Only 6.45% of cases exhibited RFC values greater than 2.25 m²°C/W, indicating that high levels of fouling are relatively rare. No occurrences were recorded at 4.25 m²°C/W, and the highest observed bin (4.75 m²°C/W) showed a frequency of just 0.12%, underscoring the infrequency of severe fouling events.

This distribution suggests that the system generally operates under low to moderate fouling conditions. However, periods of elevated fouling resistance may adversely affect performance and increase maintenance requirements.

Additionally, the correlation matrix presented in Figure 3 illustrates the relationships between fouling resistance (RFC) and key process parameters, thereby revealing variables associated with fouling behavior. RFC exhibits the strongest positive correlation with the crude oil inlet temperature (TCI), with a coefficient of 0.502, indicating that higher TCI values tend to increase fouling resistance. This observation is consistent with the variability of TCI reported in Table 3. Conversely, RFC shows a moderate negative correlation with the crude oil mass flow rate (MFC; r = −0.484), suggesting that higher flow rates may reduce fouling due to enhanced heat transfer and decreased deposition. Weaker negative correlations were also observed between RFC and the crude oil outlet temperature (TCO; r = −0.172), head reflux inlet temperature (THI; r = −0.244), and head reflux outlet temperature (THO; r = −0.120). A modest positive correlation was found between RFC and the head reflux mass flow rate (MFH; r = 0.110), indicating a limited influence of this variable on fouling resistance.

Figure 3.

Input-output variable correlations.

Results

The model analysis begins with insights from Tables 4 and 5, which evaluate model performance across varying levels of complexity. The results indicate that model fit improves with increasing complexity; however, higher-order models also carry a greater risk of overfitting. Table 4 presents the fit analysis of the reduced quartic model, which was selected by Design-Expert as the optimal model. As the model complexity increases—from linear to fifth-order—the standard deviation decreases from 0.5556 to 0.0761, and model fit metrics improve accordingly. For example, the fifth-order model achieved an R² of 0.998, compared to 0.606 for the linear model, indicating a much better explanation of data variance. However, this improvement comes at the potential cost of model simplicity and generalizability. In the quartic model, the adjusted R² and predicted R² reached 0.9928 and 0.9549, respectively, suggesting a strong balance between fit and predictive ability. The PRESS statistic, which assesses predictive performance, decreased substantially—from 202.86 in the linear model to 8.41 in the cubic model—but rose slightly to 22.54 in the quartic model, indicating a potential onset of overfitting. Since PRESS and predicted R² values are unavailable for the fifth-order model, the risk of overfitting at that level cannot be excluded. Thus, the quartic model likely offers the best trade-off between complexity and performance.

Table 4.

Model fit analysis across various complexity levels for the reduced quartic model.

Source	Std. dev.	R ²	Adjusted R²	Predicted R²	PRESS
Linear	0.5556	0.606	0.6026	0.594	202.86
2FI	0.2439	0.926	0.9235	0.916	41.96
Quadratic	0.172	0.964	0.9619	0.9571	21.43
Cubic	0.0874	0.991	0.9902	0.9832	8.41
Quartic	0.0747	0.995	0.9928	0.9549	22.54
Fifth	0.0761	0.998	0.9925

Table 5.

Evaluation of model fit and predictive accuracy for the reduced quartic model.

Statistic	Value	Statistic	Value
R ²	0.992	Std. dev.	0.086
Adjusted R²	0.991	Mean	1.06
Predicted R²	0.988	C.V. %	8.11
Adeq precision	166.010

Table 5 supports this conclusion by presenting a detailed evaluation of the reduced quartic model’s fit and predictive accuracy. The model exhibits a low standard deviation of 0.086, indicating minimal dispersion of residuals. It achieves an R² of 0.992, explaining 99.2% of the variance in the dataset. The adjusted R² (0.991) and predicted R² (0.988) are closely aligned, with a difference of only 0.003, suggesting excellent predictive performance and minimal overfitting. The model’s coefficient of variation (C.V.) is 8.11%, reflecting strong precision relative to a mean of 1.06. Furthermore, the adequate precision ratio of 166.010 significantly exceeds the recommended threshold of 4, indicating a robust signal-to-noise ratio.

Table 6 and Figure 4 provide additional insights into the characteristics of the reduced quartic model. Table 6 presents a detailed multicollinearity analysis, including estimated coefficients and associated statistical measures. The model comprises 72 terms, ranging from main effects (TCI, TCO, THI, THO, MFC, and MFH) to complex interaction and higher-order terms. A notable feature of the model is the presence of high VIF values for many terms, indicating severe multicollinearity. For example, the interaction term MFC × MFH exhibits the highest VIF of 2633.66, followed closely by MFC² at 2607.51. These extreme values suggest that these predictors are highly correlated with other variables, which may result in unstable and unreliable coefficient estimates. The main effects MFC and MFH also show substantial VIFs of 552.35 and 426.55, respectively, further highlighting the multicollinearity concern.

Table 6.

Estimated coefficients and multicollinearity analysis for the reduced quartic model.

Factor	Coefficient estimate	df	Standard error	95% CI low	95% CI high	VIF
Intercept	0.1239	1	0.0363	0.0526	0.1951
A-TCI	−0.0683	1	0.0866	−0.2385	0.1018	134.70
B-TCO	0.1965	1	0.1096	−0.0189	0.4118	85.88
C-THI	0.2669	1	0.1506	−0.0290	0.5627	89.08
D-THO	0.2640	1	0.0948	0.0778	0.4502	118.90
E-MFC	1.28	1	0.2662	0.7535	1.80	552.35
F-MFH	−1.35	1	0.2422	−1.82	−0.8722	426.55
AB	0.2765	1	0.2009	−0.1181	0.6712	71.18
AC	0.2539	1	0.1671	−0.0742	0.5821	17.89
AD	0.0507	1	0.1598	−0.2632	0.3646	65.88
AE	−1.05	1	0.3542	−1.75	−0.3545	748.34
AF	0.8579	1	0.2926	0.2833	1.43	126.29
BC	−1.24	1	0.3515	−1.93	−0.5524	68.03
BD	−0.6255	1	0.1405	−0.9014	−0.3496	39.24
BE	7.96	1	0.6962	6.59	9.32	1068.53
BF	−7.51	1	0.6733	−8.83	−6.18	452.42
CD	−0.0839	1	0.1426	−0.3640	0.1963	11.88
CE	−5.72	1	0.6357	−6.97	−4.47	629.54
CF	5.63	1	0.6528	4.35	6.92	193.60
DE	−2.59	1	0.3610	−3.30	−1.88	558.23
DF	3.02	1	0.3296	2.37	3.66	253.85
EF	−5.32	1	0.8814	−7.05	−3.59	2633.66
A²	0.0256	1	0.0511	−0.0748	0.1261	21.09
B²	0.5860	1	0.2098	0.1739	0.9981	58.98
C²	−0.0749	1	0.1184	−0.3074	0.1576	11.61
D²	−0.3551	1	0.0868	−0.5256	−0.1846	30.23
E²	4.92	1	0.7736	3.40	6.44	2607.51
F²	0.5525	1	0.1982	0.1633	0.9418	98.54
ABC	0.0156	1	0.4495	−0.8673	0.8986	18.77
ABD	−0.9409	1	0.2221	−1.38	−0.5046	25.25
ABE	0.1665	1	0.4440	−0.7055	1.04	97.56
ABF	0.1249	1	0.4005	−0.6618	0.9115	32.96
ACD	−0.0090	1	0.3312	−0.6596	0.6415	9.45
ACF	−0.6584	1	0.7250	−2.08	0.7655	38.24
ADE	−0.8230	1	0.3083	−1.43	−0.2175	77.27
AEF	2.93	1	0.5970	1.75	4.10	336.18
BCE	−0.2077	1	0.7495	−1.68	1.26	142.20
BCF	1.24	1	0.8812	−0.4899	2.97	69.77
BDF	−0.7311	1	0.3508	−1.42	−0.0421	53.33
BEF	−14.25	1	1.85	−17.89	−10.62	1700.29
CEF	12.62	1	1.46	9.74	15.49	512.79
DEF	3.94	1	0.5733	2.82	5.07	403.31
A²B	0.3025	1	0.1968	−0.0840	0.6890	34.33
A²C	−1.07	1	0.2186	−1.50	−0.6386	18.53
A²D	−0.2653	1	0.1229	−0.5067	−0.0239	26.51
A²F	0.0814	1	0.1871	−0.2860	0.4488	33.11
AE²	−2.35	1	0.4335	−3.20	−1.50	664.07
B²E	−0.8221	1	0.5367	−1.88	0.2321	119.28
B²F	−0.0018	1	0.4422	−0.8703	0.8667	57.83
BC²	0.0646	1	0.2901	−0.5051	0.6344	14.05
BD²	0.8798	1	0.1720	0.5419	1.22	41.80
BE²	4.86	1	1.14	2.62	7.10	1254.03
BF²	8.30	1	0.8672	6.60	10.00	321.15
CE²	−5.56	1	0.8384	−7.20	−3.91	534.69
CF²	−7.18	1	0.7276	−8.61	−5.75	98.95
DE²	−0.7351	1	0.5257	−1.77	0.2974	553.41
DF²	−3.58	1	0.2971	−4.16	−3.00	103.35
E²F	−5.79	1	0.8996	−7.56	−4.02	2269.18
B³	−0.7647	1	0.2662	−1.29	−0.2418	66.94
C³	−0.1401	1	0.1450	−0.4249	0.1446	10.09
E³	0.8256	1	0.6779	−0.5058	2.16	2458.58
F³	5.22	1	0.3288	4.58	5.87	214.35
ABCF	3.01	1	1.31	0.4341	5.58	23.77
ABEF	−1.29	1	0.5799	−2.43	−0.1485	30.29
BCEF	1.13	1	1.33	−1.49	3.75	85.58
A²CD	1.90	1	0.4427	1.03	2.77	7.19
A²CF	−1.55	1	0.8499	−3.22	0.1219	30.87
B³E	−4.56	1	0.9388	−6.41	−2.72	201.70
B³F	6.03	1	0.9125	4.24	7.82	122.60
BC³	0.9133	1	0.4765	−0.0225	1.85	18.76
DE³	−0.9503	1	0.3001	−1.54	−0.3609	114.20
B⁴	0.6286	1	0.2045	0.2269	1.03	23.60

Figure 4.

Statistical analysis of the reduced quartic model’s fit for the crude oil refining process.

The coefficient estimates vary considerably in both magnitude and direction. The largest positive coefficient is 12.62 for the THI × MFC × MFH interaction, while the most negative is −14.25 for TCO × MFC × MFH. Many coefficients—particularly those associated with higher-order interaction terms—exhibit large standard errors relative to their estimated values, resulting in wide 95% confidence intervals that often include zero. This observation suggests that many of these complex terms may not be statistically significant predictors in the model.

Figure 4 presents a statistical evaluation of the reduced quartic model’s fit for the crude oil refining process. The model yielded a notably low PRESS value of 6.19, indicating strong predictive performance. This result suggests that the model’s predictions are likely to be accurate when applied to new, unseen data. Additionally, the −2 log-likelihood value of −1411.37 further supports the model’s excellent fit to the observed data. However, the Bayesian Information Criterion (BIC) of −945.69 and the corrected Akaike Information Criterion (AICc) of −1248.96 are relatively high, which raises concerns. Elevated values for BIC and AICc suggest a poor balance between model fit and complexity, indicating potential overfitting or unnecessary model complexity.

To address the complexity issues identified in the reduced quartic model, Tables 7 and 8 present analyses of variance (ANOVA) for the original and improved models. These analyses aim to resolve the challenges associated with model overcomplexity while maintaining predictive accuracy. Table 7 presents the ANOVA results for the reduced quartic model, highlighting both its strengths and limitations. The model demonstrates strong explanatory power, with an F-value of 947.23 and a p-value of <0.0001, indicating that the model is statistically significant overall. Among the main effects, MFC (E) and MFH (F) exert substantial influence (p < 0.0001), whereas TCI (A) shows minimal impact (p = 0.4305). Several interaction terms—BE, BF, CE, CF, DE, and DF—are also highly significant (p < 0.0001), indicating complex interactions among variables. Higher-order terms such as F³ (p < 0.0001) and B⁴ (p = 0.0022) capture significant nonlinear effects. However, despite its high explanatory power, the model demonstrates a significant lack of fit (F = 12.21, p < 0.0001), raising concerns about its ability to fully capture variability in the data. In addition, elevated BIC and AICc values suggest the presence of unnecessary or redundant terms, contributing to excessive model complexity.

Table 7.

Analysis of variance results for the reduced quartic model.

Source	Sum of squares	df	Mean square	F-value	p-Value
Model	495.39	71	6.98	947.23	<0.0001
A-TCI	0.0046	1	0.0046	0.6223	0.4305
B-TCO	0.0236	1	0.0236	3.21	0.0737
C-THI	0.0231	1	0.0231	3.14	0.0770
D-THO	0.0571	1	0.0571	7.76	0.0055
E-MFC	0.1694	1	0.1694	22.99	<0.0001
F-MFH	0.2282	1	0.2282	30.98	<0.0001
AB	0.0140	1	0.0140	1.89	0.1693
AC	0.0170	1	0.0170	2.31	0.1291
AD	0.0007	1	0.0007	0.1006	0.7512
AE	0.0648	1	0.0648	8.79	0.0032
AF	0.0633	1	0.0633	8.60	0.0035
BC	0.0921	1	0.0921	12.50	0.0004
BD	0.1460	1	0.1460	19.82	<0.0001
BE	0.9619	1	0.9619	130.58	<0.0001
BF	0.9157	1	0.9157	124.31	<0.0001
CD	0.0025	1	0.0025	0.3458	0.5567
CE	0.5967	1	0.5967	81.01	<0.0001
CF	0.5486	1	0.5486	74.48	<0.0001
DE	0.3784	1	0.3784	51.37	<0.0001
DF	0.6166	1	0.6166	83.71	<0.0001
EF	0.2686	1	0.2686	36.47	<0.0001
A²	0.0018	1	0.0018	0.2511	0.6165
B²	0.0575	1	0.0575	7.80	0.0054
C²	0.0029	1	0.0029	0.4000	0.5274
D²	0.1232	1	0.1232	16.73	<0.0001
E²	0.2979	1	0.2979	40.45	<0.0001
F²	0.0573	1	0.0573	7.77	0.0055
ABC	8.887E-06	1	8.887E-06	0.0012	0.9723
ABD	0.1322	1	0.1322	17.94	<0.0001
ABE	0.0010	1	0.0010	0.1406	0.7078
ABF	0.0007	1	0.0007	0.0972	0.7554
ACD	5.487E-06	1	5.487E-06	0.0007	0.9782
ACF	0.0061	1	0.0061	0.8248	0.3642
ADE	0.0525	1	0.0525	7.13	0.0078
AEF	0.1770	1	0.1770	24.03	<0.0001
BCE	0.0006	1	0.0006	0.0768	0.7817
BCF	0.0146	1	0.0146	1.98	0.1596
BDF	0.0320	1	0.0320	4.34	0.0376
BEF	0.4365	1	0.4365	59.25	<0.0001
CEF	0.5470	1	0.5470	74.26	<0.0001
DEF	0.3485	1	0.3485	47.31	<0.0001
A²B	0.0174	1	0.0174	2.36	0.1248
A²C	0.1758	1	0.1758	23.87	<0.0001
A²D	0.0343	1	0.0343	4.66	0.0313
A²F	0.0014	1	0.0014	0.1895	0.6635
AE²	0.2171	1	0.2171	29.47	<0.0001
B²E	0.0173	1	0.0173	2.35	0.1261
B²F	1.277E-07	1	1.277E-07	0.0000	0.9967
BC²	0.0004	1	0.0004	0.0497	0.8237
BD²	0.1927	1	0.1927	26.16	<0.0001
BE²	0.1336	1	0.1336	18.13	<0.0001
BF²	0.6748	1	0.6748	91.62	<0.0001
CE²	0.3236	1	0.3236	43.93	<0.0001
CF²	0.7170	1	0.7170	97.33	<0.0001
DE²	0.0144	1	0.0144	1.96	0.1626
DF²	1.07	1	1.07	145.29	<0.0001
E²F	0.3049	1	0.3049	41.40	<0.0001
B³	0.0608	1	0.0608	8.25	0.0042
C³	0.0069	1	0.0069	0.9343	0.3342
E³	0.0109	1	0.0109	1.48	0.2238
F³	1.86	1	1.86	252.47	<0.0001
ABCF	0.0388	1	0.0388	5.27	0.0221
ABEF	0.0363	1	0.0363	4.93	0.0268
BCEF	0.0053	1	0.0053	0.7143	0.3984
A²CD	0.1357	1	0.1357	18.42	<0.0001
A²CF	0.0244	1	0.0244	3.32	0.0692
B³E	0.1740	1	0.1740	23.62	<0.0001
B³F	0.3218	1	0.3218	43.69	<0.0001
BC³	0.0271	1	0.0271	3.67	0.0558
DE³	0.0739	1	0.0739	10.03	0.0016
B⁴	0.0696	1	0.0696	9.45	0.0022
Residual	4.21	572	0.0074
Lack of fit	4.19	542	0.0077	12.21	<0.0001
Pure error	0.0190	30	0.0006
Cor total	499.60	643

Table 8.

Analysis of variance results for the selected model.

Source	Sum of squares	df	Mean square	F-value	p-Value
Model	480.39	19	25.28	821.09	<0.0001
A-TCI	2.69	1	2.69	87.46	<0.0001
D-THO	0.8463	1	0.8463	27.48	<0.0001
E-MFC	1.16	1	1.16	37.79	<0.0001
F-MFH	3.21	1	3.21	104.20	<0.0001
AB	0.8038	1	0.8038	26.11	<0.0001
AD	1.30	1	1.30	42.09	<0.0001
AF	2.63	1	2.63	85.47	<0.0001
BC	2.52	1	2.52	81.80	<0.0001
BE	53.38	1	53.38	1733.50	<0.0001
BF	11.00	1	11.00	357.27	<0.0001
CE	21.40	1	21.40	694.94	<0.0001
CF	5.79	1	5.79	188.04	<0.0001
DE	10.61	1	10.61	344.68	<0.0001
DF	12.04	1	12.04	391.15	<0.0001
EF	2.44	1	2.44	79.26	<0.0001
B²	2.82	1	2.82	91.57	<0.0001
D²	1.90	1	1.90	61.62	<0.0001
E²	18.21	1	18.21	591.36	<0.0001
F²	4.76	1	4.76	154.50	<0.0001
Residual	19.21	624	0.0308
Lack of fit	19.20	594	0.0323	50.99	<0.0001
Pure error	0.0190	30	0.0006
Cor total	499.60	643

Table 8 presents the ANOVA for a revised model designed to resolve these overcomplexity issues. This improved model includes only 19 terms, significantly reducing model complexity while preserving strong explanatory performance (F = 821.09, p < 0.0001). All retained terms are statistically significant (p < 0.0001). The primary effects A–TCI, D–THO, E–MFC, and F–MFH are especially influential, with F–MFH exhibiting the highest F-value (104.20). Among the interaction terms, BE (1733.50), CE (694.94), and DF (391.15) show exceptionally high F-values, highlighting strong synergistic effects. Quadratic terms (B², D², E², and F²) are also present, with E² being the most significant (F-value = 591.36), suggesting pronounced nonlinear effects. The simplified model reduces the number of terms from 71 in the reduced quartic model to 19, improving interpretability, reducing overfitting risk, and enhancing generalizability. Despite these improvements, the simplified model still exhibits a significant lack of fit (F = 50.99, p < 0.0001), implying that some systematic variation in the data remains unaccounted for.

The selected model is presented in Tables 9 and 10 and compared visually in Figure 5. Table 9 outlines the estimated coefficients for the selected model, highlighting key interactions in the crude oil refining process. This model includes several prominent interaction terms, notably BE (MFC × THI), which has a coefficient of 4.71—indicating a strong positive synergistic effect between these variables. This interaction also features a narrow confidence interval (4.49–4.94) and a low VIF of 6.76, underscoring its reliability and statistical importance. In contrast, the CE (MFC × TCO) interaction demonstrates a significant negative effect, with a coefficient of −3.47 and a confidence interval ranging from −3.73 to −3.21.

Table 9.

Coefficients for the selected model.

Factor	Coefficient estimate	df	Standard error	95% CI low	95% CI high	VIF
Intercept	−0.0020	1	0.0369	−0.0744	0.0704
A-TCI	−0.2602	1	0.0278	−0.3148	−0.2056	3.32
D-THO	0.3529	1	0.0673	0.2207	0.4851	14.35
E-MFC	0.9394	1	0.1528	0.6393	1.24	43.55
F-MFH	−1.33	1	0.1306	−1.59	−1.08	29.68
AB	−0.5018	1	0.0982	−0.6946	−0.3089	4.07
AD	0.5217	1	0.0804	0.3638	0.6796	3.99
AF	−0.9987	1	0.1080	−1.21	−0.7865	4.12
BC	−1.24	1	0.1368	−1.51	−0.9689	2.47
BE	4.71	1	0.1132	4.49	4.94	6.76
BF	−3.29	1	0.1743	−3.64	−2.95	7.25
CE	−3.47	1	0.1316	−3.73	−3.21	6.45
CF	2.86	1	0.2088	2.45	3.27	4.74
DE	−2.66	1	0.1435	−2.95	−2.38	21.09
DF	3.03	1	0.1531	2.73	3.33	13.10
EF	−2.41	1	0.2711	−2.95	−1.88	59.59
B²	0.8329	1	0.0870	0.6619	1.00	2.43
D²	−0.5669	1	0.0722	−0.7087	−0.4251	5.00
E²	4.63	1	0.1903	4.25	5.00	37.75
F²	−1.77	1	0.1424	−2.05	−1.49	12.17

Table 10.

Evaluation of the fit and predictive accuracy of the selected model.

Statistic	Value	Statistic	Value	Statistic	Value
Std. dev.	0.1755	R ²	0.9615	PRESS	21.91
Mean	1.06	Predicted R²	0.9561	BIC	−304.80
C.V. %	16.58	Adeq precision	141.8546	AICc	−392.81

Figure 5.

Comparison of actual and predicted values.

The primary effects differ significantly, with MFH (F) having the most detrimental impact (−1.33) and MFC (E) exhibiting the most favorable impact (0.9394). However, the MFC term shows a high VIF of 43.55, suggesting the presence of multicollinearity. The quadratic terms indicate notable nonlinear effects, particularly for E² (MFC²), which has a substantial positive coefficient of 4.63, suggesting significant curvature in the response surface. However, E² also has a high VIF (37.75), implying possible correlation with other predictors. Both positive and negative quadratic effects are observed: B² (THI²) and E² (MFC²) exhibit positive effects, while D² (THO²) and F² (MFH²) exhibit negative effects. Among these, E² appears to be the most influential quadratic term in the refining process, followed by F² (−1.77). Despite the elevated VIF values, the confidence intervals for these terms are relatively narrow, indicating high precision in the coefficient estimates.

The fit and predictive accuracy of the selected crude oil refining model are evaluated in Table 10. The model demonstrates a high coefficient of determination (R²) of 0.9615, accounting for 96.15% of the variance in the dataset. Its strong predictive capability—without evidence of overfitting—is confirmed by a predicted R² of 0.9561. The model’s standard deviation of 0.1755 indicates low residual dispersion, while a coefficient of variation (C.V.) of 16.58% reflects moderate relative variability. Furthermore, the model’s ability to navigate the design space effectively is supported by a high signal-to-noise ratio of 141.8546, substantially exceeding the recommended threshold of 4.

The model’s predictive strength is further reinforced by a low PRESS value of 21.91, which supports the high predicted R² and confirms the model’s reliability despite its reduced complexity. The negative values of the Bayesian Information Criterion (BIC = −304.80) and the corrected Akaike Information Criterion (AICc = −392.81) indicate a favorable balance between model fit and complexity, with lower values generally preferred. These improvements demonstrate that the simplified model performs better than the more complex reduced quartic model in terms of both parsimony and accuracy. Figure 5 illustrates that the selected model aligns more closely with the actual data points compared to the previous model, reinforcing its utility despite its more streamlined structure.

Figure 6 illustrates the relationships between key parameters and the coded response factor (RFC), complementing the statistical analysis presented in the preceding tables. Several variables are shown to influence heat exchanger fouling resistance and efficiency. In Figure 6(a), the RFC decreases from 0.2 to −0.2 as the crude oil inlet temperature (TCI) increases from 13°C to 31°C, indicating a modest negative linear trend. This observation aligns with the positive correlation between TCI and RFC shown in Figure 3, suggesting that higher TCI levels may contribute to reduced fouling resistance. As depicted in Figure 6(b), the RFC peaks between 62°C and 71.5°C, exhibiting a nonlinear, quadratic relationship with the hot fluid outlet temperature (THO). This behavior supports the presence of significant quadratic effects, as reflected in Table 9, further validating the inclusion of quadratic terms in the selected model.

Figure 6.

Effects of key operating parameters on the fouling resistance (RFC) of a heat exchanger system.

Figure 6(c) and (d) provide insights into heat exchanger system optimization, highlighting complex nonlinear interactions. Figure 6(c) illustrates a parabolic relationship between the cold fluid mass flow rate (MFC) and RFC, with a minimum fouling resistance of approximately 0 m²°C/W observed at MFC values between 58 and 63 kg/s. This suggests that a moderate increase in MFC can reduce fouling resistance, whereas excessively high flow rates may lead to increased fouling due to enhanced turbulence. Therefore, careful calibration of MFC is required to achieve optimal performance. The inclusion of both linear (MFC) and quadratic (MFC²) terms in the selected model (see Table 9) captures this complex relationship. Similarly, Figure 6(d) reveals a concave relationship between the hot fluid mass flow rate (MFH) and RFC. While moderate MFH values appear to slightly reduce fouling resistance, significantly high flow rates result in a substantial decline in RFC. However, optimizing for high MFH values must consider potential trade-offs, such as increased energy consumption and efficiency losses. As shown in Table 9, the selected model includes both negative linear (MFH) and quadratic (MFH²) terms, accurately reflecting these observed correlations.

Table 11 presents the optimization outcomes, synthesizing the insights from the preceding analysis. The results indicate that the application of the desirability function technique effectively reduced the heat exchanger fouling resistance (RFC). This outcome is attributed to the predictive capacity of the selected model and the parameter relationships illustrated in Figure 6. Out of 100 evaluated scenarios, only 37 achieved the target RFC of 0.001 m²°C/W with an ideal desirability score of 1.000. These optimal solutions exhibit consistent patterns in several parameters, most notably a hot fluid outlet temperature (THO) of 43.00°C and a hot fluid mass flow rate (MFH) of 93.00 kg/s. In contrast, other parameters—such as the cold fluid inlet temperature (TCI) and its mass flow rate (MFC)—exhibited greater variability. As shown in Figure 6(c) and (d), the interactions among MFC, MFH, and RFC are complex and nonlinear, featuring both parabolic and inverse parabolic trends. These findings underscore the importance of carefully balancing thermal and flow variables to optimize fouling resistance.

Table 11.

Optimal parameter configurations for minimizing fouling in heat exchangers.

Number	TCI	TCO	THI	THO	MFC	MFH	RFC	Desirability
1	25.74	92.19	130.88	43.00	53.10	93.00	0.001	1.000
2	29.40	98.01	123.30	43.00	25.33	93.00	0.001	1.000
3	17.13	94.17	122.98	43.00	41.18	93.00	0.001	1.000
4	17.44	106.93	133.21	43.00	28.37	93.00	0.001	1.000
5	18.51	94.56	125.42	43.00	43.55	93.00	0.001	1.000
6	27.78	95.90	120.49	43.00	26.55	93.00	0.001	1.000
7	27.11	94.28	125.70	43.00	39.00	93.00	0.001	1.000
8	27.91	96.31	133.77	43.00	46.04	93.00	0.001	1.000
9	16.30	117.24	124.11	43.00	88.46	93.00	0.001	1.000
10	23.05	104.55	131.94	43.00	27.67	93.00	0.001	1.000
11	27.46	103.65	135.91	43.00	32.38	93.00	0.001	1.000
12	17.03	110.65	135.30	43.00	23.42	93.00	0.001	1.000
13	27.76	92.30	123.74	43.00	40.11	93.00	0.001	1.000
14	30.85	92.91	133.83	43.00	52.27	93.00	0.001	1.000
15	26.06	99.47	130.68	43.00	35.33	93.00	0.001	1.000
16	21.89	101.61	132.96	43.00	36.76	93.00	0.001	1.000
17	27.69	94.53	122.02	43.00	32.18	93.00	0.001	1.000
18	29.81	92.65	127.03	43.00	43.33	93.00	0.001	1.000
19	24.99	102.60	135.58	43.00	36.03	93.00	0.001	1.000
20	30.79	92.48	114.83	43.00	23.49	93.00	0.001	1.000
21	29.43	96.11	134.12	43.00	45.99	93.00	0.001	1.000
22	30.79	103.58	135.46	43.00	29.68	93.00	0.001	1.000
23	24.55	95.69	118.37	43.00	25.29	93.00	0.001	1.000
24	20.63	92.95	129.62	43.00	53.03	93.00	0.001	1.000
25	26.75	94.19	125.02	43.00	38.37	93.00	0.001	1.000
26	16.20	92.17	130.29	43.00	60.09	93.00	0.001	1.000
27	17.00	101.44	134.77	43.00	43.90	93.00	0.001	1.000
28	23.02	92.57	125.66	43.00	45.55	93.00	0.001	1.000
29	21.42	98.13	124.19	43.00	31.23	93.00	0.001	1.000
30	22.26	102.29	125.87	43.00	23.99	93.00	0.001	1.000
31	13.75	92.77	132.46	43.00	65.44	93.00	0.001	1.000
32	14.53	98.93	135.23	43.00	53.18	92.99	0.001	1.000
33	29.69	95.15	133.30	43.02	46.93	93.00	0.001	1.000
34	14.50	116.39	135.23	43.00	90.49	92.97	0.001	1.000
35	30.54	97.90	134.37	43.00	41.39	93.00	0.003	1.000
36	13.60	101.11	114.15	43.00	92.99	93.00	0.003	1.000
37	21.32	92.23	112.44	43.00	24.52	92.96	0.002	1.000
38	22.26	92.01	114.45	43.00	28.11	92.86	0.001	0.999
39	29.03	98.77	122.72	43.00	23.00	92.78	0.001	0.999
40	27.71	97.60	126.54	43.00	32.16	93.00	0.033	0.998
41	16.29	99.82	127.59	43.00	36.87	92.36	0.001	0.997
42	18.17	92.00	136.00	43.14	69.17	92.62	0.006	0.997
43	13.00	97.29	118.88	43.02	29.28	93.00	0.068	0.995
44	22.86	118.64	130.17	43.00	93.00	91.92	0.001	0.994
45	25.90	92.00	112.23	43.02	23.00	91.73	0.001	0.993
46	13.00	104.15	108.10	43.04	89.05	91.08	0.001	0.990
47	13.00	101.34	134.01	43.00	45.91	93.00	0.144	0.989
48	25.21	92.00	127.16	44.99	50.60	93.00	0.001	0.985
49	31.00	116.91	108.64	45.08	92.85	93.00	0.001	0.984
50	14.64	97.38	110.18	43.58	93.00	90.64	0.001	0.983
51	26.39	119.00	130.40	43.00	93.00	89.47	0.001	0.982
52	21.73	106.93	122.32	43.00	93.00	89.39	0.001	0.981
53	20.73	112.44	123.54	45.51	92.85	93.00	0.001	0.981
54	24.30	104.38	107.55	43.62	93.00	89.59	0.001	0.978
55	30.74	118.92	126.41	43.00	92.99	88.17	0.001	0.975
56	31.00	96.76	133.57	46.37	47.00	92.92	0.001	0.973
57	27.00	103.92	107.00	43.00	92.94	87.90	0.001	0.973
58	16.61	92.00	136.00	46.45	76.82	92.96	0.007	0.972
59	20.22	99.45	119.50	43.00	93.00	87.45	0.001	0.971
60	13.00	109.44	136.00	45.93	34.33	93.00	0.095	0.970
61	27.03	114.47	128.96	43.00	93.00	87.30	0.001	0.970
62	25.12	101.11	108.47	43.00	92.99	87.13	0.001	0.969
63	18.46	112.83	135.63	47.03	23.21	93.00	0.001	0.969
64	13.00	97.53	120.91	43.00	35.23	86.66	0.001	0.967
65	13.00	92.00	129.73	46.36	69.50	91.54	0.001	0.967
66	20.46	101.87	114.90	43.06	93.00	88.65	0.139	0.966
67	19.49	105.68	136.00	43.07	93.00	86.20	0.002	0.964
68	23.29	96.02	115.36	43.00	23.01	86.00	0.001	0.963
69	25.23	106.31	107.02	43.00	92.95	89.40	0.257	0.962
70	28.43	92.00	109.51	48.01	23.01	93.00	0.003	0.960
71	22.81	92.00	134.14	48.05	68.08	93.00	0.001	0.960
72	14.98	105.43	136.00	48.00	93.00	92.31	0.001	0.957
73	22.27	96.74	120.12	43.00	93.00	84.90	0.001	0.957
74	21.85	95.45	125.79	43.10	92.82	83.67	0.001	0.949
75	20.83	102.99	135.92	49.35	47.05	92.96	0.023	0.947
76	23.98	97.19	125.97	43.00	93.00	82.82	0.001	0.945
77	31.00	92.00	119.93	43.00	37.42	82.42	0.001	0.943
78	24.50	110.39	136.00	43.00	23.23	81.86	0.001	0.940
79	13.00	117.47	108.72	43.00	72.80	81.33	0.001	0.937
80	26.50	119.00	136.00	43.00	87.96	80.80	0.029	0.932
81	22.24	113.71	136.00	49.14	23.00	89.07	0.001	0.931
82	13.84	104.88	130.05	43.72	36.81	80.89	0.001	0.929
83	28.34	119.00	107.08	43.00	77.93	79.59	0.001	0.927
84	31.00	119.00	134.56	43.05	89.55	80.10	0.042	0.927
85	13.00	92.61	129.70	50.21	77.71	89.86	0.001	0.926
86	13.00	100.25	113.49	49.22	23.10	88.34	0.014	0.926
87	30.99	92.00	107.97	49.81	23.52	88.32	0.031	0.920
88	13.00	111.39	107.14	43.01	73.11	78.08	0.001	0.918
89	27.49	99.90	134.92	43.00	93.00	77.42	0.001	0.914
90	21.22	116.00	136.00	53.60	23.07	93.00	0.001	0.911
91	17.77	105.53	107.00	43.00	76.55	75.64	0.001	0.903
92	31.00	98.15	135.68	54.61	56.95	93.00	0.001	0.901
93	13.05	111.93	124.98	54.64	23.00	93.00	0.015	0.900
94	13.00	92.00	121.03	43.01	82.41	76.91	0.194	0.897
95	31.00	92.00	110.83	52.34	32.23	87.15	0.001	0.894
96	18.20	94.21	135.61	55.90	93.00	93.00	0.001	0.888
97	18.76	112.67	127.57	56.05	23.05	93.00	0.001	0.887
98	31.00	117.63	107.00	43.29	74.45	73.85	0.057	0.886
99	20.87	117.74	136.00	56.29	23.00	93.00	0.001	0.884
100	30.51	92.00	131.30	46.78	93.00	76.56	0.001	0.882

The results validate the accuracy of the model presented in Table 9 for forecasting system behavior and offer valuable insights into the operational flexibility of heat exchangers. The model’s ability to capture complex system dynamics enables it to achieve optimal RFC values across a wide range of input conditions, indicating the existence of multiple efficient operating strategies. The consistent presence of high hot fluid mass flow rates (MFH) among the top-performing solutions reinforces the finding that elevated MFH values significantly reduce fouling resistance. This is supported by the negative linear and quadratic coefficients for MFH in the selected model (Table 9). These optimization results further confirm the intricate interactions among system variables, as indicated by the significant interaction terms in the prior analysis. Moreover, they provide practical guidance for operating heat exchanger systems to minimize fouling resistance while maintaining thermal efficiency under diverse operating conditions. Although the simplified model is less complex than the reduced quartic model, the optimization outcomes demonstrate its robustness and practical applicability in supporting operational decision-making (Table 9).

Discussion

The development of the quartic regression model for predicting fouling resistance in crude oil refining heat exchangers represents a significant advancement in balancing statistical rigor with operational utility. The observed progression in model performance—from linear to higher-order polynomials—demonstrates the classic trade-off between model fit and overfitting, a key consideration in thermal systems modeling. This trend aligns with the findings of Jafari et al.,³⁷ who reported optimal predictive performance using response surface methodology (RSM) for spiral heat exchangers, while also cautioning against the multicollinearity risks introduced by interaction and quadratic terms. The present study addresses this issue directly through variance inflation factor (VIF) analysis, systematically removing terms with excessive collinearity while preserving essential interactions.

The statistical validation metrics achieved by the reduced quartic model reflect significant improvements over previous approaches. The substantial reduction in predictive residual error sum of squares (PRESS = 21.91) and the high signal-to-noise ratio (141.85) demonstrate superior predictive accuracy and model robustness. These improvements surpass those reported in comparable studies, including the central composite design (CCD) model proposed by Jradi et al.²¹ (R² = 0.988) and the VIF-agnostic quadratic model of Wang et al.¹⁵ Notably, the present model retains only the most impactful nonlinear terms (e.g. MFC² and MFH²), consistent with nanofluid behavior described by Maddah et al.,¹⁶ while introducing a replicable term-elimination framework—an element often overlooked in the design of experiments literature.^19,20

From a practical standpoint, the optimization scenarios summarized in Table 11 identified operational settings—specifically, THO = 43°C and MFH = 93 kg/s—that achieved minimal fouling resistance (RFC = 0.001 m²°C/W), aligning with the 10.4% efficiency gain trends reported by Guo et al.³⁸ However, the nuanced balance between high MFH rates and turbulence-induced fouling, as illustrated in Figure 6(d), underscores the importance of strategic parameter tuning. This observation mirrors the operational complexity highlighted by Al Kumait et al.¹⁷ in nanofluid-enhanced systems and reinforces the need for a multifactorial optimization framework.

A key innovation of this study lies in the quantitative control of multicollinearity through VIF thresholds. By systematically excluding terms with VIF values exceeding 2000—such as MFC × MFH—the model avoids the coefficient instability reported in previous studies. This approach is consistent with best practices described by Jaberi and Ghassemi,³⁹ who addressed similar concerns in factorial design modeling for desalination systems. The final model, reduced by 73% (from 72 to 19 terms), maintained high predictive performance (R² = 0.9615; predicted R² = 0.9561), exceeding the R² range (0.943–0.999) typically reported for RSM-based models.^15,21

Furthermore, the use of the desirability function within Design-Expert software highlighted the role of nonlinear curvature effects, particularly in the relationship between MFC and RFC. This finding resonates with the work of Li et al.,⁴⁰ who showed that neglecting curvature in finned-tube heat exchanger optimization led to suboptimal outcomes. Likewise, the identification of an optimal THO near 71°C for fouling minimization supports the conclusions drawn by Jiang et al.⁴¹ regarding crystallization fouling thresholds in non-metallic exchangers.

The model’s interaction terms also provide deeper insights into process behavior. For instance, the parabolic response of RFC to MFC (see Figure 6(c)) aligns with the nonlinear turbulence–resistance relationships reported by Ibrahim et al.⁴² The identification of specific thresholds—such as the RFC minimum in the 58–63 kg/s MFC range—offers practical guidance not typically available from CFD-based studies,¹² reaffirming the value of empirical modeling for real-time decision-making.

Finally, the integration of diagnostic tools such as PRESS, R² disparity, and lack-of-fit analysis facilitated effective overfitting control. This approach is consistent with findings by Song et al.,⁴³ who emphasized the adverse effects of multicollinearity on model reliability in plate heat exchangers. Despite modest lack-of-fit indicators, the performance of the simplified model suggests that further improvements may be achieved by hybridizing with machine learning or physics-informed frameworks—an approach already explored by Alqahtani et al.⁴⁴ and van Veen⁴⁵ in nanofluid and steam injection modeling, respectively.

Conclusions

This study was conducted to address the critical need for effectively modeling and optimizing the performance of heat exchanger systems in crude oil refining operations, with a particular focus on minimizing fouling resistance. The objective was to develop a predictive model that balances statistical adequacy and interpretability while reliably estimating system behavior under varying operating conditions. Given the complexity and nonlinear interactions among thermal and flow parameters, the aim was to identify a model that not only provides a high-quality fit to empirical data but also avoids the typical pitfalls of overfitting and multicollinearity associated with high-complexity regression models.

To achieve this, a series of polynomial regression models—ranging from linear to fifth-order—were developed and evaluated using Design-Expert software. Performance was assessed using statistical indicators such as R², adjusted R², predicted R², PRESS, and VIF, alongside ANOVA diagnostics and multicollinearity analyses. The initial results showed that increased model complexity improved variance explanation but at the cost of overfitting, as indicated by the lack of projected R² and PRESS values for the highest-order models. The initial reduced quartic model demonstrated near-optimal fit (adjusted R² = 0.9928; predicted R² = 0.9549), but suffered from severe multicollinearity, with VIF values exceeding 2600 for some predictors—limiting its interpretability and generalizability.

A simplified model was then developed, reducing the number of terms from 71 to 19. This model retained high predictive power (R² = 0.9615; predicted R² = 0.9561), with a lower residual standard deviation (0.1755) and a significantly reduced PRESS (21.91), indicating robust predictive reliability. It also achieved improved information criteria values (BIC = −304.80; AICc = −392.81), reflecting a more balanced trade-off between model fit and simplicity. While multicollinearity was not entirely eliminated—some terms, such as MFC, still exhibited elevated VIFs—its most extreme manifestations were effectively mitigated, contributing to greater model stability.

The final model yielded meaningful insights into system behavior. It identified MFH (hot fluid mass flow rate) and MFC (cold fluid mass flow rate) as the most influential variables, both in their linear effects and through their interactions and quadratic terms. Notably, the MFC × THI and MFC × TCO interaction terms exhibited strong synergistic and antagonistic effects, respectively, as evidenced by their tight confidence intervals and high F-values. The quadratic effects of MFC² and MFH² revealed significant nonlinear trends relevant to system optimization. These findings were visually supported by response surface plots and interaction diagrams, which showed parabolic and concave patterns in RFC behavior. For instance, RFC was minimized within the 58–63 kg/s range of MFC and maximized at intermediate THO values, highlighting the complexity of the response surface and the necessity of multidimensional optimization.

In summary, this study demonstrates that simplifying a high-complexity regression model to retain only statistically significant and operationally meaningful terms can yield a robust, interpretable, and highly predictive framework for optimizing heat exchanger performance in crude oil refining. The final model preserves key linear, interaction, and nonlinear effects while mitigating overfitting and multicollinearity. Through rigorous statistical evaluation and optimization analysis, it offers actionable insights for engineers and plant operators, supporting efforts to reduce fouling resistance, enhance energy efficiency, and improve overall process reliability. Future research should explore the real-time implementation of such models, the incorporation of dynamic process variables, and the integration of machine learning techniques to increase adaptability and scalability across diverse industrial contexts.

Footnotes

Handling Editor: Sharmili Pandian

ORCID iDs

Nadjem Bailek

Jihad A. Younis

Raheem Al-Sabur

Author contributions

Conceptualization: K.I., K.D., N.B.; methodology and software validation: K.I., K.D., N.B., B.Z.; formal analysis and writing—original draft: K.I., K.D., N.B., R.A., D.A., J.A.Y., B.Z., L.M., C.A.G.S.; writing—review and editing: N.B., C.A.G.S., L.M., I.C.; visualization: K.I., N.B., I.C. All authors have read and agreed to the published version of the manuscript.

Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Data availability statement

The authors confirm that the data supporting the findings of this study are available within the article.

References

Olympios

McTigue

Farres-Antunez

, et al. Progress and prospects of thermo-mechanical energy storage—a critical review. Prog Energy 2021; 3: 22001.

Tavousi

Perera

Flynn

, et al. Heat transfer and fluid flow characteristics of the passive method in double tube heat exchangers: a critical review. Int J Thermofluids 2023; 17: 100282.

Roy

Majumder

Economic optimization and energy analysis in shell and tube heat exchanger by meta-heuristic approach. Vacuum 2019; 166: 413–418.

Hameed

Mohammed

Abdullah

MR.

Improving the thermal performance of a heat exchanger using a new passive technology. Tikrit J Eng Sci 2023; 30: 66–71.

Delrot

Guerra

Dambrine

, et al. Fouling detection in a heat exchanger by observer of Takagi–Sugeno type for systems with unknown polynomial inputs. Eng Appl Artif Intell 2012; 25: 1558–1566.

Musawel

RK.

Simulation and fouling study of propane heat exchangers. MSc Thesis, University of Basrah, Iraq, 2002.

Coletti

Hewitt

Crude oil fouling: deposit characterization, measurements, and modeling. Gulf Professional Publishing, 2014.

Diaz-Bejarano

Behranvand

Coletti

, et al. Organic and inorganic fouling in heat exchangers: industrial case study analysis of fouling rate. Ind Eng Chem Res 2018; 58: 228–246.

Al Hadad

Schick

Maillet

. Fouling detection in a shell and tube heat exchanger using variation of its thermal impulse responses: methodological approach and numerical verification. Appl Therm Eng 2019; 155: 612–619.

10.

Das

O’Connell

, et al. Assessing advances in anti-fouling membranes to improve process economics and sustainability of water treatment. ACS EST Eng 2022; 2: 2159–2173.

11.

Ikram

Djilali

Abdennasser

, et al. Comparative analysis of fouling resistance prediction in shell and tube heat exchangers using advanced machine learning techniques. Res Eng Struct Mater 2024; 10: 253–270.

12.

Berce

Zupančič

Može

, et al. Infrared thermography observations of crystallization fouling in a plate heat exchanger. Appl Therm Eng 2023; 224: 120116.

13.

Elghool

Basrawi

Ibrahim

, et al. Enhancing the performance of a thermo-electric generator through multi-objective optimisation of heat pipes-heat sink under natural convection. Energy Convers Manag 2020; 209: 112626.

14.

Sahin

A Taguchi approach for determination of optimum design parameters for a heat exchanger having circular-cross sectional pin fins. Heat Mass Transf 2007; 43: 493–502.

15.

Wang

Xiao

Wang

, et al. Application of response surface method and multi-objective genetic algorithm to configuration optimization of Shell-and-tube heat exchanger with fold helical baffles. Appl Therm Eng 2018; 129: 512–520.

16.

Maddah

Aghayari

Mirzaee

, et al. Factorial experimental design for the thermal performance of a double pipe heat exchanger using Al₂O₃-TiO₂ hybrid nanofluid. Int Commun Heat Mass Transf 2018; 97: 92–102.

17.

Al Kumait

AAR

Ibrahim

Abdullah

MA.

Experimental and numerical study of forced convection heat transfer in different internally ribbed tubes configuration using TiO₂ nanofluid. Heat Transf Res 2019; 48: 1778–1804.

18.

Wang

Shi

Fang

, et al. Heat transfer enhancement for drag-reducing surfactant fluid using photo-rheological counterion. Exp Heat Transf 2012; 25: 139–150.

19.

Bisognin

Bastos

JCSC

Meier

, et al. Influence of different parameters on the tube-to-bed heat transfer coefficient in a gas-solid fluidized bed heat exchanger. Chem Eng Process Intensif 2020; 147: 107693.

20.

Chowdhury

Borah

. Modeling and optimization of the performance parameters of a single pass shell and tube heat exchanger—an approach using response surface methodology. In: Manik

Kalia

Verma

, et al. (eds) Recent advances in mechanical engineering. Springer, 2022, pp.999–1016.

21.

Jradi

Marvillet

Jeday

MR.

Parametric study of calcium sulfate crystallization fouling in cross-flow heat exchanger using response surface methodology. Heat Mass Transf 2023; 59: 1971–1985.

22.

Hallaji

Peyghambarzadeh

Bohloul

, et al. The optimum conditions for calcium sulfate fouling rate under subcooled flow boiling using Taguchi statistical method. Int J Heat Mass Transf 2023; 204: 123859.

23.

Ferreira

SLC

Bruns

Ferreira

, et al. Box-Behnken design: an alternative for the optimization of analytical methods. Anal Chim Acta 2007; 597: 179–186.

24.

Mohapatra

Sahoo

Padhi

BN.

Analysis, prediction and multi-response optimization of heat transfer characteristics of a three fluid heat exchanger using response surface methodology and desirability function approach. Appl Therm Eng 2019; 151: 536–555.

25.

Derringer

Suich

Simultaneous optimization of several response variables. J Qual Technol 1980; 12: 214–219.

26.

Dahmani

Kouidri

Bailek

, et al. Enhanced hourly temperature prediction using advanced ensemble neural networks for energy system efficiency optimization in hyper-arid regions. AIP Adv 2025; 15: 045019.

27.

Ferkous

Guermoui

Bellaour

, et al. Enhancing photovoltaic energy forecasting: a progressive approach using wavelet packet decomposition. Clean Energy 2024; 8: 95–108.

28.

Dahmani

Ammi

Bailek

, et al. Assessing the efficacy of improved learning in hourly global irradiance prediction. Comput Mater Contin 2023; 77(2): 2579–2594.

29.

Hussein

Zerouali

Bailek

, et al. Harnessing explainable AI for sustainable agriculture: SHAP-based feature selection in multi-model evaluation of irrigation water quality indices. Water 2024; 17: 59.

30.

Difi

Heddam

Zerouali

, et al. Improved daily streamflow forecasting for semi-arid environments using hybrid machine learning and multi-scale analysis techniques. J Hydroinformatics 2024; 26: 3266–3286.

31.

Masrur Ahmed

Bailek

Abualigah

, et al. Global control of electrical supply: a variational mode decomposition-aided deep learning model for energy consumption prediction. Energy Rep 2023; 10: 2152–2165.

32.

Almorox

Arnaldo

Bailek

, et al. Adjustment of the Angstrom-Prescott equation from Campbell-Stokes and Kipp-Zonen sunshine measures at different timescales in Spain. Renew Energy 2020; 154: 337–350.

33.

Raudenbush

Bryk

AS.

Hierarchical linear models: applications and data analysis methods. Vol. 1. Sage, 2002.

34.

Allen

DM.

The relationship between variable selection and data agumentation and a method for prediction. Technometrics 1974; 16: 125–127.

35.

Bailek

Bouchouicha

El-Shimy

, et al. Improved mathematical modeling of the hourly solar diffuse fraction (HSDF)-Adrar, Algeria case study. Sol Energy 2017; 4: 8–12.

36.

Harche

Mouheb

Absi

The fouling in the tubular heat exchanger of Algiers refinery. Heat Mass Transf 2016; 52: 947–956.

37.

Jafari

Khorshidi

Sheibani

, et al. Performance investigation and optimization of spiral heat exchangers with different fluids by response surface methodology. Iran J Chem Chem Eng 2025; 44: 1222–1238.

38.

Guo

Wang

Liang

, et al. Temperature zone diagram method for designing the total site exchanger network. Therm Sci Eng Prog 2023; 43: 101979.

39.

Jaberi

Ghassemi

Optimization of electrode design for electrodialysis reversal. Report, Desalination & Water Purification Research and Development Program, USA, 2015.

40.

Wang

Chen

, et al. Fabrication and characterization of superhydrophobic Al-based surface used for finned-tube heat exchangers. Materials (Basel) 2022; 15: 3060.

41.

Jiang

Sun

, et al. Study on the effect of surface properties of non-metallic materials on the growth mechanism of crystallization fouling. Processes 2023; 11: 2232.

42.

Ibrahim

Al-Sammarraie

Al-Taha

, et al. Experimental and numerical investigation of heat transfer augmentation in heat sinks using perforation technique. Appl Therm Eng 2019; 160: 113974.

43.

Song

Lim

Yun

, et al. Composite fouling characteristics of CaCO₃ and CaSO₄ in plate heat exchangers at various operating and geometric conditions. Int J Heat Mass Transf 2019; 136: 555–562.

44.

Alqahtani

Alenazy

Pandiaraj

Optimization on heat-transfer coefficient of plate heat exchanger using MWCNT-TiO₂ nanofluid by response surface methodology. Mater Res Express 2020; 7: 84002.

45.

van Veen

THP

. Pilot plant scale experiments for direct steam injection. Master Thesis, Eindhoven University of Technology, Netherlands, 2023.