Optimal Design of Heat-Integrated Water Allocation Networks

Abstract

Industrial operations consume energy and water in large quantities without accounting for potential economic and environmental burdens on future generations. Consumption of energy (mainly in the form of high pressure steam) and water (in the form of process water and cooling water) are essential to all processes and are strongly correlated, which requires development of systematic methodologies to address their interconnectivity. To this end, the subject of heat-integrated water allocation network design has received considerable attention within the research community in recent decades while further growth is expected due to imposed national and global regulations within the context of sustainable development. The overall mathematical model of these networks has a mixed-integer nonlinear programming formulation. As discussed in this work, proposed models in the literature have two main difficulties dealing with heat–water specificities, which result in complex formulations. These difficulties are addressed in this work through proposition of a novel nonlinear hyperstructure and a sequential solution strategy. The solution strategy is to solve three sub-problems sequentially and iteratively generate a set of potential solutions through the implementation of integer cut constraints. The novel mathematical approach also lends itself to an additional innovation for proposing multiple solutions balancing various performance indicators. This is exemplified with both a literature test case and an industrial-scale problem. The proposed solutions address a variety of performance indicators which guides decision-makers toward selecting the most appropriate configuration(s) among a large number of potential possibilities. Results exhibit that despite having a sequential solution strategy, better performance can be reached compared to previous approaches with the additional benefit of providing many potential solutions for further consideration by decision-makers to select the best case-specific solution.

1. Introduction

Heat-Integrated Water Allocation Network (HIWAN) design has been extensively studied in the literature in recent decades [1,2,3,4,5]. Multiple water- and energy-related specifications of such systems have been addressed by the literature [5], including multi-contaminant problems [6] and interplant operations [7,8,9,10]. Several researchers have also proposed methodologies to incorporate other aspects of an industrial process in the pursuit of exploiting potential synergies between them. Ahmetović et al. [4] and Kermani et al. [5] have carried out systematic review of the available methodologies and have identified several gaps that have received no/little attention within the research community. These include:

• providing systematic holistic methodologies addressing non-water processes [10,11,12,13], thermal utilities [10,14,15], and other energy conversion and heat recovery technologies [15],

• considering multi-period operations and uncertainty analysis of HIWANs addressing variations of operating conditions,

• proposing Multi-Criteria Decision Making (MCDM) tools addressing conflicting set of objectives and providing a set of optimal solutions to guide decision-makers in selecting the most appropriate configuration(s) [15], and

• proposing efficient solution strategies to address the increasing difficulty of solving large-scale cases [16,17].

Kermani et al. [5] further provided a comprehensive analysis of the proposed HIWAN methodologies with particular focus on mathematical approaches. Mathematical methodologies are capable of addressing a variety of water specifications by constructing thorough superstructures that encompass all possible interconnections within the system. The overall problem comprises binary variables and nonlinear (including bilinear) constraints and hence is a non-convex Mixed-Integer Nonlinear Programming (MINLP) problem. To address the aforementioned gaps in a systematic manner using mathematical programming techniques, one of the main difficulties to consider is the interaction between the water network and the Heat Exchanger Network (HEN), that has to be well-understood and investigated. As discussed by Kermani et al. [5], this difficulty can be classified into two distinct problems of:

• selecting water streams that participate in heat exchange, and
• selecting the thermal state of these water thermal streams (i.e., whether they are hot or cold).

Depending on the origin or the destination of a water flow, a water stream can be categorized as freshwater, wastewater, inlet (to a process), outlet (of a process), recycling (among processes), or reuse (within a process). Any of these flows can be subject to thermal duties; however, considering all possibilities within the mathematical formulation introduces computational complexities. Therefore, a subset of these streams can be integrated with HEN. Kermani et al. [5] has shown that all proposed methodologies consider freshwater and wastewater interactions with HEN design while other types may or may not be addressed. Having set the water thermal streams, the second issue is to choose whether the stream is hot or cold as the HEN formulation is generally formulated by knowing these sets in advance. Within the literature, freshwater and wastewater streams are typically modeled as a succession of heat exchangers (modeled as cold streams) and splitters, and heat exchangers (modeled as hot streams) and mixers, respectively. For other types of streams, this procedure is carried out by defining two thermal streams and assigning binary variables and imposing constraints to ensure the activation of only one state. These interactions increase the complexity of the mathematical formulations which necessitates development of robust and efficient solution strategies. Kermani et al. [5] provided a detailed classification of the proposed solution strategies. The main classes can be categorized as follows:

• Simultaneous: These approaches consider all variables and constraints in a single formulation, which are mainly formulated as MINLP models. Ahmetović and Kravanja [18] proposed such an approach by placing the potential heat exchangers on freshwater, wastewater, inlet (after mixer), and outlet (before splitter) streams. They solved the water network model to provide a good initialization point for the subsequent MINLP model. Ahmetović and Kravanja [19] further extended this approach by adding potential heat exchangers on recycling and reuse water streams and essentially extending the heat exchange possibilities to all water streams. Nonetheless, they have only modelled freshwater thermal streams as cold thermal streams and wastewater as hot thermal streams.

• Decomposition: One technique to reduce the complexity of solving large MINLP models is via decomposing the model into two steps. Papalexandri and Pistikopoulos [20] and Leewongtanawit and Kim [21] proposed solution strategies by decomposing the overall MINLP model into a Mixed-Integer Linear Programming (MILP) and a Nonlinear Programming (NLP) model and solved the problem iteratively until no improvement was reached. All binary decision variables, i.e., the network configuration, were allocated to the MILP model. Despite the iterative nature of the solution strategies, only one solution was presented as the optimum one.

• Sequential: Unlike decomposition approaches, sequential approaches lack iterations and hence are considered as uni-directional solution strategies. The problem can be formulated in two or more steps, where each step potentially solves one aspect of the problem. Liao et al. [22] proposed a two-step sequential solution strategy, by optimizing the potential thermal matches and the freshwater consumption in the first step and solving the MINLP model of HEN design in the second step. Kermani et al. [14] used a more decomposed approach by proposing a three-step solution strategy, i.e., by minimizing the utility consumption in the first step and minimizing the number of thermal matches in the second step. The third step of their solution strategy is to apply the Pinch Design Method for HEN synthesis. Hong et al. [23] also proposed a three-step solution strategy by formulating an NLP model as the first step and then solving two subsequent MINLP models. The first step minimizes the freshwater consumption, however, lacks consideration of heat integration. Overall, none of the proposed sequential solution strategies addressed the possibility of having multiple solutions to the problem.

• Transformation: Another technique to reduce the complexity of solving an MINLP model is via transforming and reformulating the problem. Yan et al. [24] proposed such an approach by essentially relaxing the integrality of binary variables and hence providing an NLP formulation. Despite this transformation, the proposed superstructure still possesses the same limitation of other ones, namely, the choice of water thermal streams and their state.

Furthermore, Ibrić et al. [16] and Wang et al. [17] have proposed several heuristic to simplify the overall superstructure via eliminating impractical and/or infeasible interconnections and have shown that these techniques can lessen the burdens of solving large problems. Jeżowski [1] have also highlighted the importance of sequential, decomposition, and combination of meta-heuristic techniques to facilitate finding a global or “good” local optimum.

In the design of HIWANs, there exist several sets of conflicting objective functions (e.g., operating cost vs total investment cost). Due to presence of these types of trade-offs, no single solution can be regarded as the best alternative with respect to all objective values. This necessitates developing multi-criteria optimization methodologies to provide a set of “good” solutions to the problem. Chen and Wang [2] and Kermani et al. [5] have highlighted the importance and benefits of multi-objective optimization techniques in HIWAN synthesis problems. Very few research works have considered these techniques in their methodologies [15,25,26]. Boix et al. [26] selected four objective functions defined as minimizing the freshwater consumption, energy consumption, number of heat exchangers, and number of water interconnections. In the first stage, they solved the problem of water allocation and thermal streams selection by combining the first two objectives usingϵ -constraint while fixing the value of the latter two objectives. Although they proposed a Pareto frontier of a set of solutions, they applied the classical TOPSIS method [27,28] to choose the best alternative among all in order to proceed to the second stage. In the second stage, having fixed the water allocation network and extracting the thermal streams, they formulated an MINLP model for minimizing the total HEN cost. Their approach, however, neglects the possibility of mixing and splitting among various water streams.

The innovative aspects included in this work can therefore be summarized as the novel decomposition and sequential solution strategy introduced, a systematic generation and presentation of multiple solutions for MCDM, and the application of these aspects on an industrial case study in addition to validation using a literature test case. The proposed method thus addresses the aforementioned difficulties noted as points 1, 3, and 4, above. The overall HIWAN problem is decomposed into three sub-models. The first two models were introduced in previous work, namely: targeting (problem P1) [14] and Heat Load Distribution (HLD) (problem P2) models [29,30]. Problem P1 is formulated as MILP and its solution provides minimum targets in thermal utilities and freshwater consumption together with potential water thermal streams and their states. Adding integer variables allows all water flows to participate in heat exchange duties. Having solved problem P1, problem P2 provides a feasible set of heat exchange matches among these thermal streams. This work adapts the NLP model of HEN design [31] by including a water network model to solve the total HIWAN synthesis problem as the third model (problem P3). This model is formulated using NLP with the objective of minimizing the annualized capital investment subject to utility targets (results of problem P1, set as upper bound) and heat exchange matches (results of problem P2). To address the trade-offs between resource consumption and investment cost and to further consider non-objective performance indicators, a novel solution strategy is proposed by which the three sub-problems are solved sequentially and iteratively to generate a set of potential solutions through the implementation of integer cut constraints. Furthermore, due to the presence of many other non-objective performance indicators in the system and due to their intertwined interactions, it would be ill-advised to present only a single solution while ignoring other “promising” ones. To this end, MCDM tools which allow decision-makers and practitioners to explore the set of solutions and select the best with respect to their own criteria are promising alternatives to the classical MCDM techniques. They further allow practitioners to understand why they are selecting a specific solution and how this solution’s performance is compared to the rest. Among these tools, parallel coordinates visualization tool [32] is a powerful and practical way for evaluation and decision-making in high-dimensional spaces with applications in urban planning [33,34] and process integration [15]. This tool has been selected to visualize and further select promising alternatives similar to the approach presented by Kermani et al. [15]. The governing equations of problem P3 and the proposed solution strategy are described in Section 3 and Section 4, respectively. The methodology is validated by well-established literature case studies in Section 5.

2. Problem Statement

Two sets of water unit operations (_PinWAN,_PoutWAN) are given. Each unit j in_PinWANrequires water at temperature_Tj, with a maximum allowed inlet contamination level of_cjk,max. Each unit i in_PoutWANsupplies water at temperature_Ti, with maximum allowed outlet contamination level of_cik,max. In addition, for mass-transfer water unit operations, a mass load of_Lumust be removed from unit u. Furthermore, sets of hot and cold non-water process streams (^PH,^PC, respectively) are also considered. Thermal hot and cold utilities (^UH,^UC, respectively) are also available within the system to close the energy balance. In addition, multiple freshwater sources and wastewater sinks (_UoutWAN,_UinWAN) are given. The objective is to design a HIWAN that exhibits minimum Total Annualized Cost (TAC).

3. Mathematical Formulation

As introduced in Section 1 the overall MINLP formulation of HIWAN is decomposed in three sub-problems denoted P1, P2 and P3. In problem P1, a discretization technique is applied in which any water stream can be heated or cooled from its initial temperature to reach any other temperature in the water network [14]. By assigning binary variables to the existence of each thermal stream, the solution to problem P1 yields the minimum operating cost of the network by selecting potential thermal streams within the water network. This allows cold and hot streams to be selected for both source and sink units. Problem P2 minimizes the number of thermal matches among this set of streams. Appendix A and Appendix B provide the governing equations of problems P1 and P2, respectively. The HIWAN hyperstructure proposed for the third step (Figure 1), and the focus of this paper, is based on the NLP hyperstructure of Floudas and Ciric [31] (Figure 1a) and has been adapted to the HIWAN problem by applying the following modifications:

• Final mixers are removed from superstructures of streams associated with water source units. Hence, the modified superstructure allows for stream splitting from the outlet of any water source unit heat exchanger. (Figure 1b)

• Initial splitters are removed from superstructures of streams associated with water sink units. Hence, the modified superstructure allows for stream mixing at the heat exchanger inlet of sink water units. (Figure 1c)

• Bypass streams are added between different superstructures associated to one source or one sink (shown as dashed arrows in Figure 1b,c). This mimics_fl,^l′B,k in the HEN hyperstructure (Figure 1a), which characterizes the flow from the splitter following the match between k and^l′to the mixer preceding the match between k and l within the superstructure of stream k.

• As opposed to the HEN hyperstructure for which heat loads of each match (_Q˙lk) are fixed at their optimal values from the HLD model, the proposed HIWAN hyperstructure relaxes this constraint; hence,_Q˙lkcan become zero indicating that the match is no longer necessary. For non-water thermal streams (cold or hot), the sum over hot or cold streams, respectively, is fixed to the stream heat load. In essence,_Q˙lkfor non-water thermal streams are allowed to deviate from their optimal values determined by the HLD model.

• An important assumption in the proposed hyperstructure is the value of Heat Recovery Approach Temperature (HRAT) for which problem P1 and problem P2 are solved. This value affects the targeting values of thermal utilities as well as the set of potential thermal streams and their matches; therefore, the model is solved iteratively for several HRATs. The NLP hyperstructure is then solved for Heat Exchanger Minimum Approach Temperature (ΔT_min) values equal to HRAT.
• An initial solution to the HEN hyperstructure model is used as an initialization point for the HIWAN problem.

These modifications produce the following novelties:

• Freshwater and wastewater streams can have hot or cold thermal duties. This addition treats situations in which freshwater cooling or wastewater heating is required.
• In most superstructures proposed in the literature, only one thermal stream is assigned to the outlet of a source water unit operation. As a result, no water stream can be split at the outlet of a heat exchanger for mixing with other streams. This issue is handled in the proposed HIWAN hyperstructure by removing the final mixers of source units. This issue is also treated for sink units by removing their initial splitters.

In order to present the constraints in a coherent way, similar set definitions of Floudas and Ciric [31] have been adapted to represent stream superstructures and matches:

(1)

SetMAis the set of thermal matches(k,l):

MA=(k,l)∣k∈HSandl∈CSand_yk,lhld=1

(2)

SetHCSis the set of all thermal streams:

HCS=HS∪CS

where HS and CS are the sets of hot and cold streams, respectively.

(3)

Set_SMkis the set of streams,^k′, matching with streamk∈HCS [31], i.e., set of all streams,^k′, so that either of the two matches(k,^k′)or(k,^k′)is in the set of thermal matches,MA:

_SMk={^k′∣(k∈HS∧^k′∈CS∧(k,^k′)∈MA)∨(k∈CS∧^k′∈HS∧(^k′,k)∈MA)}

(4)

Set_SMi,kis the set of streams matching with streamk∈_Siof water uniti∈WAN, where,_Siis the set of all water thermal streams related to water unit i:

_SMi,k=^k′∣i∈WAN∧k∈_Su∧^k′∈_SMk

For mass flow variables, superscripts I, O, and B denote heat exchanger inlet, outlet, and bypass streams, respectively. 3.1. Objective Function

The objective function is the total network investment cost,^Cinv, which is the sum of investment costs over all heat exchangers:

P3:minCinv=∑(i,j)∈MAα1+α2 Qi,jUi,j LMTDi,jβ=∑(i,j)∈MAα1+α2 Ai,jβ

where,_Qi,jis the heat load of match(i,j)which, as stated before, is not fixed at its optimal value from problem P2 but rather defined as a variable. This makes the objective function non-convex. Parameters_α1,_α2, andβare fixed, area, and exponential cost parameters for heat exchangers. The_LMTDi,jis the logarithmic mean temperature difference of match(i,j)defined as:

_LMTDi,j=ΔT_1i,j−ΔT_2i,jlnΔT_1i,jΔT_2i,j

where,ΔT_1i,jandΔT_2i,jare the temperature differences at the two ends of heat exchanger of match(i,j) . Equation (6) is simplified using Chen’s approximation [35]:

_LMTDi,j=^{ΔT_1i,j·ΔT_2i,j·ΔT_1i,j+ΔT_2i,j21/3}

3.2. Constraints

Figure 1 illustrates the basics of the proposed hyperstructure for a water source, a water sink, and a non-water thermal stream which are governed by the equations presented in this sub-section.

3.2.1. Water Network Mass and Energy Balances

(1)

Mass balance at the initial splitter of source uniti∈_WANout:

_m˙i=∑k∈_Si∑l∈_{SMi,km˙i,lk,I}+∑j∈_WANinm˙i,j+∑k∈_Sj∑l∈_{SMj,km˙i,j,lk,I}

where,_m˙iis the flow rate of source unit i,_m˙i,lk,Iis the flow rate from source unit i towards the inlet of heat match(k,l)in the superstructure of stream k,_m˙i,jis the flow rate between source unit i and sink unit j, and_m˙i,j,lk,Iis the flow rate from source unit i towards the inlet of heat match(k,l)in the superstructure of stream k of sink unit j.

(2)

Mass balance at the splitter following a match of source uniti∈_WANout(k∈_Si,l∈_SMi,k):

_m˙i,lex,k=∑^k′∈_Si∑^l′∈_SMi,^k′k≠^k′⋁l≠^l′_{m˙i,^l′,lB,^k′,k}+∑j∈_{WANinm˙i,l,jk,O}+∑^k′∈_Sj∑^l′∈_{SMj,^k′m˙i,l,j,^l′k,O,^k′,I}

where,_m˙i,lex,kis the flow rate through heat exchanger of match(k,l)in the superstructure of stream k of source unit i,_{m˙i,^l′,lB,^k′,k}is the flow rate from the outlet of heat exchanger of match(k,l)in the superstructure of stream k towards the heat exchanger inlet of match(^k′,^l′)in the superstructure of stream^k′where both heat exchangers belong to source unit i,_m˙i,l,jk,Ois the flow rate towards sink unit j, and_{m˙i,l,j,^l′k,O,^k′,I}is the flow rate towards the inlet of heat exchanger of match(^k′,^l′)in the superstructure of sink unit j.

(3)

Mass balance at the mixer preceding a match of source uniti∈_WANout(k∈_Si,l∈_SMi,k):

_m˙i,lex,k=∑^k′∈_Si∑^l′∈_SMi,^k′k≠^k′⋁l≠^l′_{m˙i,l,^l′B,k,^k′}+_m˙i,lk,I

where,_m˙i,lex,kis the flow rate through heat exchanger of match(k,l)in the superstructure of stream k of source unit i,_{m˙i,l,^l′B,k,^k′}is the flow rate from the outlet of heat exchanger of match(^k′,^l′)in the superstructure of stream^k′towards the heat exchanger inlet of match(k,l)in the superstructure of stream k where both heat exchangers belong to source unit i, and_m˙i,lk,Iis the flow rate from source unit i towards the inlet of heat match(k,l)in the superstructure of stream k.

(4)

Energy balance at the mixer preceding a match of source uniti∈_WANout(k∈_Si,l∈_SMi,k):

_m˙i,lex,k·_ti,lex,k,I=∑^k′∈_Si∑^l′∈_SMi,^k′k≠^k′⋁l≠^l′_{m˙i,l,^l′B,k,^k′}·_{ti,^l′ex,^k′,O}+_m˙i,lk,I·_Ti

where_ti,lex,k,I(_{ti,^l′ex,^k′,O}) is the temperature at the inlet (outlet) of heat exchanger of match(k,l)((^k′,^l′)) in the superstructure of stream k(^k′), and_Tiis the temperature of source unit i.

(5)

Mass balance at the final mixer of sink unitj∈_WANin:

_m˙j=∑k∈_Sj∑l∈_{SMj,km˙j,lk,O}+∑i∈_WANoutm˙i,j+∑k∈_Si∑l∈_{SMi,km˙i,l,jk,O}

where,_m˙jis the flow rate of sink unit j,_m˙j,lk,Ois the flow rate from the outlet of heat match(k,l)in the superstructure of stream k towards sink unit j,_m˙i,jis the flow rate between source unit i and sink unit j, and_m˙i,l,jk,Ois the flow rate from the outlet of heat match(k,l)in the superstructure of stream k of source unit i towards sink unit j.

(6)

Mass balance at the splitter following a match of sink unitj∈_WANout(k∈_Sj,l∈_SMj,k):

_m˙j,lex,k=∑^k′∈_Sj∑^l′∈_SMj,^k′k≠^k′⋁l≠^l′_{m˙j,^l′,lB,^k′,k}+_m˙j,lk,O

where,_m˙j,lex,kis the flow rate through heat exchanger of match(k,l)in the superstructure of stream k of sink unit j and_{m˙j,^l′,lB,^k′,k}is the flow rate from the outlet of heat exchanger of match(k,l)in the superstructure of stream k towards the heat exchanger inlet of match(^k′,^l′)in the superstructure of stream^k′where both heat exchangers belong to sink unit j.

(7)

Mass balance at the mixer preceding a match of sink unitj∈_WANout(k∈_Sj,l∈_SMj,k):

_m˙j,lex,k=∑^k′∈_Sj∑^l′∈_SMj,^k′k≠^k′⋁l≠^l′_{m˙j,l,^l′B,k,^k′}+∑i∈_{WANoutm˙i,j,lk,I}+∑^k′∈_Si∑^l′∈_{SMi,^k′m˙i,^l′,j,l^k′,O,k,I}

where,_{m˙j,l,^l′B,k,^k′}is the flow rate from the outlet of heat exchanger of match(^k′,^l′)in the superstructure of stream^k′towards the heat exchanger inlet of match(k,l)in the superstructure of stream k where both heat exchangers belong to sink unit j,_m˙i,j,lk,Iis the flow rate from source unit i towards the inlet of heat match(k,l)in the superstructure of stream k of sink unit j, and_{m˙i,^l′,j,l^k′,O,k,I}is the flow rate from the outlet of heat match(^k′,^l′)in the superstructure of stream^k′of source unit i towards the inlet of heat exchanger of match(k,l)in the superstructure of stream k of sink unit j.

(8)

Energy balance at the mixer preceding a match of sink unitj∈_WANout(k∈_Sj,l∈_SMj,k):

_m˙j,lex,k·_tj,lex,k,I=∑^k′∈_Sj∑^l′∈_SMj,^k′k≠^k′⋁l≠^l′_{m˙j,l,^l′B,k,^k′}·_{tj,^l′ex,^k′,O}+∑i∈_{WANoutm˙i,j,lk,I}·_Ti+∑^k′∈_Si∑^l′∈_{SMi,^k′m˙i,^l′,j,l^k′,O,k,I}·_{ti,^l′ex,^k′,O}

where_tj,lex,k,I(_{tj,^l′ex,^k′,O}) is the temperature at the inlet (outlet) of heat exchanger of match(k,l)(^k′,^l′)in the superstructure of stream k (^k′) and_Tjis the temperature of sink unit j.

(9)

Energy balance in matches of streamk∈_Suof water unitu∈WAN:

_Qk,l=_{m˙i,lex,k cp}|_ti,lex,k,O−_ti,lex,k,I|

3.2.2. Contamination Constraints

(10)

Contamination equality constraint at the final mixer of sink unitj∈_WANinfor contaminantk∈C:

_m˙j·_cjk,in=∑i∈_WANout(_m˙i,j+∑k∈_Si∑l∈_{SMi,km˙i,l,jk,O}+∑^k′∈_Sj∑^l′∈_{SMj,^k′m˙i,j,^l′k′,I}+∑k∈_Si∑l∈_{SMi,km˙i,l,j,^l′k,O,^k′,I})·_cik,out

where_cjk,inis the contamination level of contaminant k at the inlet of sink unit j, and_cik,outis the contamination level of contaminant k at the outlet of source unit i.

(11)

Contamination at the inlet of unitj∈_WANinshould be less than a maximum allowed value:

_cjk,in−_cjk,in,max≤0

(12)

Mass transfer equality constraint for mass transfer unitu∈WAN:

∑j∈_WANu,inm˙j·_cjk,in+_Lu−∑i∈_WANu,outm˙i·_cik,out=0

where_m˙j=_m˙iand_WANu,in(_WANu,out) is the set of sinks (sources) associated to mass transfer unit u.

(13)

Outlet contamination of uniti∈_WANoutshould be less than a maximum allowed value:

_cik,out−_cik,out,max≤0

3.2.3. Mass and Energy Balances in Non-Water Streams Superstructure

The mass and energy balances of non-water thermal streams follow the same constraints of the HEN synthesis problem of Floudas and Ciric [31] (Figure 1a introduces the full set of related variables and parameters):

(14)

Mass balance at the initial splitter of streamk∈HCS:

^Fk−∑l∈_SMkflk,I=0

where,^Fkis the heat capacity flow rate of stream k and_flk,Iis the heat capacity flow rate of stream k towards the inlet of heat exchanger of match(k,l)/(l,k)in the superstructure of stream k.

(15)

Mass balance at the mixer preceding a heat exchanger of streamk∈HCS:

_flex,k−_flk,I−∑^l′∈_SMk,l≠^l′_fl,^l′B,k=0∀l∈_SMk

(16)

Mass balance at the splitter following a heat exchange match of streamk∈HCS:

_flex,k−_flk,O−∑^l′∈_SMk,l≠^l′_f^l′,lB,k=0∀l∈_SMk

(17)

Energy balance at the mixer preceding a heat exchange match of streamk∈HCS:

_flex,k·_tlex,k,I−_flk,I·_Tk−∑^l′∈_SMk,l≠^l′_fl,^l′B,k·_t^l′ex,k,O=0∀l∈_SMk

(18)

Energy balance of heat exchange match(k,l)∈MA(as mentioned previously, heat load is defined as a variable in the HIWAN model):

_Qk,l=_flex,k·(_tlex,k,I−_tlex,k,O)_Qk,l=_fkex,l·(_tkex,l,O−_tkex,l,I)

(19)

Energy balance of hot and cold non-water thermal streams:

_Qk=∑l∈_SMkQk,l∀k∈HCS|k∉⋃u∈WAN_Su

This constraint holds for all thermal streamsk∈HCSthat do not belong to any water unit operation._Qkis the heat load of non-water thermal stream k.

3.2.4. Temperature Difference Constraints

(20)

Temperature difference at the two ends of heat exchanger(k,l)∈MA:

ΔT_1k,l=_tlex,k,I−_tkex,l,OΔT_2k,l=_tlex,k,O−_tkex,l,I

(21)

Minimum approach temperature difference:

ΔT_1k,l≥Δ_Tmin=HRATΔT_2k,l≥Δ_Tmin=HRAT

4. Solution Strategy

As discussed in Section 3, the HIWAN synthesis problem is decomposed into three sub-problems:

Step 1—Targeting Model (problem P1): At this step, a value of HRAT is assumed for water thermal streams. For non-water thermal streams, this value depends on the nature of the stream and problem requirements or can be assigned based on general assumptions. A solution of this model provides the utility targets of the HIWAN synthesis problem and a list of water thermal streams with their corresponding heat loads. Moreover, a water allocation network will be generated satisfying the temperature and contamination constraints (Appendix A) [14].

Step 2—Heat Load Distribution Model (problem P2): The input of this model is the list of thermal streams with their corresponding heat loads from Step 1. A solution of this model provides a set of potential thermal matches with their corresponding heat loads which exhibits the minimum number of heat exchangers. Several solutions may exist with the minimum number of heat exchangers. This is addressed in the proposed solution strategy by applying Integer-Cut Constraints (ICCs). The water allocation network is fixed at its solution from Step 1 (Equation (A2)) for this stage.

Step 3—HIWAN Hyperstructure Model (problem P3-NLP): The model takes the water data and utility targets from Step 1, while the thermal matches and their heat loads are provided from Step 2. Since the water allocation network superstructure is embedded in this model, the water allocation network from Step 1 can be changed but provides an initial, feasible solution. Additionally, heat loads of thermal matches are relaxed at this step; hence, heat exchanger flows can reach zero, indicating that the optimized water network flows and temperatures do not require the existence of this exchanger. The variables of this problem, i.e., flows, temperatures, and heat loads of all heat exchangers, are initialized by solving the HEN hyperstructure of Floudas and Ciric [31] (model P3^init-NLP). It should be highlighted that a solution to P3^init-NLP is a feasible solution to the HIWAN synthesis problem.

Solving the above-mentioned steps in sequence provides, at most, two solutions to the HIWAN synthesis problem: the solutions of P3^init -NLP and P3-NLP. Compared to simultaneous approaches, sequential approaches are easier to solve due to reduced problem size at each step. Nevertheless, similar drawbacks observed for sequential solution strategies that are applied in HEN synthesis problems [31] are also valid in this case:

• The trade-off between operating cost and investment cost (i.e., number of matches and area) cannot be fully captured in sequential solution strategies.
• The nature of the HIWAN synthesis problem is non-convex and hence global optimality cannot be guaranteed.
• Problem P2-MILP provides only one feasible solution providing the minimum number of heat exchangers. Therefore, sequentially minimizing the number of heat exchangers followed by HEN synthesis does not guarantee a HIWAN exhibiting the lowest total cost.

To address these drawbacks, an iterative sequential solution strategy is proposed by applying an integer cut constraint on problems P1 and P2. In addition, the problem is solved for different values of HRAT to address the pinch point decomposition of temperature intervals and to get different utility targets with the goal of achieving lower investment costs. Applying integer cut constraints at Steps 1 and 2 provides a double benefit:

► By assuming the overall synthesis problem encompasses all water thermal streams and states (hot and cold), problems P1 and P2 limit the search space to a specific region, i.e., a specific set of potential thermal streams with their corresponding states and a specific set of thermal matches. This allows for a smaller problem size at Step 3 and hence provides guidance towards a good solution. Applying an integer cut constraint generates several “reduced-size” problems and increases the opportunities for identifying good solutions.

►

By applying, for example,_NiccP1and_NiccP2integer cuts at Steps 1 and 2, respectively,(_NiccP1+1)(_NiccP2+1)solutions to the HIWAN synthesis problem can be generated. These solutions can be ranked based on different key performance indicators and thus guide the decision-making process toward a “good” solution.

This approach can be regarded as iteratively optimizing the overall MINLP model within different regions of the global search space (as imposed by the formulations of problems P1 and P2 at each iteration). Figure 2 illustrates the proposed solution strategy. It should be noted that, similar to other methodologies, the choice of solver and its options affects the path toward a solution.

5. Illustrative Example

The four-water-unit-operation, single-contaminant case study of Savulescu and Smith [36] was selected to demonstrate the generation of the proposed HIWAN hyperstructure. First, problem P1 is solved providing minimum freshwater, hot, and cold utility consumption of 90 kg/s, 3780 kW, and 0 kW, respectively. Problem P1 also provides a feasible water allocation scheme. Given the list of thermal streams, problem P2 is solved which results in seven matches (Figure 3a). There are several solutions to problem P2 with an objective value of seven. For the sake of illustration, the fourth solution (by applying three integer cuts) is shown in Figure 3a and will be discussed hereafter. The CPLEX solver [37] was used for both MILP problems (Table 1). Figure 3b presents the HIWAN hyperstructure constructed based on the results provided in Figure 3a. The HIWAN hyperstructure generation follows the same approach as the HEN hyperstructure [31], as shown in Figure 3a. The inlet of unit_u2has one cold thermal stream which has one match with hot utility_uH, one match with wastewater unit_uwwand one match with its outlet; hence, one superstructure with three matches is defined for this stream. The Wastewater unit_uww has two cold streams with one and three matches, respectively. This results in a hyperstructure consisting of two stream superstructures as illustrated in Figure 3b. As presented in Figure 1 for the proposed hyperstructure, mass flows can be split from source unit heat exchanger outlets and mixed at sink unit heat exchanger inlets. This has been illustrated in Figure 3b with thicker, colored arrows.

HIWAN synthesis is solved for problem P3 as well as for the HEN model of Floudas and Ciric [31]. SNOPT solver [38] is used to solve both models (Table 1). Its algorithm is based on sequential quadratic programming (SQP) and its solutions are locally optimal. It should be noted that problem P3 is non-convex and hence the global optimality of the solutions cannot be guaranteed. The solution of the NLP hyperstructure of Floudas and Ciric [31] is shown in Figure 4a. Since this solution is based on the results of HLD, all seven heat exchangers are active. The area of heat exchangers 3 (55.4 m²) and 7 (196.4 m² ) are relatively small compared to the others (Table 2). This approach, however, does not account for the possibility of stream mixing and splitting. The HIWAN hyperstructure addresses these possibilities by embedding the water allocation network within the HEN hyperstructure. Figure 4b presents the results of the proposed hyperstructure. By considering additional possibilities for stream mixing and splitting, hence increased possibility of non-isothermal mixing, three out of the seven targeted heat exchangers are determined not to be required. The flows in bold illustrate the new possibilities embedded in the proposed HIWAN hyperstructure. Consequently the HEN cost is reduced by 22.47 % (290,799 compared to 375,089 kUSD/yr). Table 2 provides the temperatures, flows and heat loads associated to heat exchangers in both cases.

6. Validation and Discussion

Two test cases from the literature were selected to validate the mathematical formulation and the proposed iterative sequential solution strategy. The cost data are presented in Table 3. The maximum number of integer cuts was limited to 50, though one sub-problem (P1 or P2) may become infeasible before reaching its respective limit, at which point the algorithm terminates the corresponding integer cut loop. In combination with the application of integer cut constraints, small minimum sizes were introduced for all water units in problem P1 (0.05 kg/s). Similarly, for problem P2, the minimum allowed heat exchange between hot stream i and cold stream j in pinch interval l is set at a small positive number (_Qi,j,lmin=0.1 kW in Equation (A6)). These two additions prevent generation of multiple identical solutions. Furthermore, to investigate the trade-off between investment cost and operating cost, the HRAT in problem P1 was varied between 1^∘C and 10^∘C with step of 1^∘C (N_HRAT= 10). These values are passed to problem P3 to be used as the value of ΔT_min. Test case I is a threshold problem and hence the value of HRAT does not affect the utility consumption. This, however, becomes a decisive factor for test case II where both hot and cold utilities are required. In this section, only solutions of the HIWAN hyperstructure are presented. A parallel coordinate visualization tool is used to illustrate all solutions and selected key performance indicators. These key performance indicators (KPI)s are listed here:

• Resource indicators: freshwater and thermal utility consumption,_m˙fw,_Q˙uH,_Q˙uC;

• Network indicators: N_sth(number of thermal streams), N_HE(number of heat exchangers), A_HENtotal(total area of heat exchangers), N_mixer(number of mixing points), N_mixerNIM(number of non-isothermal mixing points), N_sm(number of mass streams in the water network), and Q_HENtotal(total heat load of all heat exchangers);

• Economic indicators: C^HEN(HEN cost) and C^TAC(total annualized cost).

It should be noted that test case I was intentionally designed by the scientific community to validate their mathematical formulations and solution strategies. Due to the complexity of solving non-convex MINLP problems that arise in HIWAN synthesis, and solver tendency toward local optimality in such problems, it can be argued that these test cases are not easily solved without carefully manipulating solvers and their intrinsic options to reach a single, ´´good” solution. Thus, the main goal and contribution of the proposed solution strategy in this work is to solve the HIWAN problem by first reducing the search space, then generating many solutions to the problem. To this end, many solutions can be generated for problem P3 which can be seen as approaching an optimal point from different starting points (solutions of problems P1 and P2). 6.1. Test Case I—Single-Contaminant Problem

This test case in its optimal configuration [14,40,41] exhibits 90 kg/s of freshwater consumption and 3780 kW of hot utility consumption (minimum targets, 2397.0 kUSD/yr in operating cost). Table A1 provides the assumptions and solver options for this test case. Out of 26,010 potential solutions, 6462 solutions exist to problem P2—this arises from the fact that problem P1 converged for all integer cuts, while problem P2 was found to be infeasible for low HRAT values after relatively few integer cuts. 3611 solutions exist to problem P3. They exhibit a wide spectrum of solutions with HEN cost (the objective function in problem P3) ranging from 255.1 kUSD/yr (the minimum value within the literature [5]) to 660.9 kUSD/yr. Among these, 3506 solutions exhibit HEN cost in the range reported in the literature [5] (the maximum observed cost is 413.0 kUSD/yr). To visualize this set of solutions, several filters were imposed to further reduce the number of solutions:

• ^CHEN≤300 kUSD/yr;

• N_sm ≤ 15, maximum value reported in the literature for this test case [5];

• N_mixer ≤ 13, maximum value reported in the literature for this test case [5];

• N_mixerNIM ≤ 10, maximum value reported in the literature for this test case in [5];

Figure 5 presents 2247 solutions which remain after filtering (the color spectrum is based on HEN cost). Selected results from the literature [5] with HEN cost lower than 300 kUSD/yr are also plotted in this figure in black bold. The number of heat exchangers (N_HE ) varies between two and five, where two is the lowest number of heat exchangers reported for this test case by the proposed solution strategy. Figure 5a illustrates an inverse correlation between the number of heat exchangers and their total area. The total area (A_HENtotal) varies between 3350 m²and 7090.7 m². The total area of heat exchangers reported in the literature varies between 3111 m² (Wan Alwi et al. [42]) and 6300 m² . However, it should be noted that the HEN design reported by Wan Alwi et al. [42] is infeasible due to a violation of ΔT_min.

Figure 5b presents the solutions exhibiting the lowest HEN cost in the range of 255 kUSD/yr. This figure illustrates the vast possibility of solutions with the same minimum HEN cost (375 solutions in this case). Almost all solutions in the literature are present in the list of solutions generated by the sequential solution strategy. It should also be noted that many solutions exist for problem P3 which have the same KPIs and are therefore overlapping in Figure 5a,b.

6.2. Test Case II—Simplified Industrial Case Study

The simplified industrial case study introduced by Kermani et al. [14] is revisited here. Similar to previous case studies, problems P1 and P2 are solved for 50 integer cuts. Table A1 provides all assumptions and solver options for this test case. 3588 solutions exist for problem P3. All solutions reach 80 kg/s of freshwater consumption with 15,645 kW of cold utility. To visualize this set of solutions, several filters are imposed to further reduce the number of solutions:

• ^CHEN≤ 343 kUSD/yr (Table A2);

• N_sm≤ 20;

• N_mixer ≤ 17 (Table A2);

• N_mixerNIM ≤ 11 (Table A2);

Overall, 2053 solutions satisfy the above criteria which are visualized in Figure 6a (the color spectrum is based on HEN cost). The number of heat exchangers is limited to being 8 or 9, while the area varies between 2100 and 2750 m² . It can be observed from this figure that, for a given number of mixers and mass streams, the HEN cost varies over the entire range which enunciates that minimizing HEN cost and providing a single solution would not satisfy other important criteria. The lowest total cost reached by the proposed solution strategy is 960.6 kUSD/yr which is less than the lowest reported in the literature (971.4 kUSD/yr [13]). It should also be highlighted that the solution proposed by Ibrić et al. [13] considered the optimization of tank temperatures, while these temperatures are fixed for the analysis presented herein. Moreover, all proposed solutions in this work exhibit 11% lower freshwater consumption than the solutions reported in the literature. Figure 6b illustrates all of the solutions (240) proposed by the sequential solution strategy which exhibit lower total cost than 971.4 kUSD/yr. Figure 6b shows that for a given number of heat exchangers (8 for these solutions), increasing the number of non-isothermal mixing points lowers the total cost. Figure 7 presents one of the optimal HIWAN designs. Table 4 also provides a comparison between problems P2 and P3 in terms of thermal matches and their heat loads for the selected solutions. As described in Section 5, the HLD model and consequently HEN hyperstructure model do not consider stream mixing. The proposed HIWAN model takes these possibilities into account while optimizing the temperatures and flows of heat exchangers. This increases the possibility of replacing a heat exchanger by non-isothermal mixing which led to the elimination of the two smallest heat exchangers as highlighted in Table 4. In both cases, the total heat load of each non-water thermal stream remained constant, while the flow and temperature of water thermal streams adapted to address the shift in the heat load of non-water thermal streams.

7. Conclusions Two main complexities have been observed within the HIWAN methodologies in the literature. One difficulty is selecting the state of water thermal streams in the HEN design, as this choice affects the overall mathematical formulation in MINLP models. This, in turn, introduces the second difficulty of selecting thermal matches with undefined states of water thermal streams; thus, resulting in a complex MINLP formulation. These two complexities have been addressed by solving problems P1 (targeting model) and P2 (HLD) in sequence prior to the final HIWAN design which identify the states and potential thermal matches, respectively. A novel HIWAN hyperstructure with NLP formulation has been proposed in this work for the final design of HIWANs which exhibits higher degrees of non-isothermal mixing and lower investment costs. In addition, an iterative sequential solution strategy is proposed in this work by solving the sub-problems with integer cut constraints applied in both problems P1 and P2. As discussed, this approach can be viewed as a technique of reducing the search space at each iteration to a set of potential water thermal streams and a set of potential matches followed by solving the HIWAN synthesis problem for the reduced space.

Applying the proposed methodology on two test cases from the literature illustrated that it is not only able to reach the minimum total system cost (the only objective function used in the literature), but also generates a set of alternative “good” solutions exhibiting various results with respect to other key performance indicators. The first test case, a well-known problem in the literature, yielded 375 solutions with the same (minimum) value of the single objective used. This highlights the plethora of solutions which could be selected by decision-makers to achieve cost targets but with distinguishing features in non-objective or non-quantifiable goals. The second test case yielded 240 alternatives with lower cost than the minimum reported in the literature with an additional benefit of reducing the freshwater consumption by 11%. Similar to the first test case, this reinforces the approach by proposing many “good” solutions which can then be evaluated by engineers or decision-makers based on their other merits. Complexities related to application of the proposed solution strategy and associated MCDM tool on real industrial cases, as the one in a kraft pulp process [15], is limited to characterizing the industrial plant, its specifications, and collecting relevant data.

The proposed approach is intended for the grassroot design of HIWANs and cannot be directly applied to retrofit problems; however, as shown for problem P3, the operating conditions of a water allocation network (i.e., temperatures, flow rates, and heat exchange areas) can be optimized for a given set of thermal matches. Thus, for retrofit cases in which additional heat exchangers are not envisaged, the proposed superstructure can be applied at the current level. For investigating the addition of new heat exchangers, problem P2 can be reformulated by adding new constraints. Multi-period operation of HIWANs (by considering daily/seasonal variations in process operating conditions, temperatures of freshwater, etc.) and uncertainty analysis of HIWANs (to find the most resilient networks given system uncertainties, including costs and operating conditions) should be addressed in future work. The choice of HEN cost as the main objective in the literature cases should not be the sole measure of optimality, as the basic functions for heat exchanger cost reflect average, observed correlations for overall networks and hence should not be used as the decisive factor in proposing a final design. The proposed approach in this work overcomes this limitation by generating and analyzing a set of solutions and applying multi-criteria optimization to propose competing solutions considering multiple decision criteria.

Stage	Options	Value	Description
AMPL [39]	presolve_eps	10−5	maximum difference between lower and upper bounds in constraint violations
P1-MILP	mipgap	10−4	relative improvement in integer solution below which optimization is terminated
	integrality	10−10	a variable is considered to have an integral value if it lies within that of an integer
P2-MILP	mipgap	10−4
	integrality	10−4
P3-HEN-NLP [31]	feas_tol	10−5	satisfying upper and lower bounds of variables and linear constraints
P3-NLP	feas_tol	10−4

yi,j		tjex,i,I[∘C]	tjex,i,O [∘C]	tiex,j,I [∘C]	tiex,j,O [∘C]	LMTDi,j [∘C]	Qi,j [kW]	Ai,j [m2]	Ci,jinv [USD/yr]
1	uH-u2in (2)	120.00	120.00	81.96	99.96	28.1	3780.0	269.2	42,453
2	u2out-u2in	100.00 (-)	74.61 (-)	64.62 (-)	81.62 (-)	13.8 (-)	3570.0 (-)	518.9 (-)	59,078 (-)
3	u2out-u3in	100.00 (-)	77.70 (-)	64.96 (-)	74.96 (-)	18.2 (-)	504.0 (-)	55.4 (-)	21,342 (-)
4	uww1-ufw	50.00	30.00	20.00	40.00	10.0	2352.0	470.4	56,157
		(50.24)	(30.00)	(20.00)	(40.03)	(10.1)	(3371.4)	(667.3)	(67,398)
5	uww2-u2in	75.00	30.00	20.00	65.00	10.0	9450.0	1890.0	118,932
		(91.92)	(30.00)	(20.00)	(81.93)	(10.0)	(13,004.7)	(2602.2)	(142,397)
6	uww2-u3in	75.00	50.52	40.00	65.00	10.3	1260.0	245.7	40,613
		(75.01)	(55.24)	(40.00)	(57.86)	(16.2)	(1782.5)	(220.3)	(38,550)
7	uww2-ufw	50.20 (-)	30.20 (-)	20.00 (-)	39.87 (-)	10.3 (-)	1008.0 (-)	196.4 (-)	36,514 (-)
	Total							3645.9	375,089
								(3759.0)	(290,799)

(1) i and j represent hot and cold streams, respectively. Values in parentheses are results of HIWAN hyperstructure. (2) Same values in both cases.

Parameters	Units	Values for Cases Studies
		Test Case I	Test Case II
			(Simplified Kraft)
Cost of freshwater	USD/t	0.375	0.1525
Cost of cold utility	USD/kWyr	189	18.568
Cost of hot utility (steam, 120 ∘C)	USD/kWyr	377	140.16
Temperatures of cold utility	∘C	10 → 20	10 → 35
Heat-transfer coefficient	kW/m2K	1	1
Operating hours	hr/yr	8000	8322
Interest rate	%	8	8
Plant lifetime	yr	25	25
Temperature of freshwater	∘C	20	10
Temperature of wastewater	∘ C	30	30
Specific heat capacity of water	kJ/kg.K	4.2	4.2
ΔTmin of hot and cold utilities	∘ C	10	10
Heat exchanger cost coefficients (Equation (5)) [4]
α1	USD/yr	8000
α2	USD/yr/m0.6	1200
β	-	0.6

Thermal Match	Qi,j(Problem P2) [kW]	Qi,j(Problem P3) [kW]
surface condenser →uC	5418.1	5436.8
surface condenser →ufw70∘C	2141.9	2123.2
turpentine condenser →u2,in62∘C	671.9	666.9
turpentine condenser →ufw70∘C	10,248.1	10,253.1
contaminated condensate →uC	630.0	630.0
effluent →u1,out62∘C	105.0	0
effluent →ufw70∘C	2100.0	2205.0
dryer exhaust →ufw70∘C	1050.0	1050.0
u5,in30∘C→uC	213.9	0
uww62∘C→uC	9383.1	9578.3

Author Contributions

This study was done as part of M.K.'s doctoral studies supervised by F.M. I.D.K. contributed to this paper through reviewing, editing, and discussion of results.

Funding

This research project is financially supported by the Swiss Innovation Agency, Innosuisse, and the European Union's Horizon 2020 research and innovation programme under grant agreement No 818011.

Acknowledgments

This research project is financially supported by the Swiss Innovation Agency, Innosuisse, and is part of the Swiss Competence Center for Energy Research SCCER EIP and has received funding from the European Union's Horizon 2020 research and innovation program under grant agreement No 818011.

Conflicts of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Nomenclature

Sets
PinWAN	Set of process sink (ingoing) water unit operations
PoutWAN	Set of process source (outgoing) water unit operations
PH(PC)	Set of all the hot (cold) process streams
UH(UC)	Set of all the available hot (cold) utility streams
UoutWAN	Set of freshwater sources
UinWAN	Set of wastewater sinks
WANin	Set of all sink water unit operations, i.e., process and wastewater
WANout	Set of all source water unit operations, i.e., process and freshwater
WAN	Set of all water unit operations
MA	Set of thermal matches
HS (CS)	Set of all hot (cold) streams
HCS	Set of all thermal streams
C	Set of all contaminants
Si	Set of all water thermal streams related to water unit i
SMi,k	Set of streams matching with streamk∈Siof water uniti∈WAN
SMk	Set of streams,k', matching with streamk∈HCS, i.e., set of all streams,k', so that either of the two matches(k,k')or(k,k')is in the set of thermal matches,MA
Parameters
Tj	Temperature of unit j inPinWAN[∘C]
Ti	Temperature of unit i inPoutWAN[∘C]
cjk,max	Maximum allowed inlet contamination level of unit j inPinWAN[ppm]
cik,max	Maximum allowed outlet contamination level of unit i inPoutWAN[ppm]
Lu	Mass load of unit u to be removed from mass-transfer water unit operations [mg/s]
Ui,j	Heat transfer coefficient of match(i,j)inMA[kW/m2K]
cp	Specific heat capacity of water [kJ/kg.K]
Qk	Heat load of non-water thermal stream k inHCS[kW]
yk,lhld	Heat exchange match of hot streamk∈HSand cold streaml∈CS in problem P2 [-]
α1(α2)	Fixed (Area) cost parameter for heat exchangers [USD/yr]([USD/yr/m0.6])
β	Exponential cost parameter for heat exchangers [-]
NiccP1(NiccP2 )	Number of integer cut constraints applied to problem P1 (P2) [-]
niccP1(niccP2 )	Loop counter of the number of integer cut constraints applied to problem P1 (P2) [-]
nsol	Counter of the number of solutions [-]
Variables
Cinv	Total network investment cost [USD/yr]
CHEN	Heat exchanger network cost [USD/yr]
CTAC	Total annualized cost [USD/yr]
Qi,j	Heat load of match(i,j)inMA[kW]
Ai,j	Heat exchange area of match(i,j)inMA[m2]
LMTDi,j	Logarithmic mean temperature difference of match(i,j)inMA[∘C]
ΔT1i,j(ΔT2i,j)	Temperature differences at the two ends of heat exchanger of match(i,j)inMA[∘C]
m˙i(m˙j)	Flow rate of source (sink) unit i(j) inWANout(WANin) [kg/s]
m˙i,lk,I	Flow rate from source unit i inWANouttowards the inlet of heat match(k,l)in the superstructure of stream k inSi[kg/s]
m˙i,j,lk,I	Flow rate from source unit i inWANouttowards the inlet of heat match(k,l)in the superstructure of stream k of sink unit j inWANin[kg/s]
m˙i,lex,k(m˙j,lex,k)	Flow rate through heat exchanger of match(k,l)in the superstructure of stream k of source unit i inWANout(sink unit j inWANin) [kg/s]
m˙i,l',lB,k',k(m˙j,l',lB,k',k)	Flow rate from the outlet of heat exchanger of match(k,l)in the superstructure of stream k towards the heat exchanger inlet of match(k',l')in the superstructure of streamk'where both heat exchangers belong to source unit i inWANout(sink unit j inWANin) [kg/s]
m˙i,l,jk,O	Flow rate from the outlet of heat exchanger of match(k,l)in the superstructure of stream k of source unit i inWANouttowards sink unit j inWANin[kg/s]
m˙i,l,j,l'k,O,k',I	Flow rate from the outlet of heat exchanger of match(k,l)in the superstructure of stream k of source unit i inWANouttowards the inlet of heat exchanger of match(k',l')in the superstructure of streamk'of sink unit j of source unit i inWANin[kg/s]
ti,lex,k,I(ti,l'ex,k',O)	Temperature at the inlet (outlet) of heat exchanger of match(k,l)((k',l')) in the superstructure of stream k(k') of source unit i inWANout[∘C]
tj,lex,k,I(tj,l'ex,k',O)	Temperature at the inlet (outlet) of heat exchanger of match(k,l)((k',l')) in the superstructure of stream k(k') of sink unit j inWANin[∘C]
tlex,k,I(tlex,k,O)	Temperature at the inlet (outlet) of heat exchanger of match(k,l)in the superstructure of stream k inHCS[∘C]
m˙j,lk,O	Flow rate from the outlet of heat match(k,l)in the superstructure of stream k towards sink unit j inWANin[kg/s]
m˙i,j	Flow rate from source unit i inWANouttowards sink unit j inWANin[kg/s]
cjk,in	Contamination level of contaminant k at the inlet of sink unit j inWANin[ppm]
cik,out	Contamination level of contaminant k at the outlet of source unit i inWANout[ppm]
Fk	Heat capacity flow rate of thermal stream k inHCS[kW/K]
flk,I	Heat capacity flow rate of thermal stream k towards the inlet of heat exchanger of match(k,l)/(l,k)in the superstructure of stream k inHCS[kW/K]
flex,k	Heat capacity flow rate through heat exchanger of match(k,l)/(l,k)in the superstructure of stream k inHCS[kW/K]
flk,I	Heat capacity flow rate from the initial splitter of steam k towards heat exchanger of match(k,l)/(l,k)in the superstructure of stream k inHCS[kW/K]
flk,O	Heat capacity flow rate from the outlet of heat exchanger of match(k,l)/(l,k)towards the final mixer of stream k in the superstructure of stream k inHCS[kW/K]
fl,l'B,k	Heat capacity flow rate from the splitter following the match between k andl'to the mixer preceding the match between k and l within the superstructure of stream k inHCS[kW/K]
m˙fw	Freshwater flow rate (in the optimal solution) [kg/s]
Q˙uH(Q˙uC)	Heat load of hot (cold) utility (in the optimal solution) [kW]
Nsth	Number of thermal streams (in the optimal solution) [-]
Nsm	Number of mass streams in the water network (in the optimal solution) [-]
NHE	Number of heat exchangers (in the optimal solution) [-]
AHENtotal	Total area of heat exchangers (in the optimal solution) [m2]
QHENtotal	Total heat load of all heat exchangers (in the optimal solution) [kW]
Nmixer	Total number of mixers (in the optimal solution) [-]
NmixerNIM	Total number of non-isothermal mixing points (in the optimal solutions) [-]
Superscripts
I	Denoting the inlet of a heat exchanger
O	Denoting the outlet of a heat exchanger
B	Denoting the bypass stream of a heat exchanger
H (C)	Hot (Cold) stream
ex	Denoting a heat exchanger
Subscripts
fw	freshwater
ww	wastewater
sol	solution
icc	integer cut constraint

Abbreviations

The following abbreviations are used in this manuscript:

ΔTmin	Heat Exchanger Minimum Approach Temperature
HEN	Heat Exchanger Network
HIWAN	Heat-Integrated Water Allocation Network
HLD	Heat Load Distribution
HRAT	Heat Recovery Approach Temperature
ICC	Integer-Cut Constraint
ICCs	Integer-Cut Constraints
MCDM	Multi-Criteria Decision Making
MILP	Mixed-Integer Linear Programming
MINLP	Mixed-Integer Nonlinear Programming
NLP	Nonlinear Programming
TAC	Total Annualized Cost

Appendix A. Problem P1-Targeting

Problem P1 is formulated as an MILP model following the methodology proposed by [14]. The objective function is defined as minimizing the total annualized cost,^Ctac , of the system and is given by Equation (A1): (_PinWAN,_PoutWAN) (^UH,^UC,_UinWAN)

P1:min^Ctac=^Cop+^Cinv=∑fw∈_{UoutWANm˙fw Cfw}+∑ww∈_{UinWANm˙ww Cww}+∑u∈^UH_{f˙uH quH CuH}+∑u∈^UC_{f˙uC quC CuC}·^top+irr^{(1+irr)_nyear (1+irr)_nyear} ∑u∈^UH(_{yuH If,uH}+_{f˙uH Ip,uH})+∑u∈^UC(_{yuC If,uC}+_{f˙uC Ip,uC})+irr^{(1+irr)_nyear (1+irr)_nyear} ∑i∈WAN∑T∈_TiWAN(_yi,T·_If,i,TWAN+_m˙i,T·_Ip,i,TWAN)

where_m˙fwand_m˙wware the flow rates of freshwater sourcefwand wastewater sinkww, respectively,_Cfwand_Cwware the cost of freshwaterfwand wastewaterww(for treatment), respectively,_f˙uH,_quH, and_CuH(_f˙uC,_quC, and_CuC) are the flow rate, unit heat load and cost of hot (cold) utility u, respectively. The investment cost,^Cinvis annualized given the interest rate,irrand operating lifetime,_nyear, of the system,_If,i,TWANand_Ip,i,TWAN are the linearized coefficients of the penalizing factor. The optimization is subject to mass, energy, and contamination constraints. For details of the calculations of the penalizing cost and the constraints refer to Kermani et al. [14], Maréchal and Kalitventzeff [43], Kermani et al. [44]. The results provide water and energy targets as well as a feasible water network with a set of water exchanges and heat loads of all thermal streams. Having solved problem P1, the second problem can be formulated to find the promising thermal connections.

Appendix B. Problem P2-Heat Load Distribution

The set of all hot (^HShld) and cold (^CShld ) thermal streams in the network can be listed at this stage knowing the solution of problem P1. In addition, several pinch points will be identified that divide the temperature range into several subnetworks (SN ). At this stage, the HLD model will be applied to each subnetwork to find the potential thermal matches. HLD is formulated as an MILP model, based on minimizing the total number of heat exchange matches among all hot and cold streams [29,30]. Similar notations introduced by Papoulias and Grossmann [29] are used here for familiarity purposes.

It is important to note that the result of HLD is one feasible solution among many possible solutions that exhibit the minimum number of matches. Moreover, to keep the targeted thermal utilities intact, no thermal exchange is allowed crossing the pinch point, i.e., matches are limited between thermal streams within each subnetwork. For any subnetworkl∈SN, one can readily define the set of hot and cold thermal streams within each subnetwork (_HSlhld,_CSlhld, respectively)._HSl,rhldis defined as the set of hot streams within subnetwork l that are present in temperature interval r and above, while_CSl,rhldis defined as the set of cold streams within subnetwork l that are present in temperature interval r._Qi,j,r≥0is defined as the heat exchanged between hot stream i and cold stream j in intervalr∈_SNl. Upper and lower bounds can be defined for each match (_Qi,j,lmin,_Qi,j,lmax ) to better incorporate practical infeasibilities. To avoid violating thermal utility targets, one should follow the same constraints as in problem P1. A heat cascade formulation should be completed for each hot stream (_Ri,r≥0). By defining binary variables_yi,j,l∈{0,1}to denote the existence of a match between thermal streams within each subnetwork l, the HLD problem will minimize the number of heat exchangers. Ranking factors (_pi,j,l) can be defined to favor or avoid certain matches:

P2:min_Nhld=∑l∈SN∑i∈_HSlhld∑j∈_CSlhldpi,j,l·_yi,j,l

s.t.(∀l∈SN)_Ri,r=_Ri,r-1+_Qi,r-∑j∈_{CSl,rhldQi,j,r}∀i∈_HSl,rhld∀r∈_SNl

_Qj,r=∑i∈_{HSl,khldQi,j,r}∀j∈_CSl,rhld∀r∈_SNl

∑r∈_SNlQi,j,r≤_Qi,j,lmax·_yi,j,l∀i∈_HSlhld∀j∈_CSlhld

∑r∈_SNlQi,j,r≥_Qi,j,lmin·_yi,j,l∀i∈_HSlhld∀j∈_{CSlhldQi,j,lmax}≤min∑r∈_SNlQi,r,∑r∈_SNlQj,r

Appendix C. Assumptions and Solver Options for Test Cases

Both test cases are run on a Windows machine with Intel(R) Core(TM) i7-4600U CPU 2.10CHz 2.7 GHz and 8.00 GB RAM. The models are written in AMPL programming language. Table A1 provides the assumptions and solver options for test cases.

Table

Table A1.Assumption and solver options for test cases I and II.

Table A1.Assumption and solver options for test cases I and II.

Solver Option	Value (I)	Value (II)	Description (Mainly from Respective Solver's Manual)
AMPL
reset_initial_guesses	1	1	Reset variables to their initial values between different calls to solve command (all steps)
presolve_eps	variable	variable	Maximum difference between lower and upper bounds in constraint violations
			P1: 10-6(test case I), 10-4(test case II) and 10-5 for all others
presolve	1	1	Simplifying problem prior to the solver by fixing variables and dropping redundant constraints (all steps)
P1 (CPLEX)
mipgap	10-4	10-5	Relative tolerance for optimizing integer variables: stop if abs((best bound) - (best integer)) < mipgap * (1 + abs(best bound)).
integrality	10-10	10-4	Amount by which an integer variable can differ from the nearest integer and still be considered feasible.
P2 (CPLEX)
mipgap	10-7	10-5	See above
integrality	10-10	10-8	See above
P3 (initialization) (SNOPT)
▹ Problem P3 is solved by minimizing the sum of flows between heat exchangers
meminc	1	1	Increment to minimum memory allocation
iterations	106	106	Minor iteration limit
feas_tol	10-6	10-4	Minor feasibility tolerance for all variables and linear constraints
P3 and P3 HIWAN hyperstructure (SNOPT)
meminc = 1 iterations = 106feas_tol = 10-4(see above)

HRAT = ΔTmin for each iteration. Enthalpy values of thermal stream are rounded for HEN design to 10^-6 . Optimality gap of heat balance constraints (HIWAN) at mixing points is set at 0.01%. For problem P1, sum of mass exchanges is added to the objective function (Equation (A1)) to reduce the number of mass exchanges while minimizing the total cost.

Appendix D. Benchmarking of the Simplified Kraft Pulp Case Study

Table A2 provides the benchmarking results of test case II (Section 6.2).

Table

Table A2.Key performance indicators (KPIs) for the simplified kraft pulp case study.

Table A2.Key performance indicators (KPIs) for the simplified kraft pulp case study.

		Reference	Case a [13]3	Case b [13]3	Kermani et al. [10]4
Utility indicators
Freshwater	kg/s	137	89.904	89.904	80
Hot utility	kW	3392	0	0	0
Cold utility1	kW	14,270	14,813	14,813	15,645
	kg/s	135.8	141.1	141.1	149.0
Total water consumption2	kg/s	272.8	231.0	231.0	229.0
Waste outlet temperature	∘ C	54.8	30	30	30
Network indicators
Nsth	-	15	9	9	11
NHE	-	8	6	6	10
AHENtotal	m2	2977.9	2891.8	2481.1	2444.9
Nmixer(NmixerNIM)	-	9 (7)	11 (11)	9 (8)	17 (11)
Nsm	-	5	17	16	15
Financial indicators
Cop	kUSD/yr	1366	685.8	685.8	656.0
CHEN	kUSD/yr	0	301.1	285.6	342.4
CTAC	kUSD/yr	1341	986.9	971.4	998.4

¹-For reference case, cold utility is required to cool waste streams to 30^∘C (The cold utility was assumed as freshwater from 10^∘C to 35^∘C). A counter-current heat exchanger was assumed.²-Sum of freshwater and cooling water.³ -Ibrić et al. [13]: Case a is solved by fixing the temperatures of tanks at current values. Case b was solved by optimizing these temperatures: warm water tank at 50^∘C and hot water tank at 78.15^∘C.⁴-The solution with lowest C^HENis presented here (_niccP1= 1,_niccP2 = 41).