Optimal Cleaning Cycle Scheduling under Uncertain Conditions: A Flexibility Analysis on Heat Exchanger Fouling

Abstract

Fouling is a substantial economic, energy, and safety issue for all the process industry applications, heat transfer units in particular. Although this phenomenon can be mitigated, it cannot be avoided and proper cleaning cycle scheduling is the best way to deal with it. After thorough literature research about the most reliable fouling model description, cleaning procedures have been optimized by minimizing the Time Average Losses (TAL) under nominal operating conditions according to the well-established procedure. For this purpose, different cleaning actions, namely chemical and mechanical, have been accounted for. However, this procedure is strictly related to nominal operating conditions therefore perturbations, when present, could considerably compromise the process profitability due to unexpected shutdown or extraordinary maintenance operations. After a preliminary sensitivity analysis, the uncertain variables and the corresponding disturbance likelihood were estimated. Hence, cleaning cycles were rescheduled on the basis of a stochastic flexibility index for different probability distributions to show how the uncertainty characterization affects the optimal time and economic losses. A decisional algorithm was finally conceived in order to assess the best number of chemical cleaning cycles included in a cleaning supercycle. In conclusion, this study highlights how optimal scheduling is affected by external perturbations and provides an important tool to the decision-maker in order to make a more conscious design choice based on a robust multi-criteria optimization.

1. Introduction

Due to the higher CAPital EXpenses (CAPEX) related to equipment oversizing and OPerating EXpenses (OPEX) related to energy and production losses as well as frequent maintenance, fouling still represents, nowadays, a relevant issue for the process industry.

The total heat exchanger fouling costs for highly industrialized countries were about 0.25% of the countries’ Gross National Product (GNP) in 1992 [1,2].

Table 1. Fouling vs. Gross National Product (GNP) (2018).

Country	Costs (M$/a)	GNP 2018 (M$/a)	Costs/GNP (%)
US	14,175	20,891,000	0.13%
Germany	4875	4,356,353	0.21%
France	2400	2,962,799	0.15%
Japan	10,000	5,594,452	0.33%
Australia	463	1,318,153	0.06%
New Zealand	64.5	197,827	0.06%

shows the annual costs of fouling in some different countries based on the actualization of the money value of the 1992 estimation, considering 2018 GNP.

As it can be noticed, additional energy duties required to compensate for an ineffective heat transfer and frequent unit cleaning and maintenance result in highly relevant expenses. A considerable part of this work is then focused on the OPEX optimization for already designed heat exchanger systems.

The overall fouling process is the net result of two simultaneous sub-processes, namely a deposition and a removal process [3,4]. The combination of these basic phenomena affects the growth of the deposit on the surface, mathematically defined as the rate of deposit growth (fouling resistance or fouling factor, $R_{\mathrm{f}}$).

The two best ways to describe fouling in a model suitable way have been found in literature as the “Two-layer model” and “Distributed model”. The latter still does not result as reliable as the first one for optimization purposes, due to its computational intensiveness [5].

In order to keep the energetic and economic efficiencies of the whole heat transfer process acceptable, heat exchangers cleaning should be often performed. Cleaning methodologies can be classified into two main groups, namely chemical cleaning and mechanical cleaning.

Mechanical cleaning methods completely restore the heat transfer surface from fouling but may damage the equipment. On the contrary, chemical cleaning techniques, which are not actually able to completely remove the fouling layers, do not cause stress or damage to the heat exchanger internals. As will be later discussed, in the developed model mechanical cleaning removes both the layers, while chemical cleaning only removes the one exposed to the fluid.

This physical description of the fouling process is useful to assess the system heat transfer performances and thus to find the optimal operation scheduling. Batch process scheduling is indeed a major research topic in process system engineering and its optimization algorithms were widely studied during the last decades [6,7,8].

With the advent of major concerns related to the sustainability issue, the environmental considerations concerning energy and waste analysis of batch processes have become part of the scheduling optimization domain thanks to several studies carried out by the most influential exponents of the topic [9,10].

Until the end of the 20th century, the vast majority of the studies concerning batch processes design were mainly referred to as nominal operating conditions, i.e., no uncertainty was accounted for. Moreover, the scheduling optimization and the control design were performed in two different steps of the process design procedure even though they substantially affect each other.

With the advances in design under uncertainty methodologies and the spread of its applications to thermodynamics [11], unit operations [12,13,14,15,16], reacting systems [17], and other fields of process engineering, the way process scheduling was conceived has started changing even if—differently from process units—a considerably smaller amount of publications about this topic is available in the literature. A pioneering work under this perspective was performed by Balasubramanian and Grossmann that analyzed the scheduling optimization under uncertain processing times with a branch and bound [18] and with a fuzzy programming [19] approach. Bonfill et al. (2005) [20] later discussed the scheduling optimization accounting for variable market demand with a two-stage stochastic optimization strategy.

Two papers, with similar research topics to that presented in this publication, were proposed by Smaïli et al. (1999) [21] and Diaby et al. (2016) [22]. The first one deals with the optimal scheduling of cleaning cycles in heat exchanger networks subject to fouling and the uncertainty analyzed refers to the measurements used to derive the fouling model. On the other hand, Diaby et al. propose the use of the genetic algorithm in order to optimize heat exchanger networks as well.

Uncertainty in batch process scheduling may refer as well to the control system performances to optimize the operation times. Further details concerning the current state of the art and the new challenges of simultaneous scheduling design and control under uncertainties can be found in the interesting literature review recently carried out by Dias and Ierapetritou (2016) [23].

The main innovative aspect of this article with respect to the literature presented here above is that, given a heat exchanger unit undergoing fouling, the optimal cleaning scheduling is thoroughly analyzed and optimized accounting for uncertain operating conditions, i.e., the so-called aleatory uncertainty related to the process input parameters is discussed under a flexibility point of view. The energy analysis of the system is coupled indeed with the economic one in order to provide a reliable estimation of each operating cycle time as well as of the associated cleaning costs.

In order to have a clear insight into this research study, the case study followed by further details and hypotheses concerning the fouling process, the flexibility assessment, and the optimization algorithm is presented in the following sections.

2. The Heat Exchanger Case Study

This research work aims at the application of the newly proposed procedure starting from the most elementary unit in order to set the basis for the possible system scale up to the several possible arrangements of multiple heat exchanger to form a more complex network. Therefore, the selected case study is a simple heat exchanger unit undergoing fouling as shown in Figure 1. The purpose of this unit is the effective heat recovery between hot and cold hydrocarbon process streams. The former passes on the tube sides while the latter on the shell side of the unit.

Two additional external duty exchangers are present to make the streams achieve the desired temperature specifications in case the heat recovery was not sufficient. The value of these heat duties expressed in Watts are defined as Q_cold and Q_hot respectively.

The process parameters for this case study were accurately selected in order to be compliant with the industrial practice. Based on the same case study proposed by Ishiyama et al. 2011 [24] and later used by Pogiatzis et al. (2012) [5], some process parameters, in particular specific heat capacities and heat transfer coefficient, were adjusted to more reliable values for the given mixtures according to the most reputable literature in the process engineering domain [25,26,27,28].

The equipment sizing was performed as well, in the case of maximum heat recovery for a fixed minimum temperature approach under nominal operating conditions and a heat transfer surface area equal to 140 m² was found.

The complete list of the obtained process parameters used for this research work is reported in

Table 2. Process parameters.

	Symbols	Quantities	Value	Unit
Heat exchanger	$U$	Overall heat transfer coefficient	350	W/(m2 K)
	$A$	Surface area	140	m2
Cold stream	$q_{\mathrm{m},\mathrm{c}}$	Mass flow rate	135	kg/s
	$T_{\mathrm{c}}^{\mathrm{in}}$	Inlet temperature	525	K
	$c_{\mathrm{P}}^{\mathrm{c}}$	Specific heat capacity	3125	J/(kg K)
Hot stream	$q_{\mathrm{m},\mathrm{h}}$	Mass flow rate	68	kg/s
	$T_{\mathrm{h}}^{\mathrm{in}}$	Inlet temperature	575	K
	$c_{\mathrm{P}}^{\mathrm{h}}$	Specific heat capacity	2200	J/(kg K)
Fouling	$K_{\mathrm{g}}$	Deposition rate	5 × 10−6	(m2 K)/(W d)
	$K_{\mathrm{c}}$	Ageing rate	2.5 × 10−7	(m2 K)/(W d)
	${\lambda }_{\mathrm{g}}$	Gel thermal conductivity	0.1	W/(m K)
	${\lambda }_{\mathrm{c}}$	Coke thermal conductivity	0.8	W/(m K)

. The details concerning the cleaning cycle scheduling and optimization procedure as well as the flexibility analysis applied to the case study here above are presented and discussed in detail in the following section.

3. Methodology

The methodology section deals with the three main aspects related to cleaning cycle scheduling under uncertain operating conditions.

The first part to be defined concerns the fouling kinetic model used to describe the fouling physical phenomenon that will affect heat transfer effectiveness dynamics and will be used to assess the optimal cycle duration. Moreover, the two different cleaning techniques are presented in order to highlight the effect of the cleaning process on the fouling parameters.

After that, the definition of uncertain operating conditions via the flexibility index is carried out. The analysis is focused both on deterministic and stochastic indexes in order to provide a complete overview of flexibility and the way it can be quantified and associated with process variables.

Finally, the scheduling optimization and the corresponding decisional algorithm are discussed and compared to those already employed in the available literature. The definition of the algorithm is indeed required to outline a thorough and general procedure to be used no matter the case study under analysis.

3.1. Fouling Kinetics and Cleaning Techniques

3.1.1. Fouling Kinetics

The fouling model employed in this work is the so-called “two-layer model” proposed in similar publications by Ishiyama et al. (2011) [24] and Pogiatzis et al. (2012) [5].

As stated by its name, this model assumes a fouling deposition distributed over two layers characterized by different physical properties and defined as:

• Gel layer, located at the interface between the hard solid deposition and the process fluid;
• Coke layer, formed between the Gel layer and the exchanger heat transfer surface.

Please note that “coke” as defined here above is used only for fouling formed starting from an organic compound, otherwise, it should be referred to as “crust”.

The zeroth-order kinetic model assumed to describe the behavior of the growing layers are:

• Coke layer:

$$\{\begin{array}{c}\frac{\mathrm{d}{\delta }_{\mathrm{c}}}{\mathrm{d}t}={\lambda }_{\mathrm{c}}\cdot K_{\mathrm{c}}, {\delta }_{\mathrm{g}}>0\\ \frac{\mathrm{d}{\delta }_{\mathrm{c}}}{\mathrm{d}t}=0, {\delta }_{\mathrm{g}}=0\end{array}$$

(1)

• Gel layer:

  $$\frac{\mathrm{d}{\delta }_{\mathrm{g}}}{\mathrm{d}t}={\lambda }_{\mathrm{g}}\cdot K_{\mathrm{g}}-\frac{\mathrm{d}{\delta }_{\mathrm{c}}}{\mathrm{d}t}$$

  (2)

  where δ_g and δ_c are the thicknesses and where λ_g and λ_c are the thermal conductivities of the gel and coke layer, respectively.

The sum of these basic components represents the deposit growth on the heat exchanger surface. The thermal fouling resistance (or fouling factor) of the fouling layer R_f, also defined as the rate of deposit growth, is evaluated by treating the layers as a pair of thin slabs of insulating material and mathematically described as:

$$R_{\mathrm{f}}={\mathsf{\Phi }}_{\mathrm{d}}-{\mathsf{\Phi }}_{\mathrm{r}}=\frac{{\delta }_{\mathrm{g}}}{{\lambda }_{\mathrm{g}}}+\frac{{\delta }_{\mathrm{c}}}{{\lambda }_{\mathrm{c}}}$$

(3)

where Φ_d and Φ_r are the rates of deposition and removal respectively. R_f, as well as the deposition and the removal rate, can be expressed in the units of thermal resistance as [(m² K)/W] or in the units of the rate of thickness change as [m/s] or units of mass change as [kg/(m² s)] [3].

Finally, in the fouling model adopted, the overall heat transfer coefficient U is given by:

$$\frac{1}{U}=\frac{1}{U_{\mathrm{clean}}}+R_{\mathrm{f}}\iff U=\frac{U_{\mathrm{clean}}}{1+U_{\mathrm{clean}}\cdot R_{\mathrm{f}}}$$

(4)

where U_clean indicates its value at t = 0.

3.1.2. Cleaning Techniques

Both the two different cleaning techniques presented by Ishiyama et al. (2011) [24] have been studied in this research work. They can be classified and qualitatively described as follows:

• Chemical Cleaning: the heat exchanger is shut down in order to perform the cleaning. Only the gel layer is removed with this technique by means of a proper solvent. This cleaning procedure is generally shorter and has lower fixed costs with respect to mechanical cleaning;
• Mechanical Cleaning: the heat exchanger is shut down in order to perform the cleaning procedure. Deposits are completely removed by means of mechanical strength and the heat exchanger recovers its original heat transfer capacity. This cleaning procedure is generally longer and has higher fixed costs with respect to chemical cleaning.

The cost related to the two techniques are those suggested in the referenced publications. Although the results are affected by those values, the procedure and the related code are of general validity whatever the cost function is.

The effect of these cleaning techniques on the model equations concerns the terms δ_g and δ_c. In particular, the chemical cleaning process resets the gel layer δ_g to zero, restoring the heat transfer coefficient defined by the Equation (4) to the value:

$$\frac{1}{U}=\frac{1}{U_{\mathrm{clean}}}+\frac{{\delta }_{\mathrm{c}}}{{\lambda }_{\mathrm{c}}}\iff U=\frac{U_{\mathrm{clean}}\cdot {\lambda }_{\mathrm{c}}}{{\lambda }_{\mathrm{c}}+U_{\mathrm{clean}}\cdot {\delta }_{\mathrm{c}}}$$

(5)

It is worth remarking that, once the chemical cleaning is performed, the coke layer keeps on increasing starting from the thickness achieved before the unit shut down.

On the other hand, the mechanical cleaning process is able to remove the entire fouling layer so that the heat transfer coefficient is equal to U_clean once the unit is restarted.

The analysis of more cleaning cycle techniques with respect to the usual scheduling based on the Mechanical Cleaning only implies the possibility of multiple cleaning cycles configurations obtained by combining them. The implications related to what will be called “supercycle” and the way we should treat it during the design phase are better discussed in Section 3.3.

3.2. Flexibility Indices

Flexibility, defined as the ability of a process to accommodate a set of uncertain parameters [29], can be quantified by the mean of several indices proposed in the literature. As already observed by Di Pretoro et al. (2019) [13], flexibility indices can be classified into deterministic and stochastic according to the way uncertainty is treated. In this research work, both a deterministic index and a stochastic one have been employed in order to compare their behavior on the applied case study.

The main definitions and mathematical formulations are discussed here below in order to provide a complete overview before their application to the cleaning cycle scheduling problem.

The first to be defined and the most widespread flexibility index is the Swaney and Grossmann (1985) [30] one.

Let us define the general hyperrectangle:

$$T\left(\delta \right)=\{\theta :{\theta }^N-\delta \cdot \mathsf{\Delta }{\theta }^{-}<\theta <{\theta }^N+\delta \cdot \mathsf{\Delta }{\theta }^{+}\}$$

(6)

With δ non-negative scalar variable for which

$$\{\begin{array}{c}\delta <1,\hspace{1em}\hspace{1em}T\left(\delta \right)\in T\\ \delta >1,\hspace{1em}\hspace{1em}T\left(\delta \right)\ni T\end{array}$$

(7)

The flexibility index F_SG, proposed by Swaney and Grossmann, is then the solution to the constrained optimization flexibility problem:

$$F_{\mathrm{SG}}=\underset{}{\max }\delta$$

(8)

$$\mathrm{s}.\mathrm{t}.\underset{\theta \in T\left(\delta \right)}{\max }\underset{z}{\min }\underset{j\in J}{\max }{f}_j\left(d,z,\theta \right)\le 0$$

(9)

where d and z represent the design and the manipulated variables respectively. By referring to Figure 2, F_SG graphically represents the maximum scale factor so that the corresponding hyperrectangle is bounded by the feasible zone. Moreover, for constraints jointly quasi-convex in z and 1D quasi-convex in θ, the solution lies at a vertex of the hyperrectangle allowing to solve the optimization problem by evaluating the feasibility of the design at each vertex. On the contrary, certain types of non-convex domains may lead to nonvertex solutions.

The index presented here above can be referred to as “deterministic” since it is based on the assessment of a “perturbation magnitude” in the uncertain domain without taking into account the perturbation likelihood. Thus, this index definition results in a rather conservative estimation of the system flexibility.

For these reasons, a few years later, Pistikopoulos and Mazzuchi (1990) [31] introduced a stochastic flexibility index based on the perturbation likelihood. Given the uncertain parameters Probability Distribution Function (hereafter PDF) P(θ) and the feasible region

$$\mathsf{\Psi }\left(d,z,\theta \right)\le 0$$

(10)

the stochastic flexibility index SF can be defined as:

$$SF={{\displaystyle \int }}_{\mathsf{\Psi }}^{}P\left(\theta \right)\cdot d\theta$$

(11)

Unlike the deterministic indices, the stochastic flexibility one does not quantify a perturbation magnitude but assesses the percentage of the possible operating conditions for which the system keeps being feasible. However, the higher cost to pay for such an accurate result does not lie on the higher need for data only. As shown in Figure 3 indeed, more than one subregion in the uncertain domain corresponding to the same value of integral (11) could exist.

Therefore, the design choice and the related Time Average Losses (TAL, hereafter denoted by Θ) associated with each value of SF is the result of an economic optimization problem stated as:

$$\Theta \left(SF\right)=\underset{D}{\max }\Theta$$

(12)

where

$$D=D^k|\underset{D^k}{\max }\Theta \le \underset{D^j}{\max }\Theta \forall j :{{\displaystyle \int }}_{D^k}^{}P\left(\theta \right)\cdot d\theta ={{\displaystyle \int }}_{D^j}^{}P\left(\theta \right)\cdot d\theta =SF$$

(13)

This formulation states that, for each value SF of the stochastic flexibility index, the associated subregion D, represented by the maximum value of Θ in it, is the cheapest one among all the subregions satisfying the Equation (11). Thus, for each step of the stochastic flexibility analysis an additional optimization problem should be solved, considerably increasing the required computational effort.

The particular definition of this index makes it different from the deterministic ones by two main properties as pointed out by Di Pretoro et al. (2019) [13]:

• The stochastic flexibility index SF has always a value between 0 and 1;
• Under nominal operating conditions, SF can be higher than 0.

The importance of including flexibility when assessing operating costs has already been proved by Di Pretoro et al. (2020) [32]. However, in order to do that the decisional algorithm for the operation scheduling to which the flexibility indices will be applied should be thoroughly defined. For this reason, a detailed explanation will be discussed in the following section.

3.3. Decisional Algorithm for Operation Scheduling

Given the heat exchanger case study under analysis with the related cleaning techniques and the associated fouling kinetic model described in Section 2 and Section 3.1 respectively, the decisional algorithm for the cleaning operation schedule should give the answers to two main questions:

• WHEN the cleaning operation should be performed?
• WHICH cleaning operation should be performed?

In order to do that, some cost functions are required. The fact that the cycle is optimized with respect to an economic criterion implies that, once questions 1 and 2 are answered, a third question “what are the expected operating costs?” is answered as well.

The first cost function to use is called the “Energy Loss” function and it is defined as:

$$EL=C_{\mathrm{E}}\cdot {{\displaystyle \int }}_0^t\left(Q_{\mathrm{cl}}-Q\right)\cdot dt$$

(14)

where C_E is the cost per energy [€/J], t is the operating time, Q_cl is the heat exchanged under clean conditions, i.e., at t = 0, and Q is the actual heat recovery at t = t_op. In particular, it can be observed that the integrated function represents the external heat duty demand, i.e., the sum between Q_cold and Q_hot. This cost function describes operating costs to be afforded when the heat recovery is still ongoing.

Two other cost functions should be defined in order to quantify the expenses incurring during the unit shutdown. The first one refers to the cost of the energy to be provided by the external duty when no heat recovery is performed. This value is equal to:

$$E{L}_{\mathrm{C}}=C_{\mathrm{E}}\cdot Q_{\mathrm{cl}}\cdot {\tau }_{\mathrm{C}}$$

(15)

for solvent cleaning or

$$E{L}_{\mathrm{M}}=C_{\mathrm{E}}\cdot Q_{\mathrm{cl}}\cdot {\tau }_{\mathrm{M}}$$

(16)

for mechanical cleaning where ${\tau }_i$ is the time required for the corresponding cleaning operation.

Finally, the last cost item is represented by the cost of the cleaning operation itself, abbreviated with $C_{\mathrm{C}}^{\mathrm{Cl}}$ and $C_{\mathrm{M}}^{\mathrm{Cl}}$ for the chemical and mechanical cleaning techniques respectively. The costs of cleaning have been assumed as constant in time. Although the code does not account for inflation and market trends’ influence, it is still possible to implement within it these aspects for real applications. However, this is not the purpose of this work.

Once all these cost items are calculated, the optimization problem concerning the optimal cycle duration needs to be solved. Its general formulation is:

$$\underset{t, {\tau }_{\mathrm{M}},{\tau }_{\mathrm{C}}}{\min }\Theta \left(t,{\tau }_{\mathrm{M}},{\tau }_{\mathrm{C}}\right)$$

(17)

$$\Theta \left(t,{\tau }_{\mathrm{M}},{\tau }_{\mathrm{C}}\right)={\Theta }_{\mathrm{M}}\left(t,{\tau }_{\mathrm{M}}\right)+{\Theta }_{\mathrm{C}}\left(t,{\tau }_{\mathrm{C}}\right)$$

(18)

This formulation is used to solve the so-called supercycle, i.e., the sequence of multiple chemical cleaning operations (N_C) required before performing the mechanical one. In particular, whether a cleaning operation is required, i.e., when the energy losses are higher than the cost of cleaning, the minimum value between Θ_M and Θ_C determines which kind of cleaning technique should be employed. In fact, each time the solvent cleaning is preferred, the coke layer keeps increasing until the chemical operation becomes more expensive than the mechanical one and the so-called supercycle is closed.

In the case where only mechanical cleaning is considered, the optimization problem is simplified by removing all the chemical cleaning related terms.

The decisional algorithm described here above has been graphically summarized in the flowchart reported in Figure 4.

The following section presents the obtained results both for the conventional procedure and for the stochastic flexibility based one. In the latter case, the economic optimization carried out to select the optimal scheduling has been coupled with the optimization required by the stochastic flexibility assessment. This innovative procedure allows outlining the trend of the cycle duration and the associated costs obtained from the first optimization as a function of the flexibility index resulting from the second one.

The heat exchanger and fouling models, the flexibility analysis, and the decisional algorithm for scheduling optimization calculations were all performed by means of MatLab^® codes.

4. Results

The following sections present in detail the results obtained. In particular, they have been classified into three main parts. The first one shows the results of the scheduling algorithm applied to our case study under nominal operating conditions both in case of a single cleaning process and of the combination of the two.

A sensitivity analysis with respect to two possible uncertain variables will follow in order to detect the most critical one and use it to perform the flexibility assessment.

Finally, the third part presents and discusses in detail the outcome of the stochastic flexibility assessment for the single heat exchanger case study by coupling economic and operational aspects.

4.1. Nominal Operating Conditions

This section presents the results obtained with the process parameters discussed in Section 2 accounting both for the mechanical cleaning methodology only and for both chemical and mechanical ones.

4.1.1. Mechanical Cleaning Only

The scheduling optimization has been performed first according to the conventional methodology accounting for mechanical cleaning only. The time average losses have been calculated for a range of cycle time and the t_op corresponding to its minimum value has been detected.

The overall Energy Losses and the TAL trend as a function of the operation time are shown in Figure 5a,b respectively. As it can be noticed, for low t_op values, the TAL are very high since the cleaning costs are much higher than the cost of the external duty required to compensate for the lower performance of the heat exchanger due to the fouling phenomenon. This value decreases until a minimum of 306.16 $/d in correspondence of 195 days and starts increasing again for the longer t_op because of the higher impact of the EL with respect to the cleaning expenses.

4.1.2. Supercycle

The decisional algorithm for scheduling optimization has been then applied accounting for both the gel and the coke layers according to the methodology presented in Section 3.3. The TAL cost function was evaluated for different numbers of chemical cleaning cycles N_C before performing the mechanical operation.

Figure 6a shows the optimal TAL value for N_C ranging from 0 to 5. The function in 0 represents the result obtained when no chemical cleaning is performed before the mechanical one. This value is exactly the one obtained for a simple cycle scheduling in the previous section, i.e., 306.16 $/d. The TAL trend shows a minimum of 234 $/d in correspondence of N_C equal to 2. This means that the optimal supercycle consists of two chemical cleanings and a mechanical one. For lower values, i.e., 1 and 0, the energy losses due to the lower efficiency of the heat recovery do not justify the cost of a mechanical operation. On the other hand, for higher values, the impact of the increasing coke layer on the heat transfer implies substantial energy losses that would compromise the profitability of the operation more than the cost of complete cleaning.

In order to ease the understanding of this phenomenon, the gel and coke layer thicknesses as a function of time have been represented in Figure 6b for the optimal operation scheduling.

The time interval when the layer thickness has a stationary trend refers to cleaning operation time, i.e., the unit maintenance shutdown time. In particular, it can be noticed that the last shutdown required for the mechanical cleaning is longer than the previous ones (5 days vs. 1 day) due to the higher complexity related to the coke layer removal.

The most relevant outcome of this section is that the combination of two different cleaning strategies allows to save up to 24% of the time average losses thanks to the lower complexity and duration of the solvent cleaning subcycle.

4.2. Sensitivity Analysis

The sensitivity analysis was performed on the hot side inlet temperature $T_{\mathrm{h}}^{\mathrm{in}}$ and the heat transfer surface area A. The choice of perturbing both an input parameter and a design variable has the purpose to verify which one of them has a higher impact on the operation scheduling. The most critical one indeed will be selected as an uncertain parameter to perform the stochastic flexibility assessment.

Figure 7 shows the results of the sensitivity analysis for wide ranges of the perturbed variables. In particular, Figure 7a refers to the optimal cycle TAL while Figure 7b shows the corresponding optimal cycle time.

On the one hand, in the former Figure, an almost linear trend of the Time Average Losses as a function of the two perturbed variables can be noticed. More precisely, for the same deviation in terms of percentage, the TAL shows a slightly higher sensitivity with respect to $T_{\mathrm{h}}^{\mathrm{in}}$ than to A. In particular, the sensitivity of the cost function with respect to one of the two variables becomes more relevant as the value of the other one gets higher.

On the other hand, in the latter Figure, the higher sensitivity of the optimal cycle time with respect to the temperature results to be much more relevant. The perturbation of t_op with respect to the heat transfer surface area is almost negligible for low $T_{\mathrm{h}}^{\mathrm{in}}$ values and shows a moderate response for very high-temperature perturbations. On the contrary, the optimal cycle time exhibits an almost exponential growth with respect to the inlet temperature for any value of the heat transfer surface area. This behavior can be reconducted to the fact that this parameter directly acts on the external duty demand required to achieve the hot side specification while the heat exchanger was already properly sized and affects the optimal scheduling solution only in case of very ineffective heat transfer.

Thus, in the light of these results, the intermediate value A = 140 m² was finally kept as the nominal one while the hot side inlet temperature has been selected as uncertain parameters for the flexibility analysis presented in detail in the next section.

4.3. Flexibility Assessment

As already discussed in detail by previous studies dealing with the aleatory uncertainty and the stochastic flexibility index [13,31,32], the uncertainty characterization, i.e., the probability distribution function used to describe the uncertain parameter likelihood, plays a main role in the final result. For this reason, in order to have a more complete overview, two PDFs, namely the Gaussian and the Beta distributions, are employed in this study. The former describes an uncertain variable perturbation with a symmetric behavior with respect to the central point (i.e., its mode) usually corresponding to the nominal operating conditions. On the other hand, the Beta distribution defines an aleatory parameter for which either the positive or the negative deviation is more likely to occur. Further details about their mathematical formulation are provided in the Appendix A.

Both the Gaussian (or Normal) and the Beta distributions belong to the class of two-parameter PDFs, therefore two conditions should be imposed in order to uniquely determine their trends. In the case of real applications, the distribution of the uncertain parameter can be outlined according to the available frequency data or accounting for its expected behavior. For this paper, since the study is not related to any existing unit, the maximum likelihood of the PDF was set equal to the nominal operating conditions point and the variance was selected so that the entire uncertain space is covered with a residual probability lower than the 0.01%.

The trends resulting from these hypotheses have been outlined in Figure 8. As it can be noticed, the Gaussian distribution results particularly narrow with respect to the entire temperature range since there is a minimum inlet temperature below which the heat transfer is not feasible due to a minimum temperature approach lower than 10 °C. This value was then set as the lower boundary of the uncertain domain.

Once the uncertain parameter PDF has been properly defined, the stochastic flexibility assessment can be finally carried out according to the methodology presented in Section 3.2. Each probability distribution function was then integrated over the uncertain domain and the maximum TAL and cycle time have been assessed at each value of the SF index. This procedure allowed to estimate the expected costs and operation time to be respectively afforded and waited for each simple cycle.

The cost trends obtained this way are shown in Figure 9. As it can be noticed, the expected TAL grows fast for very low flexibility index values until a lower slope increased is attained over a wide range of the uncertain parameter. As the residual probability gets lower, much higher losses should be expected in order to almost entirely cover the uncertain domain.

The interval width reflects the PDF shape as expected and has already been discussed in the previous sections. In particular, the “steady” region is much wider in the case of the Gaussian one. The TAL for lower values resulted higher for the Beta PDF while the cost increase in the residual part of the stochastic flexibility index is more relevant when using the Gaussian one.

In regards to the optimal cycle time, the obtained curves show an opposite trend with respect to the corresponding TAL functions. This behavior is due to the fact that higher inlet temperature requires cleaning cycles to be more frequently performed as already pointed out by the sensitivity analysis. Analogously to the economic remarks discussed in Section 4.1.2 about the costs, the use of a combined chemical and mechanical cleaning methodology allows extending the operating time before a complete shutdown is required.

In any case, the two trends are qualitatively similar. It is worth remarking that the presence of some edges in these trends is due to the discretization of the uncertain domain that was not kept too dense in order to avoid excessive computational time.

The same procedure was then carried out for the supercycle accounting for both chemical and mechanical cleaning operations.

The optimal number of solvent cleaning cycles over the uncertain domain always resulted as 2. However, the corresponding costs and optimal cycle time considerably vary over the inlet temperature interval as shown by the plots in Figure 10a,b for the two probability distributions. Given that the optimal supercycle configuration does nott change, relatively smooth trends can be observed for both the output variables and both the PDFs without any particular discontinuous point. As for the simple cycle, a fast increase in the cost function can be observed for very low and very high SF values. However, for the supercycle, the Beta distribution exhibits a higher slope in the intermediate range with respect to the Gaussian one.

Even in this case, an opposite trend can be detected for the optimal cycle time as for the previous analysis. In general, higher costs, and thus lower optimal cycle time, can be pointed out for the Beta distribution due to its higher variance.

Therefore, in the light of these results, it can be also concluded that, in the case of a narrow probability distribution, the scheduling parameters obtained accounting for nominal operating conditions only are much closer to the real estimation. This is due to the fact that the trends TAL and t_op vs. SF exhibit a lower slope over a wider flexibility range and just the very last part of residual probability significantly affects their values.

5. Conclusions

The purpose of this research work was to include the uncertainty of the operating condition in the heat cleaning cycles scheduling optimization algorithm. The outcome of the study was successful and allowed for the definition of a thorough procedure to be applied in the case of heat exchanger cleaning cycle optimal scheduling.

The general validity of the innovative algorithm with respect to the uncertainty characterization was proved by testing both a symmetric (i.e., Gaussian) and a skewed (i.e., Beta) probability distribution function. Moreover, the procedure does not depend on the fouling kinetic model that is used to describe the physical phenomena related to fouling.

Furthermore, the employment of a preliminary sensitivity analysis is proposed in this study with the purpose of restraining the number of uncertain parameters to be accounted for by identifying the most critical ones from a flexibility perspective and thus decreasing the dimension of the uncertain domain.

A further result obtained during this study was the algorithm aimed at the optimally combined solvent and mechanical cleaning cycle under uncertainty. As for nominal operating conditions, the possibility to exploit different cleaning operations allows reducing the costs related to fouling removal. However, the optimization of the so-called “supercycle” implies the introduction of an additional decisional variable, i.e., chemical vs. mechanical, requiring a non-negligible computational effort already in the case of a single unit.

In conclusion, a general procedure for operation scheduling under uncertainty was outlined by means of a stochastic flexibility indicator and, in particular, its application to a single heat exchanger cleaning cycle was analyzed in deep. The aleatory uncertainty characterization by means of a probability distribution function combined with the proposed algorithm provides a tool to the decision-maker to have a more reliable cost estimation and to define the best scheduling solution according to the specific operating conditions.

This work sets the basis for future applications on more complex heat exchanger networks as well as additional perspectives on other kinds of operations different than unit maintenance. Moreover, besides the use of the methodology proposed in this research study to deal with an undesired disturbance, i.e., fouling, the scheduling algorithm could be also actively exploited to take advantage of certain operational parameters in order to reduce operating costs and enhance the process performances.