A Novel Sequential Approach for the Design of Heat Exchanger Networks

This paper presents a new algorithm for the design of heat exchanger networks (HEN) that tries to take advantage of the strengths of the sequential and simultaneous approaches. It is divided into two sequential parts. The first one is an adaptation of the transportation model (TransHEN). It maintains the concept of temperature intervals and considers the possibility of heat transfer between all the hot and cold streams inside those intervals, and at the same time it allows the a priori calculation of the logarithmic mean temperature difference between all possible heat exchanges, and therefore it maintains the area estimation linear in the model. The second step (HENDesign model), uses a superstructure that contains all the possible alternatives in which the matches predicted by the first stage model can exchange heat to design the final heat exchanger network. Unlike the sequential approach, in this model, all heat flows, temperatures, areas, etc. are reoptimized maintaining the set of matches predicted in the first stage. The model is highly nonlinear and nonconvex, however, it is relatively easy to get good results, because the model starts with the values predicted by the TransHEN model. The algorithm has been tested using fifteen benchmark problems commonly used in literature to compare the performance of heat exchanger network algorithms. In eleven out of the fifteen cases present better or equal results than the best ones reported in the open literature. In three the results presented only marginal differences in total annualized cost (lower than 0.5%) and only a difference of 2.4% in the largest one.


INTRODUCTION
In a worldwide ever-growing energy consumption scenario, heat integration has, since the 1970's, risen as a topic of major interest in process engineering. Efficient energy management in processing plants may lead to more profitability and more environmentally friendly industry. To achieve those goals, facilities with fluid material streams must consider heat exchange among these streams. Such a task is performed indirectly by heat exchangers. Proper allocation of these units includes selecting what streams are passing through them for a hot one to provide heat to a cold one. This mitigates the need for external utilities, which are also performed in heat exchangers identified as heaters or coolers. In heaters, a hot utility (e.g.,: steam from a boiler) is used to provide heat to a cold process stream, while in coolers, a cold utility (e.g.,: cold water from cooling towers) is employed as a cooling mean. Note that utilities production yield operating costs associated, for instance, to fuel consumption in boilers or electricity to run cooling towers or chiller cycles.
The set of heat exchangers implemented in a plant is named heat exchanger network (HEN). HEN synthesis is an intriguing area under the scope of industrial energy management. That is since in general the mathematics involved in describing mass/ energy balances, heat transfer, and thermodynamic phenomena in HEN, even when simplified, yield models that are non-linear and non-convex. Moreover, as the number of streams in a plant grows, the number of possible heat exchanger matches becomes inconveniently large and turns out that exhaustively evaluating these combinations is computationally impossible with current technology. For large-scale cases, a globally optimal solution is currently impossible. As pointed out by Furman and Sahinidis (2001), these aspects result in an NP-hard problem in the strong sense.
The development of strategies to overcome these complicating characteristics methodologically is thus of great value. A pioneering method is the noteworthy pinch analysis (Linnhoff and Flower, 1978;Flower and Linnhoff, 1980). Pinch analysis proposed systematic procedures to find energy consumption targets for a given set of process streams based on thermodynamic insights. Moreover, simple recommendations were also given for the allocation of units by dividing the network into two according to temperature intervals. Later, apart from energy-related targets, Linnhoff's research group was able to include goals concerning total heat transfer area. This means predicting the minimum area that a network would have in near-optimal solutions. According to the proposing works, these energy and area predictions could lead to a priori knowledge of the total annual costs (capital + operating costs) of a network with a 5% accuracy (Linnhoff and Ahmad, 1990). Concepts proposed by the time were expanded and are used to the day in several works. A detailed review on pinch analysis was performed by Klemeš et al. (2018) where more than 300 works that to some degree are based on pinch are cited.
Pinch analysis has been extended to other areas besides heat integration. It is remarkable the parallelism in concepts and methodology between heat integration and mass integration (MEN -mass exchanger networks-) (El-Halwagi and Manousiouthakis, 1989;Manousiouthakis and Martin, 2004). An excellent recent review about MEN can be found in Short and Isafiade (2021). Especial interest has the simultaneous heat and water integration because it takes simultaneously into account two of the major challenges that process engineering must face: the efficient use of energy and water -See for example Ibrić et al. (2021), Kamat and Bandyopadhyay (2021a), Kamat and Bandyopadhyay (2021b). Even though in this paper we only focus on heat integration, the similarity between the algorithms in MEN and HEN, indicates that likely the proposed methodology can be extended to the design of mass exchanger networks as well. In any case, if we consider two interconnected networks: HEN and water reuse and regeneration, the heat exchanger network can be designed using the current approach.
With considerable evolution in computer processing capacities in the 80's and 90's, more automated HEN synthesis schemes were developed. In special, mathematical programming models achieved considerable success in obtaining good solutions to the problem. These models were presented first in a sequential fashion. Optimization models were developed to mimic crucial algorithmic concepts from the pinch method and some are worth highlighting here. In an analogy to the transshipment model, the optimization linear programming (LP) mathematical formulation of Papoulias and Grossmann (1983) was able to obtain energy targets as in the pinch methodology. In the same work, the authors proposed a second step where another optimization model enabled obtaining the minimum number of heat exchange units. From that solution, a HEN structure could be straightforwardly derived. The second step optimization model has a mixed-integer linear programming (MILP) formulation, solvable via branch and bound methods. Steps comprising models for minimum utility use and a minimum number of matches were also employed by Floudas et al. (1986) with some additional steps. The first was the derivation of a superstructure from the minimum units solution where each stream could have its possible interconnections. In the second step, non-linear programming (NLP) optimization model was applied so that flow rates and heat duties were optimized with the objective of area minimization, therefore minimizing total annual costs. Other important sequential approaches include the vertical heat exchange MILP model of Gundersen and Grossmann (1990) where approach temperatures were treated independently for the whole heat recovery project and each heat exchanger; the improved version of those techniques where an NLP model was included to optimize areas (Gundersen et al., 1997); the modification of the minimum units number sub-problem by Anantharaman et al. (2010) using physical insights/heuristics; Chen et al. (2015) weighted objective function for utility/capital costs and reformulations for MILP transshipment models; the two-step MILP/MINLP method of Nemet et al. (2018) which was able to obtain promising results, including a network for a 173streams problem.
A disadvantage of sequential approaches is that decisions taken in a given step are carried throughout the following ones. Especially, in HEN synthesis, we can assess that a match selection step may be a source of uncertainty. Even though when methodologically performed, globally optimal match choice is not possible given the large number of possible combinations (exhaustive evaluation is infeasible in large problems). To overcome such uncertainties, simultaneous (one-step) mathematical programming models were developed. A very prominent one is Floudas and Ciric (1989) hyperstructure. It contained a wide range of possibilities for HEN configurations including stream splitting, by-passes, etc. The model has a mixed integer non-linear programming (MINLP) problem formulation and was solved with a decomposition approach. From a theoretical point of view a mathematical programming-based approach that simultaneously takes into account the trade-off between area (investment cost) and energy (operating costs), and at the same time, offers a rich space of alternatives for the connectivity between heat exchangers (heat exchange in series, parallel, or the combination of parallel and series, bypasses, etc.) should provide the best possible solution. However, the complexity of the resulting model, usually a non-convex MINLP model, constraints the applicability of those models to short to medium-scale problems. Likely the most successful simultaneous model was presented by Yee and Grossmann (1990), the model known as SYNHEAT is a stage-based approach in which some feasible alternatives are sacrificed to maintain the model simple and robust. A good number of modifications/improvements have been presented to the original version to deal for example with non-isothermal mixing (Björk and Westerlund, 2002;Huang et al., 2012) multiple utilities (Ponce-Ortega et al., 2009), etc. Recently Chang et al. (2020a, Chang et al. (2020b) proposed a smart enumeration algorithm that authors conjecture is the global optimal solution and can tackle even large-scale problems.
Moreover, within simultaneous approaches for HEN synthesis, a noteworthy research trend focuses on investigating alternative solution methods from those that are deterministic. Namely, these approaches are based on meta-heuristics. These methods are derivative-free methods and, to some degree, involve randomized search techniques to improve solutions. These methods were able to provide promising solutions to simultaneous models, in special to SYNHEAT (Yee and Grossmann, 1990) and variants. The following are noteworthy contributions under this scope. Lewin et al. (1998) employed genetic algorithm (GA) techniques to solve the SYNHEAT model obtaining good near-optimal results. That work proposed HEN structures with GA while continuous variables were optimized with a deterministic approach (sequential quadratic programming), giving rise to a hybrid method. GA was also employed by Ravagnani et al. (2005) in a hybrid approach that also included pinch-based techniques. Peng and Cui (2015) applied a GA approach for the HEN structure and a continuous SA (CSA) adaptation to the continuous variables (heat duties) in a superstructure based on SYNHEAT but with no stream splits. Pavão et al. (2017a) employed SA to the structure level and a hybrid CSA/particle swarm optimization (PSO) approach to the continuous level, which was called rocket fireworks optimization (RFO). The superstructure used was the stage-wise superstructure from SYNHEAT but without the isothermal mixing assumption. More recent efforts include the improved GA by Rathjens and Fieg (2020) and Random Walk algorithms with non-structural models (Xiao et al., 2021).
Simultaneous approaches are in general considered superior to sequential ones because in the sequential approach the decisions in a stage are maintained in the rest of the stages and a costly iterative process is needed if it is desired to modify some of the decisions of previous stages. However, some numerical tests showed that in some cases, the sequential approach can get equal or even better results than the simultaneous one (i.e., SYNHEAT) even in small scale problems at least for two reasons: In large scale problems simultaneous models get trapped in poor local solutions and global optimization deterministic solvers can, in the best case, provide some local solutions but most of the times, they are not able to achieve proof of global optimality. This problem is mitigated in the sequential approach-although it is not avoided-because models are easier to solve. Besides, in the sequential approach, the superstructures for the final design of the network include alternatives that are not being considered in the most common simultaneous approaches.
In this work, we present an alternative two-stage sequential model that tries to combine some of the advantages of the simultaneous approach and the richer space of alternatives of the last stage of the sequential one. The first stage is an adaptation of the TransHEN model (Nemet et al., 2018).
The TransHEN model relies on the concept of temperature intervals and requires, a priori fixing a minimum heat recovery approach temperature (HRAT). In the "classical" approach for calculating the minimum utilities using a Tableau (Linnhoff et al., 1982), the intervals of temperature are formed only by the inlet temperatures of the process streams. The reason is that the pinch point is always induced by the inlet temperature of a process stream. However, Nemet et al. (2018) pointed out that other alternatives like including both the inlet and outlet temperatures or even using some intermediate equally spaced intervals could help to improve the estimation of area. In this work, we consider both the inlet and outlet temperatures, that although could eventually create some extra temperature intervals facilitate the implementation, maintains a good accuracy in area estimation, and only marginally complicates the model.

METHODOLOGY
The TransHEN model combines in a single model the calculation of utilities and the cost (or area) estimation. The model takes the form of a transportation problem in which a given hot stream (either a process or a utility) can exchange (transport) heat between the temperature intervals in which that stream is present and the temperature intervals where a cold stream is present. The temperature of the cold streams in the pair of intervals that exchange heat must be lower than the temperatures of the hot streams. The great advantage of this transportation model is that it is possible to calculate a priori the logarithmic mean temperature difference when heat is exchanged between interval "k" and interval "kk". In that way, the area needed to exchange heat between a hot stream in interval k and a cold stream in interval kk is maintained linear in the model.
A comprehensive description of this model is as follows: Let us assume that the following data is known: A set of hot process streams with their inlet and outlet temperatures ( Cost coefficients for the investment cost estimation of any possible heat exchanger. It is assumed that the cost of a heat exchanger can be correlated with an equation of the type: Given a minimum approach temperature (ΔT min or HRAT) it must be generated a shifting scale of temperatures in such a way that a hot stream and a cold stream with the same temperature T* in the shifted scale are separated in the actual scale by a difference ofΔT min . To that end, we subtract ΔT min /2 to the inlet and outlet temperatures of the hot process or utility streams and add the same quantity to the cold process or utility streams. All the temperatures (inlet and outlet, hot and cold, and process and utility streams) in the shifted scale are sorted. Two consecutive temperatures define a temperature interval. Figure 1 shows a sketch of the main elements of the model.
Formally, the following index sets must be defined: Note that the set INT is an ordered set. The lower the position of k the higher the temperature.
It is also necessary to define some sets that explicitly take into account to/from which temperature intervals a given stream provide/demand heat: HK i,k {The hot process stream i gives heat to interval k } CK j,k {The cold process stream j demands heat from interval k } HUK m,k {The hot utility m can release heat to interval k} CUK n,k {The cold utility n can remove heat from interval k} Once the temperature intervals have been established the heat provided/demanded by the process streams in each temperature interval can be calculated: T k , T k refer to the higher and lower temperature in the interval k respectively. The logarithmic mean temperature difference for all the feasible heat exchanges between a hot stream in interval k with a cold stream in interval kk is as follows; Finally, it is possible a priori calculate the global heat transfer coefficients for each possible match between any hot or cold stream.
We have now all the data required by the model. It is also necessary to define the following positive variables: Q HU m Total heat flow provided by hot utility m. Q CU n Total heat flow removed by cold utility n. Q i,j,k,kk Heat exchanged between hot stream i in interval k and cold stream j in interval kk. Q I−CU i,n,k,kk Heat exchanged between hot stream i in interval k and cold utility n in interval kk.
Q HU−J m,j,k,kk Heat exchanged between hot utility m in interval k and cold stream j in interval kk.
Area i,j Area of the heat exchanger between process streams i and j.
Area CU i,n Area of the heat exchanger between process hot stream i and cold utility n.
Area HU m,j Area of the heat exchanger between hot utility m and process cold stream j.
CostHE i,j Cost of the heat exchanger between process stream i and j.
CostHE CU i,n Cost of the heat exchanger between hot process stream i and cold utility n.
CostHE HU m,j Cost of the heat exchanger between hot utility m and process cold stream j.
And the following binary variables: y i,j It takes the value of 1 if streams i and j exchange heat, and 0 otherwise. y CU i,n It takes the value of 1 if stream i exchanges heat with utility n, and 0 otherwise. y HU m,j It takes the value of 1 if utility m exchanges heat with stream j, and 0 otherwise.
The objective is to minimize the total cost.
In the previous equation, it is worth remarking that investment costs must be in the same units as operational costs (i.e., we consider annualized costs).
The rest of the model requires energy balances in points 1, 2, 3, 4 of Figure 1.
The heat released by the hot stream i into interval k can be exchanged with the cold process or utility streams in intervals kk at lower or equal temperatures than interval k -remember that in the shifted scale of temperatures the actual difference for the same value for a hot and a cold stream is actually ΔT min -(point 1 in Figure 1): The heat demanded by cold stream j from interval kk can be exchanged with the hot process or utility streams in intervals k at higher or equal temperatures than the interval kk (point 2 in Figure 1) Energy balance in points (3) and (4) in Figure 1. Heat exchanged by the hot and cold utilities: If a heat exchanger does not exist then the heat flow must be forced to be zero.
Q UP i,j , QCU UP i,n , QHU UP m,j are the maximum heat flow that hot and cold streams can exchange: The area of a heat exchanger can be calculated by adding the areas of the intervals in which two streams exchange heat: Finally, the cost of a heat exchanger is usually a non-linear function in terms of the area. To maintain the linearity of the model the actual cost is substituted by a piecewise linear approximation. To decrease the error in the piecewise optimization we calculate the optimal distribution of intervals to minimize the error between the actual cost and the linear approximation. Figure 2 shows an example of the optimal distribution for the case in which the cost function is approximated by 5 linear intervals.
Assuming that the cost function is approximated by «t» area intervals, to model the piecewise approximation it is necessary to introduce a new index set: And the following new Boolean variables: YD i,j,t Takes the value true if there is a heat exchanger between streams i and j in the area interval t.
YD HU m,j,t Takes the value true if there is a heat exchanger between utility stream m and process stream j in the area interval t.
YD CU i,n,t Takes the value true if there is a heat exchanger between process streams i and utility stream n in the area interval t.
The cost functions can be approximated by the following disjunctions: The following disjunctions can be reformulated in terms of binary variables using the big M or Hull reformulations (Trespalacios and Grossmann, 2014). The hull reformulation, used in this work, is as follows: y HU m,j t∈IT yd HU m,j,t ∀m ∈ HU, j ∈ COLD y CU i,n t∈IT yd CU i,n,t ∀i ∈ HOT, n ∈ CU where yd i,j,t , yd HU m,j,t , yd CU i,n,t are binary variables that take the value 1 if the corresponding Boolean variable takes the value of True.
Area LO m,j,t yd HU m,j,t ≤ Aread HU m,j,t ≤ Area UP m,j,t yd HU m,j,t ∀m ∈ HU, j ∈ COLD, t ∈ IT (31) i,n,t + b CU i,n,t Aread CU i,n,t ∀i ∈ HOT, n ∈ CU, t ∈ IT Finally, those variables in the temperature intervals in which is not possible a heat exchange, either because the stream does not exist in that interval or there is a violation in the driving force must be fixed to zero: The final model is a MILP model formed by Eqs 7-38.
In the case that in the model isothermal streams (either process or utility streams) appear, there are two alternatives: 1) introduce a small fictitious temperature variation (i.e., 1°C). In this case, further modifications are needed, or 2) assume a discontinuity in the set of intervals (or an interval with 0°C difference). In this case, it is necessary to substitute the calculation of heat in this interval by the product of latent heat by mass flowrate, the rest of the model is not further modified.
Even though it is explicitly included in the model, it is worth pointing out that considering multiple utilities at different temperature levels with different costs is straightforward. The model is also very flexible; it allows to forbid matches between streams or to avoid matches with very low (or too large) areas.

The HENDesign Model
The TransHEN model gives us information on the matches between different streams, the heat exchanged, and a good estimation of the area and costs of all heat exchangers. With all those information it is possible to generate a valid heat exchanger network. However, those results were obtained assuming a fixed heat recovery approach temperature. Besides, There are information about which streams should exchange heat but if a given stream exchanges heat with more than one stream, we do not have information on the particular arrangement (series, parallel, or a combination of both) with better performance in a performance index, for example, total cost. Unlike the classical sequential approach, in which in the final stage the heat exchanged by the process and utility streams are fixed, in the HENDesign model, the only parameter that must be fixed is the set of matches predicted by the TransHEN model. But all the heats exchanged between process and utility streams are reoptimized.
The superstructure is inspired by that proposed by Floudas et al. (1986). The superstructure and the basic notation are sketched in Figure 3.
For the HENDesign model, when possible, we use the same notation as in the TransHEN model, However, there are differences in some parameters, and variables. Therefore, to avoid confusion, in the following paragraphs, the complete description of index sets, parameters, and variables of the model are included.
Let us define the following index sets: HS { i | i is a hot stream. It includes both process and utility streams} CS { j | j is a cold stream. It includes both process and utility streams} FH i Heat flow of hot process streams (i ∈ HS). FC j Heat flow of cold process streams (j ∈ CS). U i,j Overall heat transfer coefficient between hot stream i and cold stream j.
Cf i,j Coefficient in cost estimation for exchanger i−j. Cv i,j Coefficient in cost estimation for exchanger i−j. β i,j Exponent in cost estimation for exchanger i−j. Cost HU i Unitary cost of the hot utility i (i ∈ HU). Cost CU j Unitary cost of the cold utility j (j ∈ CU). The continuous variables: FH i Heat capacity flowrate. If i is a utility FH i is a variable, otherwise, its value is fixed (it is a parameter).
FC j Heat capacity flowrate. If i is a utility FC j is a variable, otherwise, its value is fixed (it is a parameter).
Note that the cost of heat exchangers must be in the same units as the cost of utilities (i.e., annualized cost). When a "classical" sequential approach is solved in which first the energy costs are calculated (cost of utilities) and then those values are fixed to determine the net with minimum capital cost, the same network is obtained independently of whether the capital costs are annualized or not. However, in simultaneous approaches, this is something to take into account. In most of the benchmark problems used to test algorithms to synthesize heat exchanger networks, especially those that do not represent an actual process, it is not clear whether the costs are annualized or not. However, usually, the capital and operational costs are simply added without any comment. In this work, we assume that all capital costs are coherent with the operational cost.
The constraints of the model are mass and energy balances in the mixers and splitters of the superstructure, heat exchanger design equations, and cost estimation: Material balance in the initial mixer for hot and cold streams (see Figure 3).
Note that FH i , FC j are fixed for the process streams (they are known data) but are variables for utility streams.
Mass and energy balance in the final mixer for hot and cold streams (see Figure 3) Material and energy balance in the mixer previous to each heat exchanger (see Figure 3) Material and energy balance in the splitter at the exit of each heat exchanger (see Figure 3) Heat exchangers design equations The logarithmic mean temperature difference (LMTD) is approximated by the Chen equation (Chen, 1987) to avoid numerical problems if the temperatures of both ends in a heat exchanger are equal.
Cost estimation: If a heat exchanger does not exist, its heat flow, and consequently the area, must be zero.
Bounds on variables: The previous model takes the form of a Mixed Integer NonLinear Programming (MINLP) problem. However, considering that the set of matches is fixed the model could be just a Nonlinear programming problem (NLP). Two reasons justify using an MINLP model instead of simply an NLP: 1) The optimization of the heat exchanged between the different process and utility streams could eventually remove some of the heat exchangers initially assumed. If this were the case, the constant term in the cost must be forced to be zero. If we want to solve an NLP model either we check that in the final solution all the heat exchangers initially assumed have areas different from zero or a priori fix a minimum area in all heat exchangers. 2) If the optimal solution-set of matches-is very dependent on the assumed HRAT, it is also possible to extend the superstructure with the union of matches obtained at different HRATs. In this case, the MINLP is mandatory because the optimal solution will only contain a subset of those matches. Nevertheless, extending the superstructure produces more complex models, that are considerably more difficult to converge to good solutions.
An alternative to extending the MINLP consists of solving the TransHEN-HENDesign problem for different values of HRAT. Taking into account that the HENDesign problem indirectly "reoptimize" the minimum approach temperature, it is not necessary to check too many HRAT values. There is a tradeoff between solving more simpler problems versus just a single large problem. In this work, we followed the first approach (test different HRAT values) with very good performance even in large-scale problems.

Examples: Comparison With Benchmark Problems
To test the new model we have selected a set of fifteen wellknown problems from open literature that has been widely used to test heat exchanger network algorithms and range from small problems with three process streams to a large scale one with 22 hot streams and 17 cold streams. The models were implemented in GAMS (Bussieck and Meeraus, 2004) (version 34.3) using CPLEX (CPLEX IBM ILOG, 2009) as MILP solver and the global solver BARON (Sahinidis, 1996) as MINLP solver on a PC machine working under windows (i7 2.90 GHz, 32.0 GB). The MILP problems were solved to global optimality (optimality GAP fixed to zero) except for the two largest ones that were stopped with relative optimality gaps of 1.2 and 3.5% after 2 h of CPU time. In the case of MINLP problems BARON is not able to close the optimality gap, however, it is able to find very good solutions in the first minutes of execution and then typically is stacked for hours. We stopped the BARON execution after 1,000 or 3,600 s of CPU time depending on the problem size.
The TransHEN model provides a good estimation of the structure of the matches between hot and cold streams (either process or utilities), the heat exchanged by each couple of streams, the area, and the cost of the network. The major limitation of the TransHEN model is the discrete nature of the temperature intervals: the logarithmic mean temperature difference is calculated between the extreme temperatures of the intervals, and even though a given stream can exchange only a fraction of its heat in an interval the TransHEN model assumes that the inlet and outlet temperatures are those of the interval. This problem can be mitigated by adding new temperatures to create new intervals. Notwithstanding, if the TransHEN model can capture the structure of heat exchanges, the HENDesign model corrects the inaccuracy due to the incorrect estimation of logarithmic mean temperature differences.
In the TransHEN model, it is possible to get solutions in which two streams exchange heat in nonconsecutive temperature intervals. That means that there should be two heat exchangers with the same hot and cold streams. In those cases, the area estimation is correct, but the model assumes a single heat exchanger, and therefore the cost is underestimated (except in those cases in which the cost varies linearly with the area). It is possible to modify the TransHEN model to either, forbid that situation, or take into account this fact, but it requires a more complex model with new binary variables. The good results obtained suggest that those modifications are not necessary.
The determination of the correct number and distribution of temperature intervals is case-dependent. We use the following heuristic: we solve twice the TransHEN model. First with the original set of temperature intervals and then adding new intervals. By default, we simply double the number of intervals, but other strategies like dividing the larger ones or those with larger heat flows are also valid alternatives. If the results are similar we go to the HENDesign step, if not we increment again the number of intervals. Note that the larger the number of intervals, the lower the objective function. Of course, increasing the number of intervals increases the size of the problem and the time to solve it, but numerical tests have shown that usually, it is not necessary to create new intervals.
Interestingly, in large problems, there is, in general, an important overlapping of streams which results in a large number of  (2020) a The digits after the "H" and "C" indicate the number of hot and cold process streams respectively.
In bold appears the best known solution after this work. temperature intervals. The consequence is that adding new temperature intervals has a small effect on the final result. The HENDesign model assumes that two streams can only exchange heat between them once. Therefore, although it is possible to modify the TransHEN model to consider the possibility of more than a single exchange between two streams, that change will be not reflected in the HENDesign model. It is possible to modify the superstructure in the HENDesign model to allow that situation, however, the good results obtained in the test problems suggest that the modification is likely to have a small impact on the final solution.
Unlike the usual approach in the sequential algorithm, in HENDesign we re-optimize the network maintaining the structure of heat exchanges predicted in TransHEN. While the model is highly non-linear and nonconvex, the good initial values provided by the transportation model and the reduced set of heat exchangers allow global deterministic solvers to obtain good solutions in some minutes of CPU time. Even though global optimality cannot be guaranteed, we can guarantee a very good solution as prove the excellent results obtained in the test problems.
In the Supplementary Material, there is a comprehensive description of the fifteen examples, a comparison with available solutions, and the best solution obtained following the TransHEN-HENDesign approach. Table 1 shows a summary of the main results, and

CONCLUSION
In this work, we have presented a new algorithm that tries to get the advantages of the sequential and simultaneous approaches for the design of heat exchanger networks. The first part is an adaptation of the TransHEN model presented by Nemet et al. (2018). The model maintains the concept of temperature intervals and postulates all the possible heat transfers between the hot and cold stream in each of the temperature intervals. The explicit consideration of the heat transportation of individual streams between temperature intervals allows the a priori determination of the logarithmic mean temperature difference between all possible heat exchanges, which permits maintaining the area estimation linear in the model. To maintain the linearity the cost of a heat exchanger is approximated by a piecewise linear approximation.
In the second step, a superstructure inspired by the works by Floudas et al. (1986) that contains all the possible alternatives in which two streams can exchange heat allow the final design of the heat exchanger network. Unlike the sequential approach, in this model, all heat flows, temperatures, areas, etc. are reoptimized maintaining the set of matches predicted in the first stage.
Although the model is nonlinear and nonconvex it is relatively easy to get good results, because the model starts with the values predicted by the TransHEN model.
The TransHEN model maintains the concept of temperature intervals and therefore requires de a priori definition of a minimum approach temperature (HRAT). While this is generally an important drawback-especially if the assumed value is far away from the best one-the reoptimization in the HENDesign model minimizes that problem. Of course, it is still needed to approximately know an estimation of the HRAT, but in general, a precise value is not required and any of the iterative pinch-based approaches to estimate the area (or cost) of the network changing the HRAT provide good initial estimations.
Fifteen examples commonly used in the literature to evaluate algorithms for heat exchangers network design were tested. In most cases, we got better or equal solutions than the best-known solutions published in the open literature. For the rest, the optimal solution obtained is only marginally worse than the best-reported solution (difference lower than 0.5%) except for the larger case in which the difference is around 2.4%.
A comparison between the best obtained solution and the best-known solution until that moment does not show any kind of correlation with problem size or number of hot or cold process streams. The differences besides are really small, the largest one is around 4%. This is at indication that, at least inside the problem size of the tested benchmark problems there are no an apparent degradation in performance of the proposed model. However, although there is a correlation between the number of process streams and the difficulty in solving the problem (measured in CPU time) this is not the only factor. For example, the degree of overlapping between streams is directly related with the number of intervals generated and consequently the problem size (i.e., discrete variables) in the TransHEN model. For example, the CPU time for solving problem 7 (H5C5) was around 33 s while in problems 1 to 10 (they are sorted by the total number of process stream) is below 1 s. Even more, the solution of the TransHEN model has a large influence on the HENDesign model. For example, if a given hot/cold stream exchange heat with a large number of other cold/hot streams the number of non-convex terms (and the size of the model) could be larger than another problem with more process streams in which the solution of the TransHEN model predicts that all the streams exchange heat with a reduced number of other streams.
Of course, we cannot know if the result is close or not the global optimal solution (except for comparison with other published models). We would say that, with the actual software and our available computer capacity, we are not too far away of the maximum problem size we can solve without any kind of reformulation, but it is difficult to know. But again we have not information about how our solution deteriorates with the problem size and or problem features.
The algorithm has proved to be robust, and efficient and the test has shown good performance in medium to large heat exchanger network design problems.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.