A Modeling Framework to Investigate the Inﬂuence of Fouling on the Dynamic Characteristics of PID-Controlled Heat Exchangers and Their Networks

: The present work is an extension of the authors’ previous research, where changes in the dynamic behavior of heat exchangers induced by fouling build-up were studied. In the present work, the authors used the previously elaborated and validated mathematical model of transient heat exchange with the inﬂuence of thermal resistance of fouling taken into account. The behavior of speciﬁc Heat Exchanger Networks (HENs) coupled with a Crude Distillation Unit together with their control loops is simulated using Simulink/MATLAB and the inﬂuence of fouling build-up on speciﬁc indices of quality of operation is investigated. According to the presented results, the higher the number of heat exchangers in the PID control loop and the greater the number of heat exchangers interacting in the network, the smaller the inﬂuence of fouling on the control quality indices, and in the extreme case, this inﬂuence may be negligible. This might be caused by the compensation of the negative effects of fouling build-up when the heat exchangers are interacting in the HEN. Nevertheless, potential adverse effects of fouling on HEN operation can be prevented by periodic adjustments of the optimal values of PID gains.


Introduction
Heat Exchanger Network (HEN) fouling is a chronical problem encountered in many process industries. The operation of a Heat Exchanger (HE) may be affected by fouling which builds up on the heat transfer surface. For example, fouling of HENs in the oil refining industry results in increased energy consumption (burning extra fuel to compensate for reduced heat recovery), reduced plant throughput when the exchangers are cleaned, and induced costs of cleaning interventions [1,2]. In recent years, various approaches to the mitigation of fouling effects in industrial HEs and HENs have been reported in the literature [3,4]. Fouling leads to the reduction of steady-state heat recovery [5], but also to changes in the transient states of HEs [6] and inefficient control of HEs [7,8] that may have an adverse effect on the overall performance of the HEN [9,10]. In the literature, publications devoted to the effect of fouling on the dynamic behavior of HEs and the role of fouling in control issues are rare and limited in scope [11,12]. For a more regular introduction please refer to the authors' previous work [6], where the relevant research field has been reviewed and key publications have been cited. Figure 1 shows an example of dynamic characteristics of a HE (output signal y as a function of time t, in response to a step change in an input signal that occurred at t = 0) and indicates possible Figure 1 shows an example of dynamic characteristics of a HE (output signal y as a function of time t, in response to a step change in an input signal that occurred at t=0) and indicates possible deformation of the characteristics induced by the deterioration of control quality. The dynamic characteristics are described by the following parameters: gain K0, delay time td and time constant t1. The gain is a coefficient determining the change of the output signal with respect to the change in the input control signal; the larger the gain value, the stronger output signal's response to the input signal. The delay time defines the waiting period between the step change in the input signal (in Figure 1, t=0) and the change of the output signal (that is, delay time describes the speed of response on given control). The time constant is a measure of the capacity of the process and determines how long after the cessation of the delay does the output signal reach 63.2% of its final value. When fouling builds up on the heat transfer surface, the dynamic characteristics of the HEs operated in the network may be changed. In a previous work [6], the present authors studied the influence of fouling deposition in the individual HE units on their dynamic characteristics and on the quality of their control, with the aim of describing fouling effects quantitatively. Moreover, four examples of control loops with single HE and without significant interactions from the other HEs were investigated. In appendix A, as a complement to the mentioned closed-loop characteristics, the corresponding values of control-quality indices are shown in Table A1. In the present paper, using the validated dynamic model, the behavior of specific crude preheat trains (branches of HENs interacting with each other -Figures 2 and 3) together with their PIDcontrol loops is simulated and the influence of fouling build-up on the specific indices of the quality of operation is investigated. Figures 2 and 3 show two considered examples of real-life PID-controlled HENs coupled with a Crude Distillation Unit (CDU). The graphical form of Figures 2 and 3 was developed by the present authors on the basis of schemes and process data made available by the owner of two different CDUs that operated in a Polish oil refinery.
The conclusions of the previous publication [6] seemed to suggest that fouling build-up in heat exchangers usually leads to the significant deterioration of control-quality indices of PID-based control loops. However, from the results of the present research in which more complex PIDcontrolled loops in large HENs have been investigated, different conclusions can be drawn. It was found that the larger the number of heat exchangers in PID-control loops and the larger the number of interacting heat exchangers in the HEN, the less pronounced is the influence of fouling on the indices of control quality. This observation, being new in the pertinent literature, is supported by the presented simulation results and the discussion of a case study. It may be ascribed to the compensation of adverse effects of fouling build-up in HENs characterized by strong interactions between the heat exchangers. In the present paper, using the validated dynamic model, the behavior of specific crude preheat trains (branches of HENs interacting with each other-Figures 2 and 3) together with their PID-control loops is simulated and the influence of fouling build-up on the specific indices of the quality of operation is investigated. Figures 2 and 3 show two considered examples of real-life PID-controlled HENs coupled with a Crude Distillation Unit (CDU). The graphical form of Figures 2 and 3 was developed by the present authors on the basis of schemes and process data made available by the owner of two different CDUs that operated in a Polish oil refinery.
The conclusions of the previous publication [6] seemed to suggest that fouling build-up in heat exchangers usually leads to the significant deterioration of control-quality indices of PID-based control loops. However, from the results of the present research in which more complex PID-controlled loops in large HENs have been investigated, different conclusions can be drawn. It was found that the larger the number of heat exchangers in PID-control loops and the larger the number of interacting heat exchangers in the HEN, the less pronounced is the influence of fouling on the indices of control quality. This observation, being new in the pertinent literature, is supported by the presented simulation results and the discussion of a case study. It may be ascribed to the compensation of adverse effects of fouling build-up in HENs characterized by strong interactions between the heat exchangers.
The generated knowledge could lead to the development of methods and techniques to prevent heat-recovery reduction that may occur when HEN control is affected by fouling build-up in the exchangers.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 3 of 24 The generated knowledge could lead to the development of methods and techniques to prevent heat-recovery reduction that may occur when HEN control is affected by fouling build-up in the exchangers.

Materials and Methods
Using multi-cell steady-state model of a HE, a control-theory based approach was proposed for the identification and evaluation of the influence of fouling on the dynamic behavior of the heat exchangers and on the quality of their control. A prerequisite for reliable monitoring of the quality of HEN control under fouling conditions is that acquisition and processing of operation data are well organized. Figure 4 illustrates a simplified scheme of the necessary data flow. It is believed that the

Materials and Methods
Using multi-cell steady-state model of a HE, a control-theory based approach was proposed for the identification and evaluation of the influence of fouling on the dynamic behavior of the heat exchangers and on the quality of their control. A prerequisite for reliable monitoring of the quality of HEN control under fouling conditions is that acquisition and processing of operation data are well organized. Figure 4 illustrates a simplified scheme of the necessary data flow. It is believed that the proposed approach can be applied to the HENs used in continuously operated process plants of oil refining, chemical, food processing and other industries. This is illustrated by a case study in which HENs coupled with crude distillation units are investigated.
The scheme of data flow and its details needed for the minimization of uncertainty margins of the monitoring of HEN control are presented below in four stages.

Stage 1 -Acquisition and Pre-processing of the Data Base
A prerequisite for identification of the effect of HE fouling on the quality of HEN control is the availability of geometrical data of all the relevant HEs -both those included in the control loop and other ones that may interact with loop components. Equally important is continuous availability of the data on HEs operation, that is, mass flows, temperatures and chemical compositions of the involved process streams, that are necessary for the determination of physico-chemical properties of

Stage 1-Acquisition and Pre-processing of the Data Base
A prerequisite for identification of the effect of HE fouling on the quality of HEN control is the availability of geometrical data of all the relevant HEs-both those included in the control loop and other ones that may interact with loop components. Equally important is continuous availability of the data on HEs operation, that is, mass flows, temperatures and chemical compositions of the involved process streams, that are necessary for the determination of physico-chemical properties of the media flowing through HEN components. However, as raw process data may also reflect inaccurate measurements, errors in data transmission and recording, as well transient states of the HEN, it is necessary to apply data pre-processing by filtering, averaging and reconciliation. While filtering is aimed at the elimination of gross errors in the recorded data, averaging (over representative time intervals) is needed for the determination of parameter values that enable application of mathematical models of steady-state heat transfer. Finally, data reconciliation makes it possible to minimize uncertainties induced by measurement errors and deviations from steady state of the HEN. Appropriate methods of data filtering, averaging and reconciliation are presented elsewhere [13,14].

Stage 2-Identification of fouling effects in HEs in the steady state
For each heat exchanger in the studied HEN, the data base established in Stage 1 is used for determining the characteristics of fouling. The existence of fouling and its time behavior in the HE (see example in Figure 5) are represented by the evolution of coefficient R f of the total thermal resistance of fouling layers on both sides of the heat transfer surface. Parameter R f is calculated as the difference between thermal resistances of fouled and clean heat transfer surface.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 6 of 24 application of mathematical models of steady-state heat transfer. Finally, data reconciliation makes it possible to minimize uncertainties induced by measurement errors and deviations from steady state of the HEN. Appropriate methods of data filtering, averaging and reconciliation are presented elsewhere [13,14].

Stage 2 -Identification of fouling effects in HEs in the steady state
For each heat exchanger in the studied HEN, the data base established in Stage 1 is used for determining the characteristics of fouling. The existence of fouling and its time behavior in the HE (see example in Figure 5) are represented by the evolution of coefficient Rf of the total thermal resistance of fouling layers on both sides of the heat transfer surface. Parameter Rf is calculated as the difference between thermal resistances of fouled and clean heat transfer surface. The mathematical model includes widely known relationships describing heat transfer phenomena and energy balance in the heat exchanger, as well as multi-cell representation of steadystate operation of shell-and-tube HE [15]. It also includes the algorithm of least-squares based determination of the heat transfer coefficient that has been presented in earlier publications by the present authors [16,17]. Details of the elaborated model are valid for shell-and-tube HEs only, but by adapting the relationships describing the heat transfer phenomena and energy balance for other types of heat exchangers (e.g. spiral or plate HEs), the model can be generalized.

Stage 3 -Development and Validation of the Dynamic HE Model
The planning of efficient use of HEs under changing operating conditions (e.g., conditions resulting from fouling build-up with time) requires the application of adequate dynamic models. Stage 3 is based on the mathematical model proposed by Trafczynski et al. [6], of transient heat exchange with the influence of thermal resistance of fouling taken into account.
According to the scheme shown in Figure 4, the three main steps of Stage 3 are: • Determination of a cell-based dynamic HE model based on the operator transmittances. By solving the equations of the mathematical model, relationships employing operator transmittances can be obtained between disturbances occurring at cell inlet and changes in temperature at the cell outlet. Operator transmittance G(s) is a widely used tool for describing a dynamic system. This step is extensively described in Section 2 of the authors' previous work The mathematical model includes widely known relationships describing heat transfer phenomena and energy balance in the heat exchanger, as well as multi-cell representation of steady-state operation of shell-and-tube HE [15]. It also includes the algorithm of least-squares based determination of the heat transfer coefficient that has been presented in earlier publications by the present authors [16,17]. Details of the elaborated model are valid for shell-and-tube HEs only, but by adapting the relationships describing the heat transfer phenomena and energy balance for other types of heat exchangers (e.g. spiral or plate HEs), the model can be generalized.

Stage 3-Development and Validation of the Dynamic HE Model
The planning of efficient use of HEs under changing operating conditions (e.g., conditions resulting from fouling build-up with time) requires the application of adequate dynamic models. Stage 3 is based on the mathematical model proposed by Trafczynski et al. [6], of transient heat exchange with the influence of thermal resistance of fouling taken into account.
According to the scheme shown in Figure 4, the three main steps of Stage 3 are: • Determination of a cell-based dynamic HE model based on the operator transmittances. By solving the equations of the mathematical model, relationships employing operator transmittances can be obtained between disturbances occurring at cell inlet and changes in temperature at the cell outlet.
Operator transmittance G(s) is a widely used tool for describing a dynamic system. This step is extensively described in Section 2 of the authors' previous work [6]. • Implementation of the dynamic HE model. Starting from HEN block diagram in which the role of operator transmittances was visualized and using MATLAB/Simulink program package, a software module was developed to simulate the performance of HEN control. In order to make simulation possible, a database is needed for providing the values of relevant parameters in all the cells at steady state (from stage 2), of all the HEs in the HEN. This step is extensively described in Section 2.2 of the authors' previous work [6].

•
Validation of the dynamic HE model using operational data of a real-life HEN coupled with a CDU. The values of simulated and real temperature at heat exchanger outlet in transient states of the exchanger were compared and found to be in close agreement. This step is in more detail described in Section 3 of the authors' previous work [6].

Stage 4-Identification and Evaluation of the Influence of Fouling on the Dynamic Behavior of PID-Controlled HEs and on the Control-Quality Indices
The first step is to study the open-loop unit step responses simulated at different periods of fouling build-up for all controlled HEs.
When parameter R f is increased, the thermal inertia of the HE is changed, leading to changes in its dynamic behavior. Such changes can be detected by studying the open-loop unit step responses simulated at different stages of fouling build-up. A typical response of a HE system is illustrated in Figure 6a. In Figure 6b, open-loop responses of the HE model are plotted for a step upset +5% in the shell-side flowrate M s . As can be seen in these responses, fouling build-up on the exchanger's heat transfer surface induces changes in the values of gain K o , delay time t d and time constant t 1 .
Appl. Sci. 2019, 9, x FOR PEER REVIEW 7 of 24 • Validation of the dynamic HE model using operational data of a real-life HEN coupled with a CDU. The values of simulated and real temperature at heat exchanger outlet in transient states of the exchanger were compared and found to be in close agreement. This step is in more detail described in Section 3 of the authors' previous work [6].

Stage 4 -Identification and Evaluation of the Influence of Fouling on the Dynamic Behavior of PID-Controlled HEs and on the Control-Quality Indices
The first step is to study the open-loop unit step responses simulated at different periods of fouling build-up for all controlled HEs.
When parameter Rf is increased, the thermal inertia of the HE is changed, leading to changes in its dynamic behavior. Such changes can be detected by studying the open-loop unit step responses simulated at different stages of fouling build-up. A typical response of a HE system is illustrated in Figure 6a. In Figure 6b, open-loop responses of the HE model are plotted for a step upset +5% in the shell-side flowrate Ms. As can be seen in these responses, fouling build-up on the exchanger's heat transfer surface induces changes in the values of gain Ko, delay time td and time constant t1. The second step is modeling of a control-loop unit including the dynamic exchanger model together with PID control configurations in Matlab-Simulink (see Figure 7a). In the next two steps, assuming PID control of the exchanger unit, three gain coefficients (proportional Kp, integral Ki and derivative Kd) are needed to determine its closed-loop characteristics. The values of gain coefficients can be determined by the Ziegler-Nichols method [18], which is commonly used in industry. The second step is modeling of a control-loop unit including the dynamic exchanger model together with PID control configurations in Matlab-Simulink (see Figure 7a). In the next two steps, assuming PID control of the exchanger unit, three gain coefficients (proportional K p , integral K i and derivative K d ) are needed to determine its closed-loop characteristics. The values of gain coefficients can be determined by the Ziegler-Nichols method [18], which is commonly used in industry.
The second step is modeling of a control-loop unit including the dynamic exchanger model together with PID control configurations in Matlab-Simulink (see Figure 7a). In the next two steps, assuming PID control of the exchanger unit, three gain coefficients (proportional Kp, integral Ki and derivative Kd) are needed to determine its closed-loop characteristics. The values of gain coefficients can be determined by the Ziegler-Nichols method [18], which is commonly used in industry. For specific dynamic characteristics similar to that shown in Figure 7b, description and evaluation of the control quality can be based on the quality indices including:

•
Overshoot Mp -percentage of the maximum deviation of step response y(t) from its steady-state value: where: ymax -y(t) at its maximum, yss -y(t) at steady state (yss ≤ ymax). The value of overshoot is determined during control-system design and may be used as a measure of system stability; large overshoot values are not recommended.
• Peak time tp -time interval to the maximum of y(t), that is, y(tp) = ymax. • Delay time td -time interval to the step response reaching 50% of its value at steady state, that is, y(td) = 0.5 yss.

•
Rise time tr -time interval to the step response reaching 80% of its value at steady state, that is, y(tr) = 0.8 yss.

•
Settling time ts -time interval to the step response staying within the tolerance margin of its steady-state value, usually yss ± 5% (see Figure 7b).
Being easy to determine, the abovementioned quality indices can be used to evaluate the characteristics of the control system on the basis of its response to step changes of process variables. The dynamic characteristics were simulated at different stages of fouling build-up, that is, after one, two and three years of the continuous operation of the HE unit. Initially, transient responses of an exemplary HE unit were simulated assuming PID control with constant values of the gain coefficients that were determined for clean heat exchange surface (that is, without fouling). As can be seen in the responses obtained for the consecutive periods of HE operation (Figure 7c), the build-up of fouling and the increased thermal resistance would lead to oscillations of the controlled temperature, a too slow response to set-point changes and the risk of significant temperature overshoot that may be dangerous especially during the execution of start-up procedures. However, adverse changes in control quality can be prevented by periodic adjustments of the gain coefficients. This can be seen in Figure 7d, which depicts simulated step responses of an exemplary HE with controller tuning parameters adjusted for the consecutive periods of fouling build-up. These characteristics indicate that if the real-life controller tuning was adjusted to fit the requirements of efficient control, then despite increased values of the thermal resistance, the indices of control quality would not be For specific dynamic characteristics similar to that shown in Figure 7b, description and evaluation of the control quality can be based on the quality indices including:

•
Overshoot M p -percentage of the maximum deviation of step response y(t) from its steady-state value: where: y max − y(t) at its maximum, y ss − y(t) at steady state (y ss ≤ y max ). The value of overshoot is determined during control-system design and may be used as a measure of system stability; large overshoot values are not recommended. • Peak time t p -time interval to the maximum of y(t), that is, y(t p ) = y max . • Delay time t d -time interval to the step response reaching 50% of its value at steady state, that is, y(t d ) = 0.5 y ss . • Rise time t r -time interval to the step response reaching 80% of its value at steady state, that is, y(t r ) = 0.8 y ss . • Settling time t s -time interval to the step response staying within the tolerance margin of its steady-state value, usually y ss ± 5% (see Figure 7b).
Being easy to determine, the abovementioned quality indices can be used to evaluate the characteristics of the control system on the basis of its response to step changes of process variables. The dynamic characteristics were simulated at different stages of fouling build-up, that is, after one, two and three years of the continuous operation of the HE unit. Initially, transient responses of an exemplary HE unit were simulated assuming PID control with constant values of the gain coefficients that were determined for clean heat exchange surface (that is, without fouling). As can be seen in the responses obtained for the consecutive periods of HE operation (Figure 7c), the build-up of fouling and the increased thermal resistance would lead to oscillations of the controlled temperature, a too slow response to set-point changes and the risk of significant temperature overshoot that may be dangerous especially during the execution of start-up procedures. However, adverse changes in control quality can be prevented by periodic adjustments of the gain coefficients. This can be seen in Figure 7d, which depicts simulated step responses of an exemplary HE with controller tuning parameters adjusted for the consecutive periods of fouling build-up. These characteristics indicate that if the real-life controller tuning was adjusted to fit the requirements of efficient control, then despite increased values of the thermal resistance, the indices of control quality would not be adversely affected.
Overall, the presented results for control loop of a HE ( Figure 7) indicate that if the thermal resistance of fouling is increased, unchanged parameters of controller tuning could lead to the deterioration of the indices of control quality.
By adjusting the values of proportional-integral-derivative gains K p , K i , K d , these adverse effects of fouling could be prevented. For a given value of the thermal resistance of fouling, appropriate gain values could be determined using the dynamic model of the heat exchanger and the suitability of these values can be tested by simulation -which is the last step in the Stage 4 of the proposed procedure (see Figure 4).

Case Study-Results
In order to investigate the influence of fouling build-up on the dynamic behavior of the HENs and on the quality of their control, two cases (the real-life HENs coupled with a CDU plant-see Figures 2 and 3) were considered. Fractional distillation of crude oil is a highly energy-intensive process that requires the crude to be heated from ambient temperature to around 370 • C. The required heat is provided through a set of HEs in which heat from the distillation products and pump-around streams of the distillation columns is recovered, and a furnace fuelled by heavy fuel oil. The crude is pumped through the first part of the HEN to a desalting unit where it is washed with water to remove inorganic water-soluble impurities. After that, the crude flows through the second HEN part, and further to the furnace where it is heated up to the temperature needed for entering the fractional distillation column.
Using operational data available from the period of three years of continuous HEN operation, exchanger characteristics were studied at different stages of fouling build-up, that is, after 1, 2 and 3 years (passed from operation start-up when HE surfaces had been clean).

Case No. 1
In case no. 1 four branches ABCD (the crude preheat trains) were selected from a real-life HEN coupled with a CDU rated 110 kg/s of crude oil. Twenty-six shell-and-tube, two-pass HEs with straight tubes and floating heads are connected as schematically shown in Figures 2 and 8.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 9 of 24 model of the heat exchanger and the suitability of these values can be tested by simulation -which is the last step in the Stage 4 of the proposed procedure (see Figure 4).

Case study -Results
In order to investigate the influence of fouling build-up on the dynamic behavior of the HENs and on the quality of their control, two cases (the real-life HENs coupled with a CDU plant -see Figures 2 and 3) were considered. Fractional distillation of crude oil is a highly energy-intensive process that requires the crude to be heated from ambient temperature to around 370°C. The required heat is provided through a set of HEs in which heat from the distillation products and pump-around streams of the distillation columns is recovered, and a furnace fuelled by heavy fuel oil. The crude is pumped through the first part of the HEN to a desalting unit where it is washed with water to remove inorganic water-soluble impurities. After that, the crude flows through the second HEN part, and further to the furnace where it is heated up to the temperature needed for entering the fractional distillation column.
Using operational data available from the period of three years of continuous HEN operation, exchanger characteristics were studied at different stages of fouling build-up, that is, after 1, 2 and 3 years (passed from operation start-up when HE surfaces had been clean).

Case No. 1
In case no. 1 four branches ABCD (the crude preheat trains) were selected from a real-life HEN coupled with a CDU rated 110 kg/s of crude oil. Twenty-six shell-and-tube, two-pass HEs with straight tubes and floating heads are connected as schematically shown in Figures 2 and 8. Owing to limited measurement data, it was not possible to determine the relationship between the thermal resistance of fouling Rf (fouling factor) and time t, for each HE. The measurements of temperature and mass flow were performed only at the inlet and outlet of the studied HEN but no temperature measurements were available between the HEs. In order to resolve this issue, the Rf values of HEs that were used in the simulation studies had been postulated by the authors on the Owing to limited measurement data, it was not possible to determine the relationship between the thermal resistance of fouling R f (fouling factor) and time t, for each HE. The measurements of temperature and mass flow were performed only at the inlet and outlet of the studied HEN but no temperature measurements were available between the HEs. In order to resolve this issue, the R f values of HEs that were used in the simulation studies had been postulated by the authors on the basis of values recommended by TEMA standards [19] (see Table 1). As demonstrated in reference [6], such thermal resistance values may significantly affect the performance of HE control. As shown in the HEN scheme in Figure 2, the crude-oil feed stream is split in parallel branches A and B before the desalting unit and in parallel branches C and D after the desalting unit. For the two control loops with PID controllers 1 and 2 as indicated in Figure 2, the split ratios in branch pairs AB and CD are adopted as manipulated variables (which can be changed by the action of control valves 1AB and 2AB), while the controlled variables are defined as the differences between the studied outlet temperatures: The control objective is to maximize heat recovery, understood as total heat flow Q transferred in the HEN, and the setpoint values of the controlled variables should be CV 1 = CV 2 = 0. In other words, when process disturbances occur, the controllers installed in the HEN are required to adjust the split ratios in network branches AB and CD to ensure that the values of the controlled variables return to zero.

Case No. 2
In case no. 2 the crude preheat trains were selected from another real-life HEN coupled with a CDU rated 220 kg/s of crude oil. Fourteen shell-and-tube, two-pass HEs with straight tubes and floating heads are connected as schematically shown in Figures 3 and 9.
For each HE, the relationship between the fouling factor R f and time t, was determined using method described in the work [16]. In this case, all measurements of temperature and mass flow at the inlet and outlet of the studied HEs were available. Obtained R f values of HEs that were used in the simulation studies are presented in Table 2.
As shown in the HEN scheme in Figure 3, the desalted crude-oil stream is split in parallel branches and after exchangers E3-14AB and E3-16ABC the branches are connected again into the one preheat train. There is one simple control setup with PID controller 3. In the control loop 3 with exchanger E3-18AB as indicated in Figure 3, the controlled variable is the tube-side outlet temperature T 3 before the preflash column, while the manipulated variable is shell-side by-pass mass flow rate (which can be changed by the action of control valve 3). The other process variables are the disturbances.   For the different periods of HEN operation during which fouling was building up, simulations have been carried out in Simulink. According to the obtained results, when the thermal resistance of fouling is increased, the thermal inertia of every HE is changed leading to changes in the dynamic behavior of the interacting A, B, C, D branches shown in Figure 2. Such changes can be detected by studying the open-loop step responses of the end temperatures (after parallel branches TABend, TCDend)   For the different periods of HEN operation during which fouling was building up, simulations have been carried out in Simulink. According to the obtained results, when the thermal resistance of fouling is increased, the thermal inertia of every HE is changed leading to changes in the dynamic behavior of the interacting A, B, C, D branches shown in Figure 2. Such changes can be detected by studying the open-loop step responses of the end temperatures (after parallel branches T ABend , T CDend ) simulated at the different stages of fouling build-up. The features of a typical response of a heat exchanger system (network branch) are illustrated in Figure 6a.
In Figure 10a,b, simulated open-loop responses of the studied HEN models are plotted for a step upset +1 • C in the tube-side inlet temperatures of the branches A (exchanger E1-11AB) and B (exchanger E1- 21). The open-loop responses of the branches to a step upset +1 • C in the shell-side inlet temperatures of the exchangers (E1-13, 14, 27 and E2-14, 15AB, 25AB, 26, 27), are presented in Figure 11a,b. Next, Figure 12a        PID controllers for loops 1 and 2 (see Figure 2) were separately tuned according to the Skogestad tuning rules [20], by assuming step (10%) increases in the crude oil mass flows (+6.11 kg/s in MAt and MCt) in each of the branches A and C. The control variable responses for each of the selected operation periods are shown in Figure 14a,b and the resulting values of the tuning parameters for PID controllers 1 and 2 are presented in Table 3.  PID controllers for loops 1 and 2 (see Figure 2) were separately tuned according to the Skogestad tuning rules [20], by assuming step (10%) increases in the crude oil mass flows (+6.11 kg/s in MAt and MCt) in each of the branches A and C. The control variable responses for each of the selected operation periods are shown in Figure 14a,b and the resulting values of the tuning parameters for PID controllers 1 and 2 are presented in Table 3. In order to prevent that from happening, it is advisable to investigate all the three components of the tuning of each PID controller (K p , K i , K d ) that is, gain values in the proportional, integral and derivative components) and to check the resulting transient responses.
PID controllers for loops 1 and 2 (see Figure 2) were separately tuned according to the Skogestad tuning rules [20], by assuming step (10%) increases in the crude oil mass flows (+6.11 kg/s in M At and M Ct ) in each of the branches A and C. The control variable responses for each of the selected operation periods are shown in Figure 14a,b and the resulting values of the tuning parameters for PID controllers 1 and 2 are presented in Table 3.  In case no. 2, the changes in the dynamic behavior of the E3-18AB HE unit operated in HEN (see Figure 3) can be detected by studying the open-loop step responses of the outlet temperature T3 simulated at the different stages of fouling build-up. Figure 15a shows the open-loop responses of the studied outlet temperature T3 to +1% step change in the shell-side flowrate (+0.47 kg/s) of the HE. Next, Figure 15b depicts the open-loop responses of the studied outlet temperature T3 to -1% step change in the tube-side flowrate (-0.73 kg/s) of the HE. In Figure 15c, simulated open-loop responses of the studied outlet temperature T3 are plotted for a step upset +1°C in the shell-side inlet temperature of the HE. Finally, the open-loop responses of the studied outlet temperature T3 to a step upset +1°C in the tube-side inlet temperature of the exchanger E3-18AB, are presented in Figure 15d.  In case no. 2, the changes in the dynamic behavior of the E3-18AB HE unit operated in HEN (see Figure 3) can be detected by studying the open-loop step responses of the outlet temperature T 3 simulated at the different stages of fouling build-up. Figure 15a shows the open-loop responses of the studied outlet temperature T 3 to +1% step change in the shell-side flowrate (+0.47 kg/s) of the HE. Next, Figure 15b depicts the open-loop responses of the studied outlet temperature T 3 to −1% step change in the tube-side flowrate (−0.73 kg/s) of the HE. In Figure 15c, simulated open-loop responses of the studied outlet temperature T 3 are plotted for a step upset +1 • C in the shell-side inlet temperature of the HE. Finally, the open-loop responses of the studied outlet temperature T 3 to a step upset +1 • C in the tube-side inlet temperature of the exchanger E3-18AB, are presented in Figure 15d. PID controller for loop 3 (see Figure 3) was tuned according to the Ziegler-Nichols method [18]. The values of the parameters of the open-loop characteristics (shown in Figure 15a) for each of the selected operation periods and the resulting values of the tuning parameters for PID controller 3 are presented in Table 4.  Because of the fouling build-up on the exchangers' heat transfer surfaces, the visible changes in the delay time t d , time constant t 1 and gain K 0 values in the studied open-loop step responses, may impair the quality of PID control considerably. In order to prevent that from happening, it is also advisable to investigate the components of the tuning PID controller 3 (K p , K i , K d ) and to check the resulting transient responses.
PID controller for loop 3 (see Figure 3) was tuned according to the Ziegler-Nichols method [18]. The values of the parameters of the open-loop characteristics (shown in Figure 15a) for each of the selected operation periods and the resulting values of the tuning parameters for PID controller 3 are presented in Table 4.

Closed-loop Control Analysis
Using the dynamic HE model outlined in Section 2.3, the entire HENs together with control configurations were modelled employing Simulink software; the block diagram of the HEN model is presented in Figures 8 and 9.
For case no. 1, the simulations of transient responses were carried out and their results shown in Figure 16a,b (for control loops 1 and 2) demonstrate that fouling build-up induces insignificant changes in CV rise time t r and settling time t s . This can be seen as an indication that in the studied case, no adjustments of K p and K i values are needed and the indices of control quality would not be adversely affected by fouling of heat-exchanger surfaces-see Table 5.

Closed-loop Control Analysis
Using the dynamic HE model outlined in Section 2.3, the entire HENs together with control configurations were modelled employing Simulink software; the block diagram of the HEN model is presented in Figures 8 and 9.
For case no. 1, the simulations of transient responses were carried out and their results shown in Figure 16a,b (for control loops 1 and 2) demonstrate that fouling build-up induces insignificant changes in CV rise time tr and settling time ts. This can be seen as an indication that in the studied case, no adjustments of Kp and Ki values are needed and the indices of control quality would not be adversely affected by fouling of heat-exchanger surfaces -see Table 5.

Closed-Loop Responses for
Case No. 1

R f after Period of Iperation (years)
Base PID Parameters 1 (K p /K i /K d )  Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters.   Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters.  For case no. 2 (HEN with control loop 3), the closed-loop step responses under fouling conditions were simulated with three different sets of the PID parameters:

1.
Assuming constant values of the base PID parameters obtained using Ziegler-Nichols method [18] for clean HEN- Figure 17a  2.
With the adjusted PID parameters, for the consecutive periods of fouling build-up, in accordance with the data shown in Table 4- Figure 17b 3.
With the optimal PID parameters obtained using Signal Constraint toolbox in SIMULINK [21] under fouling conditions- Figure 17c In case no. 2, control loop 3 comprises heat exchangers E3-18AB whose operation is affected by the interactions with the remaining exchangers in the studied HEN (see Figure 3). Qualitative evaluation of the obtained dynamic closed-loop characteristics (Figure 17a-c) can be complemented by the values of quality indices- Table 6.      Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. 1 characteristics were simulated assuming constant values of the PID parameters obtained using Ziegler-Nichols method [18] for clean HEN 2 characteristics with the adjusted PID parameters, for the consecutive periods of fouling build-up, in accordance with the data shown in Table 4 3 characteristics with the optimal PID parameters obtained under fouling conditions using Signal Constraint toolbox in Simulink [21] In this case, judging from the information presented in Figure 17 and Table 6, the effect of fouling that builds up during HEN operation is more pronounced than that observed in case no. 1. At unchanged PID-controller settings, settling time t s is increased from 124.6 s for the clean exchanger to 173.2 s for the exchanger fouled after 3 years of HEN operation. While some quality indices including peak time t p , delay time t d and rise time t r , are changed, overshoot M p remains nearly constant (Figure 17a). If adjustments of PID-controller settings were applied in reaction to a fouling build-up, then the resulting dynamic characteristics and control-quality indices would not deteriorate (Figure 17b and Table 6). Using the Ziegler-Nichols method [18] or the Signal Constraint toolbox in Simulink [21], optimal controller settings can be determined to eliminate excessively large values of overshoot M p (Figure 17c). In this context, a more advanced approach recently introduced by Oravec [12] in cooperation with Trafczynski and Markowski can be mentioned. In their work, robust Model Predictive Control -MPC with integral action is used for optimizing the control performance when the operation of heat exchangers has been affected by fouling that induces changes of the exchangers' parameters.

Discussion
Using Simulink software, a validated multi-cell dynamic model of a shell-and-tube HE was applied in simulating the operation of PID-controlled HEs (see Figures 8 and 9). A control-theory based approach was proposed for the identification and evaluation of the influence of fouling on the dynamic behavior of the HEN and on the quality of its control (Figure 4). The dynamic model was applied to a case study on the HEs and HENs operated in the crude distillation unit under fouling conditions.
In case no. 1, control loops no. 1 and 2 include all the HEs operated in the HEN and interactions between the HEs are significant (Figure 2). The simulated step responses prove that as fouling was building up, the quality indices of network control remained nearly unchanged even if the tuning of PI controllers was not adjusted (see Figure 16 and Table 5).
In case no. 2, control loop no. 3 includes a set of HEs that interact with other exchangers present in the HEN (Figure 3). From the qualitative and quantitative estimates presented above, it can be inferred that the effect of fouling on HEN operation is more pronounced than that observed in case no. 1. Although most indices of control quality remain unchanged as the fouling increases, the settling time becomes longer. Periodic adjustments of PID controller tuning are required in the consecutive stages of fouling build-up because the value of parameter t s is significantly increased (see Figure 17 and Table 6).
In previous publications [6,21], the present authors evaluated dynamic characteristics of four different sets of heat exchangers operated in simple control loops according to the scheme shown in Figure 7a, that is, without meaningful interactions with other exchangers in the HEN. The performance of those control loops have now been simulated and their closed-loop step responses are presented in Appendix A, Figures A1-A4, while the corresponding values of control quality indices are shown in Table A1. In each of the studied PID-controlled HEs, it was found that increased fouling led to the deterioration of all the indices of control quality so that periodic adjustments of PID-controller tuning appeared necessary.
It can be mentioned that the dynamic model of shell-and-tube heat exchangers developed by the present authors found application in the work done by Borges de Carvalho et al. [22], who performed the dynamic analysis of fouling build-up in the HEs designed according to TEMA standards.
The same author team also tested several tuning strategies for the PID-controlled HEs under fouling conditions [23] and arrived at results that appear to be consistent with those of the present authors.
Overall, according to the results of the mentioned case studies, the higher the number of heat exchangers in the PID control loop and the more interactions occur between heat exchangers in the network, the weaker the influence of fouling on the control quality indices (see Table 7). This observation may be attributed to underestimated values of R f (calculated according to TEMA standards) and/or to the compensation of the negative impacts of fouling on the heat transfer in the HEs. Such a compensation is possible only in the network where significant interactions occur between the HEs (that is, if antecedent exchangers are operated on the process streams- Figure 18). As previously observed by the present authors [17,24], the larger the number of interacting exchangers, the better the compensation of the detrimental effects of fouling. Fouling on the heat transfer surface of a HE operated in the HEN brings about a change in the exchanger capacity as well as changes in the outlet temperatures of process streams. However, the operation of the HE can be affected by other exchangers serving the same process streams (antecedent exchangers); examples of such exchangers in the HENs can be found in Figures 2 and 3. As fouling builds up on the heat transfer surfaces of the antecedent exchangers, temperatures of process streams at HE inlet are increased. Due to that, although heat transfer intensity has been reduced by fouling, the thermal power of the HE may remain unchanged. Similarly, the indices of control quality in the associated control loop may also remain unchanged (Figure 16a,b). Table 7. The number of HEs in a control loop and the number of antecedent HEs compared to the level of fouling influence on the control-quality indices.

The Number of Antecedent HEs
The Level of Fouling Influence on the Control-quality Indices. on Hot Stream on Cold Stream Case no. 1 (Figure 2) (see Table 5 Table 6) for control loop 3 1 7 7 low (Figure 17a) Cases in work [6] ( Figure 18) (see Table A1) for HE11 control loop (see Table 6) low (Figure 17a Table A1

Conclusions
In conclusion, insufficient quality of HEN control may lead to excessive oscillations (increased settling time and overshoot) of process parameters, as well as to excessive consumption of energy and raw materials, resulting in increased production costs. It may also generate the risk of dangerous process perturbations such as, exceeding safety margins of temperature values. In order to prevent such situations from occurring, various approaches to the determination of controller-tuning parameters can be applied such as the trial-and-error procedure, the Ziegler-Nichols method, or the MPC methodology. Appropriate PID-gain values can be determined using the dynamic model of the heat exchanger network and the suitability of these values can be tested by simulation. In order to ensure a satisfactory performance of PID control when fouling layers build up on the heat-transfer surfaces of the exchangers in the HEN, periodic adjustments of PID-controller tuning are needed. A more costly alternative is to apply periodic cleaning of the exchangers. Where the rate of fouling build up is very high and therefore exchanger cleaning cannot be avoided, the adjustments of controller tuning may help to reduce the frequency of cleaning interventions, thus lowering their total cost.

Conclusions
In conclusion, insufficient quality of HEN control may lead to excessive oscillations (increased settling time and overshoot) of process parameters, as well as to excessive consumption of energy and raw materials, resulting in increased production costs. It may also generate the risk of dangerous process perturbations such as, exceeding safety margins of temperature values. In order to prevent such situations from occurring, various approaches to the determination of controller-tuning parameters can be applied such as the trial-and-error procedure, the Ziegler-Nichols method, or the MPC methodology. Appropriate PID-gain values can be determined using the dynamic model of the heat exchanger network and the suitability of these values can be tested by simulation. In order to ensure a satisfactory performance of PID control when fouling layers build up on the heat-transfer surfaces of the exchangers in the HEN, periodic adjustments of PID-controller tuning are needed. A more costly alternative is to apply periodic cleaning of the exchangers. Where the rate of fouling build up is very high and therefore exchanger cleaning cannot be avoided, the adjustments of controller tuning may help to reduce the frequency of cleaning interventions, thus lowering their total cost.

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

Appendix A
It can be mentioned that in previous publications [6,21], the present authors qualitatively evaluated dynamic characteristics of some other components of the studied HEN, namely heat exchangers E11AB, E15AB, E30AB and E35AB-see Figures A1-A4. These HEs were assumed to operate in the control loops similar to that shown in Figure 7a in Section 2.4., that is, in the absence of meaningful interactions with other HEN components. As a complement to the mentioned characteristics, the corresponding values of control-quality indices are shown in Table A1.
Regarding exchanger E15AB, the evaluation of its control performance is similar to that discussed for case no. 2 (exchanger E3-18AB) in Section 3.4. As a consequence of fouling build-up at unchanged PID-controller settings ( Figure A2a), settling time t s increases from 70.9 s to 145 s, for the clean HE and fouled HE after 3-year operation (Table A1). Concurrently, peak time t p , delay time t d and rise time t r are slightly changed, while overshoot M p remains nearly constant.
The indices of control quality of the other HEs affected by fouling build-up are generally deteriorated, however the extent of change is differentiated. At constant PID-controller settings, the indices of E11AB are changed as follows: overshoot M p is increased from 47.3% to 73.8%, peak time t p from 18.7 s do 31.3 s, and settling time t s from 61.3 s to >300 s, for the clean HE and fouled HE after 3-year operation, respectively. The remaining indices, that is, delay time t d and rise time t r are changed insignificantly (see Table A1 and Figure A1a). Analogous changes determined for E30AB are: overshoot M p from 45.5% to 54.3%, and settling time t s from 76.4 s to 152 s, while peak time t p , delay time t d and rise time t r are nearly unchanged (Table A1 and Figure A3a). The situation of heat exchanger E35AB is different because at constant controller settings, fouling build-up may lead to drastic deterioration of control-quality indices and unstable control performance (see Figure A4a and Table A1). However, if periodic adjustment of the settings of PID-controllers were applied for all the mentioned heat exchangers (E11AB, E15AB, E30AB and E35AB), then increasing thermal resistances of the fouling layers would not induce deterioration of control characteristics (Figures A1b, A2b, A3b and A4b) and their corresponding quality indices (Table A1). The adjustments of controller settings, optimized using Ziegler-Nichols method [18] or Signal Constraint toolbox in Simulink [21],    Figure A2. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E15AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A2. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E15AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A2. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E15AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A3. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E30AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A3. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E30AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A2. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E15AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A3. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E30AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Table A1. Values of the control-quality indices for used sets of the PID parameters in different periods of operation (cases in work [6]). Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Table A1. Values of the control-quality indices for used sets of the PID parameters in different periods of operation (cases in work [6]).  Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters.    Figure A4. Closed-loop step responses (case in work [6]) under fouling conditions of the HE E35AB models: without adjustment (a) after adjustment (b) and after adjustment of optimal (c) PID tuning parameters. Table A1. Values of the control-quality indices for used sets of the PID parameters in different periods of operation (cases in work [6]).