Systematic Selection of Green Solvents and Process Optimization for the Hydroformylation of Long-Chain Olefines

Abstract

Including ecologic and environmental aspects in chemical engineering requires new methods for process design and optimization. In this work, a hydroformylation process of long-chain olefines is investigated. A thermomorphic multiphase system is employed that is homogeneous at reaction conditions and biphasic at lower temperatures for catalyst recycling. In an attempt to replace the toxic polar solvent N,N-dimethylformamide (DMF), ecologically benign alternatives are selected using a screening approach. Economic process optimization is conducted for DMF and two candidate solvents. It is found that one of the green candidates performs similarly well as the standard benchmark solvent DMF, without being toxic. Therefore, the candidate has the potential to replace it.

1. Introduction

Since the 12 principles of “Green Chemistry” were postulated in 1998 [1], the concept has gained an evergrowing interest [2]. Nowadays, chemical engineering not only focuses on improving process efficiency by means of economic factors, but also focuses on ecological aspects. In the short term, ecologic and economic objectives are often contradictory [3]. However, in the long term, these objectives could very well coincide due to stricter regulations, e.g., on hazardous chemicals, social pressure [4], and limited resources. Solvents are a class of chemicals that finds frequent use in the chemical industry for diluting or dissolving other chemical species. Because they pose a potential health risk for personnel due to the many ways they may enter the body [5], the search for non-toxic alternatives for well-established toxic solvents is important.

Solvent properties may directly influence the process economics, for example due to a higher boiling point, which increases the cost of distillation, or a lower extraction efficiency. Hence, identifying economically well-performing solvents meeting additional non-economical requirements, like non-toxicity, is non-trivial. For this reason, computer-aided and optimization-based methods for solvent selection have been developed. Choosing a new solvent and simultaneously optimizing the process is called integrated process and solvent design. These design problems are usually formulated as complex mixed-integer non-linear programs (MINLPs). Owing to the difficulty in solving MINLPs, hierarchical approaches emerged and two different strategies are commonly used.

The first strategy relies on computer-aided molecular design (CAMD), where so-called group contribution methods [6] are used to estimate the properties of a molecule during optimization. Burger et al. [7] used a hierarchical approach named “HiOpt”, where a reduced process model is optimized using multi-objective optimization (MOO) to shrink the molecular search space for a subsequent MINLP optimization of the full process model. There are several other examples where a MOO was first used to select solvent candidates and followed by process optimization to select the best among those candidates [8,9]. Karunanithi et al. [10] and Cignitti et al. [11] decomposed the problem into easy-to-solve sub-problems by a step-wise reduction of the search space in conjunction with a step-wise increase in problem complexity. Eden et al. [12] optimized a process with the molecule properties left as additional degrees of freedom and afterward found the best fitting molecule using CAMD. A three-step approach for the solvent selection was proposed by Kossack et al. [13] for extractive distillation based on an initial screening, a short-cut model evaluation, and a rigorous process optimization. Gopinath et al. [14] developed a tailor-made outer approximation algorithm for CAMPD problems and applied it to the separation of methane and CO${}_2$ at high pressures. First et al. [15] performed a database screening for zeolites used in natural gas separation followed by a rigorous process optimization and determined the economically most beneficial process. In order to chemically store the fluctuating energy produced from renewable sources, Jens et al. [16] identified an optimal storage molecule from six candidates and an optimal solvent for a this process using a combination of screening, rigorous process models, and experimental data. Overall, several decomposition methods were developed with the objective to reduce the solvent space and make the solution of the MINLP tractable. An alternative to solving the MINLP deterministically is to use stochastic algorithms, such as genetic algorithms [17].

The second strategy relies on data-based screening methods for the selection of an appropriate solvent. First, unsuitable candidates from selected molecule databases are removed using MOO approaches. A subsequent process optimization is performed for the remaining candidates to find the economically best-performing solvent [18,19,20]. Scheffczyk et al. [21] introduced a method that evaluates a process using simplified pinch-based short-cut models with a large solvent set. They used very basic screening criteria, such as the existence of a miscibility gap, to avoid excessive exclusion of potentially interesting molecules. This, coupled with the simplified process models, allowed them to evaluate thousands of possible solvents and to identify the best candidates that minimize energy demand. This framework was applied to find novel solvents for the use in several extraction-distillation processes [21,22] and adsorption solvents in an integrated CO${}_2$ capture and utilization process producing CO [23]. The latter was expanded by Fleitmann et al. [24] to the initial adsorption of the CO${}_2$ from a natural gas stream. By including this first adsorption step, new and different optimal solvents were identified. This again emphasizes the importance of considering the overall process in solvent design.

Another approach to solve the integrated solvent and process design problem is continuous-molecular targeting. Here, a physically based equation of state (EOS) is used without specifying a particular solvent molecule. The EOS parameters for the solvent are left as optimization variables and determined by the solver. In a following step, a Taylor series expansion of the objective function is used to find real-world molecules with the best matching EOS parameters [25,26,27]. Lampe et al. [28] increased the applicability of this approach by using GPC-SAFT to design a molecule with EOS parameters that are nearest to those optimally found instead of depending on a database of given values.

In the present work, the second approach, database screening, is used. A modified version of a recently developed screening approach [29] is used to select solvent candidates. A subsequent process optimization is conducted for the ecologically benign candidates as well as for a toxic and commonly used benchmark solvent.

The case study considered here is the hydroformylation of n-decene. The two major challenges of the process are achieving a high space–time yield of the catalyzed reaction and the recovery of the rhodium based catalyst. To meet these challenges, a thermomorphic multiphase system (TMS) [30] consisting of the polar solvent N,N-dimethylformamide (DMF) and the non-polar solvent n-dodecane is employed. In the course of the collaborative research center “InPROMPT”, this process was already investigated extensively, albeit with a longer olefin chain length [31,32,33]. Several optimization studies were conducted [34,35,36] and the process was successfully operated continuously and optimized in real-time with a TMS based separation [37,38] and with a surfactant-based separation [39,40]. Some aspects of green chemistry are already fulfilled. Long-chain olefines are a class of raw materials that may be synthesized from renewable resources using the Fischer–Tropsch process with biogas or directly derived from unsaturated oleochemicals from renewable resources [41,42,43]. The reactant used in this case study, n-decene, is a representative of this group. Furthermore, the process is operated efficiently to save energy and thus reduce the CO${}_2$ footprint. However, one major drawback is the use of the polar solvent DMF. Despite its wide usage throughout the chemical industry, it is developmentally toxic [44] and also damaging to the liver [45]. DMF fulfills the criteria of Article 57 of the European chemicals ordinance REACH and is on the list of substances of very high concern (SVHC) [46]. Therefore, DMF should be replaced by a safer, ecologically benign alternative in order to obtain a green process.

The article is structured as follows: Section 2 will give an overview of the screening approach and the process model. Section 3 will present the results of our hierarchical process and solvent design. Finally, Section 4 will give a summary of the results and a short outlook for further research.

2. Materials and Methods

Our approach to the integrated process and solvent design can be split into two parts. The first part is the screening of databases to exclude unsuitable solvents from further considerations and the second part is a process optimization with the remaining solvent candidates and the benchmark solvent DMF. Both steps and their necessary prerequisites are described in this section.

2.1. Solvent Screening

The screening methodology identifies promising catalyst solvents by searching a database for desired properties. The approach was published recently as a conference paper [29]. In summary, the screening starts by taking into account all components within the COSMObase database. This database provides geometry optimized molecular structures in the COSMO state, which can be used for thermodynamic property predictions in the liquid state using the continuum solvation model for realistic solvation (COSMO-RS) [47,48]. The first screening step takes physical properties into account. Among others, these are structural constraints to guarantee stability during reaction (e.g., carbon double bonds are not allowed), a reasonable molecular weight, and an appropriate boiling point in order to avoid azeotropes in solvent recovery distillation. The second screening step looks at environmental, health, and safety (EHS) properties. This is important, because [49] identified the developmentally toxic dimethylformamide (DMF) as the best performing catalyst solvent using a systematic solvent screening approach. The current work expands upon the former screening approach by including six EHS criteria in order to identify an environmentally more benign solvent to replace DMF. The evaluated EHS properties are persistence, bioaccumulation, toxicity, carcinogenicity, mutagenicity, and developmental toxicity. QSPR models published in VEGA [50] and EPISuite [51] are used to predict these properties for the candidate solvents. All solvent candidates, which are considered hazardous with respect to these properties, are rejected, and those remaining are considered to be environmental friendly. This list of candidates is also augmented by green solvents identified in the literature [52] in order to validate the screening.

Then, the thermodynamic performance of each of the green solvent candidates is evaluated. The relative solubility is a measure for the solvent’s ability to extract the catalyst, and potential solvents were rejected if their relative solubility is below a certain threshold. Finally, the quality of the miscibility gap is investigated at separation conditions (T = 25 °C) using n-decane as the non-polar solvent, as done in [49]. The ternary liquid–liquid equilibrium (LLE) for the system composed of the two solvents and the product n-undecanal were then predicted. Systems with less than one LLE point were excluded. In order to guarantee desirable TMS phase behavior, the miscibility gap between the solvents and the reactant 1-decene was calculated at the reaction temperature of 100°C as a final screening step. All systems showing a heterogeneous phase domain are removed.

As a result, three potential polar solvents were identified. Two could be purchased commercially and their LLE behaviors and reaction feasibility were investigated experimentally. Both candidates, dimethyl succinate (DSUC) and tetrahydro-4H-pyran-4-one (THPO), performed successfully as reaction solvents and showed the desired phase behavior. The systematic process-wide comparison of both candidates and the state of the art solvent DMF are now investigated for separation.

2.2. Process Model

The flowsheet of the considered process is shown in Figure 1 and can be roughly split into three parts: The reactor, the decanter cascade, and the black-box representing additional process parts. These parts will be presented in the following. The chemical species involved in the considered process can be found in

Table 1. Chemical species. Species highlighted in grey are explicitly modeled in the reactor only and lumped otherwise. Only one of the solvents highlighted in green is present at any given time.

Index	Name	Purpose	Produced by
C10en	n-decene	reactant
iC10en	iso-decene	side product	isomerisation of n-decene (reversible)
C11al	n-undecanal	desired product	hydroformylation of n-decene (irreversible)
iC11al	iso-undecanal	side product	hydroformylation of iso-decene (irreversible)
C12an	n-dodecane	non-polar solvent
DMF	dimethylformamide	polar solvent
DSUC	dimethyl succinate	polar solvent
THPO	tetrahydropyranone	polar solvent
C10an	n-decane	side product	hydrogenation of n-/iso-decene (irreversible)
RHO	rhodium	catalyst
BIP	biphephos	ligand
H${}_2$	hydrogen	reactant
CO	carbon monoxide	reactant

.

2.2.1. Reactor

The reactor model is taken from [32,33]. The model describes a hydroformylation reaction in a continuous stirred tank reactor. CO and H${}_2$ are provided to the homogeneous liquid reaction phase from a second, gaseous phase in the reactor. Other reactants and additives are provided in liquid form. Deviating from the aforementioned publications, we only consider the stationary case, and the equilibrium between the liquid and the gas phase is described by an artificial neural network (ANN) trained on the PC-SAFT equation of state [34]. The activation function of the ANN, a hyperbolic tangent, is replaced by a piecewise linear approximation [34]. Furthermore, the chemical species are different. The chain length of the olefin compared to the previous work is reduced by 2 and the new solvents are introduced. A list of all species involved in the reaction can be found in Table 1. However, the original reaction kinetics were used because the deviations due to the different chain lengths are small and the influence of the polar solvent on the reaction is negligible [49].

The mass balance is defined as

$$\begin{array}{cc}\hfill F_{\mathrm{reac},i}^{\mathrm{out}}& =F_{\mathrm{reac},i}^{\mathrm{in}}+V_{\mathrm{reac}}{c}_{\mathrm{cat}}{M}_{\mathrm{cat}}\sum _{j\in \mathrm{RCT}}{\nu }_{i,j}{r}_j,\hfill \end{array}$$

(1)

$$\begin{array}{cc}\hfill F_{\mathrm{reac},i}^{\mathrm{in}}& =Z_i+F_{\mathrm{product},i}^{\neg \mathrm{C}11\mathrm{al}}(1-\varpi )+F_{\mathrm{recycle},i},\hfill \end{array}$$

(2)

where $c_{\mathrm{cat}}$ is the active catalyst concentration inside the reactor, $M_{\mathrm{cat}}$ is the molar mass of the catalyst, $\mathit{\nu }$ and $\mathit{r}$ are the stoichiometric coefficients and reaction rates, respectively, for the eight reactions, RCT, taking place. The general form of the reaction rates $\mathit{r}$ of all eight reactions in the reactor is given by

$$\begin{array}{c}\hfill r_j=k_j·\left(\prod _i\left(c_i^{{\nu }_{j,i}}\right)\right)·{\left(1+\sum _i\left(K_{i,1}·c_i^{K_{i,2}}\right)\right)}^{-1}.\end{array}$$

(3)

The reaction rate constants $\mathit{k}$ are calculated using Arrhenius law and the parameters $\mathit{K}$ are taken from the literature and can be found in [33]. To calculate the concentration of the active catalyst $c_{\mathrm{cat}}$ the following equation is used

$$\begin{array}{c}\hfill c_{\mathrm{cat}}={\displaystyle \frac{c_{\mathrm{cat},\mathrm{tot}}}{1+K_{\mathrm{cat},1}·c_{\mathrm{CO}}^{K_{\mathrm{cat},2}}}},\end{array}$$

(4)

where $c_{\mathrm{cat},\mathrm{tot}}$ is the total catalyst concentration.

The reactor inlet ${\mathit{F}}_{\mathrm{reac}}^{\mathrm{in}}$ is composed of three elements. The first of these is the feed stream, $\mathit{Z}$, which contains the make-up streams for the catalyst, the solvents, and the reactants. The second is the recycle stream, ${\mathit{F}}_{\mathrm{recycle}}$, coming from the decanter cascade, which recycles the polar solvent (thereby the catalyst), unused reactants, and some of the products. The final stream is the recycle coming from the black box, ${\mathit{F}}_{\mathrm{product}}^{\neg \mathrm{C}11\mathrm{al}}$, minus the purged stream.

2.2.2. Downstream Processing

The reactor outlet, ${\mathit{F}}_{\mathrm{reac}}^{\mathrm{out}}$, is fed into the downstream section consisting of four decanters and a distillation column.

The distillation column is required to recover the extraction solvent. It is assumed to be operated at a pressure of 60 mbar, with a linear pressure drop of 30 mbar along the column. The requirement of the column is that it recovers all of the circulated extraction solvent $F^{\left(psol\right)}$ at the top. Instead of a rigorous distillation column model, a separately calculated surrogate function is employed for each investigated polar solvent $psol$. The procedure to obtain the surrogate function can be split into two parts.

First, the relative volatilities $\mathit{\alpha }$ are calculated. An ensemble of 810 distillation columns with a different number of stages, reflux, and reboil ratios were simulated. Each distillation column is assumed to have constant molar flows, a constant molar overflow, and a total condenser. Furthermore, each is assumed to behave ideally, thus the gas phase is described by Dalton’s law ${\mathit{P}}_{\mathrm{vap}}=\mathit{\chi }·P$ and the liquid phase is described by Raoult’s law ${\mathit{P}}_{\mathrm{vap}}=\mathsf{{\rm Y}}·{\mathit{P}}_{\mathrm{vap}}^{*}$, where $\mathit{\chi }$ and $\mathsf{{\rm Y}}$ are the liquid and vapor mole fractions, respectively. For the calculation of the vapor pressures the vapor pressure correlation

$$\begin{array}{c}\hfill {log}_{10}\left({\mathit{P}}_{\mathrm{vap}}^{*}\right)={\mathit{a}}_0+{\mathit{a}}_1·T^{-1}+{\mathit{a}}_2·{log}_{10}\left(T\right)+{\mathit{a}}_3·T+{\mathit{a}}_4·T^2\end{array}$$

(5)

is used, with T in K and ${\mathit{P}}_{\mathrm{vap}}^{*}$ in mmHg. The parameters ${\mathit{a}}_i$ can be found in

Table A1. Vapor pressure parameters (${10}^{a_0+a_1{T}^{-1}+a_2{log}_{10}\left(T\right)+a_{3}T+a_4{T}^2}$) [mmHg]. Vapor pressure correlations for DSUC and THPO are fitted with the method presented in [62], the other values are taken from [63].

	a${}_0$	a${}_1$	a${}_2$	a${}_3$	a${}_4$
DMF	−47.9857	$-2.385$ × ${10}^3$	28.8	$-5.8596$ × ${10}^{-2}$	$3.1386$ × ${10}^{-5}$
DSUC	117.8014	$-6.3944$ × ${10}^3$	−42.5731	$3.0869$ × ${10}^{-2}$	$-9.2995$ × ${10}^{-6}$
THPO	74.2227	$-4.2846$ × ${10}^3$	−25.9627	$1.8373$ × ${10}^{-2}$	$-5.4143$ × ${10}^{-6}$
C10en	2.2678	$-3.12$ × ${10}^3$	5.43	$-2.01$ × ${10}^{-2}$	$1.12$ × ${10}^{-5}$
C12an	−8.5899	$-3.5241$ × ${10}^3$	$10.806$	$-2.8161$ × ${10}^{-2}$	$1.4267$ × ${10}^{-5}$
C11al	−31.8129	$-3.14$ × ${10}^3$	20.4	$-3.73$ × ${10}^{-2}$	$1.75$ × ${10}^{-5}$

. Using these 810 simulations, mean temperatures for the top of the column, ${\overline{T}}_D$, and the bottom of the column, ${\overline{T}}_B$, can be calculated. These mean temperatures are used to calculate the relative volatilities using Equation (5) and

$$\begin{array}{c}\hfill {\alpha }_i=\sqrt{{\displaystyle \frac{{\overline{P}}_{\mathrm{vap},\mathrm{B},i}^{*}}{{\overline{P}}_{\mathrm{vap},\mathrm{B},\mathrm{HK}}^{*}}}·{\displaystyle \frac{{\overline{P}}_{\mathrm{vap},\mathrm{D},i}^{*}}{{\overline{P}}_{\mathrm{vap},\mathrm{D},\mathrm{HK}}^{*}}}},\end{array}$$

(6)

where $\mathrm{HK}$ denotes the highest boiling component, i.e., C11al. The obtained values for $\mathit{\alpha }$ can be found in

Table A3. Relative volatilities.

	${\mathit{\alpha }}_{\mathit{DMF}}$	${\mathit{\alpha }}_{\mathit{DMS}}$	${\mathit{\alpha }}_{\mathit{THPO}}$
DMF	31.9942	[-]	[-]
DMS	[-]	4.0208	[-]
THPO	[-]	[-]	22.4438
C12an	3.1307	2.7766	3.0696
C10en	16.8669	11.7308	15.8739
C11al	1	1	1

.

The optimal costs for the distillation are calculated using the constant $\mathit{\alpha }$. Global optimizations using BARON [53], with a Fenske–Underwood–Gilliland (FUG) shortcut model [54], were conducted with 1000 varied feeds for each solvent. For these global optimizations the constraints are as follows:

• DMF: 99% purity of DMF in the distillate
• DSUC: 99% purity of DSUC + C10en in the distillate
• THPO: 99% purity of THPO in the distillate

The varied purity requirements for the different polar solvent candidates are necessary because of their different separation properties. These properties will be explained later.

The total annualized costs (TAC) in $\$/\mathrm{a}$ for the distillation column are calculated with a lumped cost function [55]

$$\begin{array}{cc}\hfill J_{\mathrm{col}}& ={\kappa }_1·V_{\mathrm{col}}+{\kappa }_2·{(l_{\mathrm{col}}+{\kappa }_3)}^{{\kappa }_4}·{\left(V_{\mathrm{col}}·{\kappa }_5\right)}^{{\kappa }_6}+{\kappa }_7·{(l_{\mathrm{col}}+{\kappa }_8)}^{{\kappa }_9}·{\left(V_{\mathrm{col}}·{\kappa }_{10}\right)}^{{\kappa }_{11}},\hfill \end{array}$$

(7)

using the parameters from

Table A2. Cost function parameters for each investigated solvent.

Parameter	DMF	DSUC	THPO
${\kappa }_1$	$17,764$	15,071	17,783
${\kappa }_2$	$2463.6$	2988.6	$2202.3$
${\kappa }_3$	$7.2$	$7.2$	$7.2$
${\kappa }_4$	$0.81$	$0.81$	$0.81$
${\kappa }_5$	1	1	1
${\kappa }_6$	$0.525$	$0.525$	$0.525$
${\kappa }_7$	$155.4357$	202.9541	$133.1394$
${\kappa }_8$	$-0.8$	$-0.8$	$-0.8$
${\kappa }_9$	$0.97$	$0.97$	$0.97$
${\kappa }_{10}$	1	1	1
${\kappa }_{11}$	$0.725$	$0.725$	$0.725$
${\lambda }_1$	64,191	168,833	51,628
${\lambda }_2$	20,638	19,656	20,695
${\lambda }_3$	−56,984	−182,134	−53,368
${\lambda }_4$	−9.6219	−4.2671	−12.5046
${\lambda }_5$	98,841	439,840	71,678
${\eta }_{\mathrm{C}}10\mathrm{en}$	$0.21041$ × ${10}^{-3}$	$0.21041$ × ${10}^{-3}$	$0.21041$ × ${10}^{-3}$
${\eta }_{psol}$	$0.07915$ × ${10}^{-3}$	$0.12395$ × ${10}^{-3}$	$0.08491$ × ${10}^{-3}$
${\eta }_{\mathrm{C}}12\mathrm{an}$	$0.233642$ × ${10}^{-3}$	$0.233642$ × ${10}^{-3}$	$0.233642$ × ${10}^{-3}$
${\eta }_{\mathrm{C}}11\mathrm{al}$	$0.195184$ × ${10}^{-3}$	$0.195184$ × ${10}^{-3}$	$0.195184$ × ${10}^{-3}$
${\theta }_1$	504,155	504,155	504,155
${\theta }_2$	0.586667	0.586667	0.586667
${\theta }_3$	$2.9661$ × ${10}^{12}$	$2.9661$ × ${10}^{12}$	$2.9661$ × ${10}^{12}$
${\theta }_4$	85,150	85,150	85,150
${\theta }_5$	17,248.42	17,248.42	17,248.42
${\theta }_6$	0.62	0.62	0.62
${\theta }_7$	$2.8512$ × ${10}^7$	$2.8512$ × ${10}^7$	$2.8512$ × ${10}^7$
${\theta }_{8,\mathrm{C}10\mathrm{en}}$	0.661	0.661	0.661
${\theta }_{8,psol}$	0.0731	0.0731	0.0731
${\theta }_{8,\mathrm{C}12\mathrm{an}}$	0.0714	0.0714	0.0714

, where $V_{\mathrm{col}}$ is the vapor flow rate and $l_{\mathrm{col}}$ is the length of the column.

The values of $J_{\mathrm{col}}$ are used for a polynomial regression of the costs with respect to the circulated extraction solvent stream $F^{\left(psol\right)}$ and the recovery rate $rec=F^{\left(psol\right)}/F_{\mathrm{col},psol}^{\mathrm{in}}$,

$$\begin{array}{c}\hfill TA{C}_{\mathrm{col}}={\lambda }_{psol,1}+{\lambda }_{psol,2}·F^{\left(psol\right)}+{\lambda }_{psol,3}·rec+{\lambda }_{psol,4}·{\left(F^{\left(psol\right)}\right)}^2+{\lambda }_{psol,5}·re{c}^2,\end{array}$$

(8)

with parameters $\mathit{\lambda }$ given in Table A2. Equation (8) is used as a surrogate for the distillation column.

The number of decanters is chosen because it was cost optimal in a previous study [56]. Following the results of [56] once again, the decanters are operated at a fixed temperature of $298.15$ K. The decrease in temperature compared to the reaction conditions leads to a liquid–liquid phase split. The mass balances are defined as follows,

$$\begin{array}{c}\hfill {\mathit{F}}_{\mathrm{dec},n}^{\mathrm{in}}=\left\{\begin{array}{cc}{\mathit{V}}_2+{\mathit{F}}_{\mathrm{col}}^{\mathrm{out}},\mathrm{D},\hfill & \mathrm{for}n=1,\hfill \\ {\mathit{L}}_{n-1}+{\mathit{V}}_{n+1},\hfill & \mathrm{for}1<n<n_{\max },\hfill \\ {\mathit{F}}_{\mathrm{reac}}^{\mathrm{out}},\hfill & \mathrm{for}n=n_{\max },\hfill \end{array}\right.\end{array}$$

(9)

where the decanters are numbered from right to left, i.e., $[n_{\max },\cdots ,1]$. The upper outlet of the decanters $\mathit{V}$ is the less polar, product-rich phase with composition $\mathit{y}$ and the lower outlet of the decanters $\mathit{L}$ is the more polar, catalyst-rich phase with composition $\mathit{x}$.

The boiling points of the three solvents vary, with DMF being the easiest to separate and DSUC being the hardest to separate (see the relative volatilities $\mathit{\alpha }$ in Table A3). Therefore, the separation properties of the three solvents regarding distillation vary. With DMF and THPO, it is possible to solely recover the solvent in the distillate stream. With DSUC, C10en becomes the component with the highest relative volatility. Therefore, we assume that all of the C10en contained in the column feed will leave through the distillate.

Because the components contained in the distillate flow differ depending on the separation properties of the solvents, the value of $F_{\mathrm{col},i}^{\mathrm{out}},\mathrm{D}$ needs to be calculated as

$$\begin{array}{c}\hfill F_{\mathrm{col},i}^{\mathrm{out}},\mathrm{D}={\beta }_i^{\left(psol\right)}{F}_{\mathrm{col},i}^{\mathrm{in}},\end{array}$$

(10)

because the column itself is replaced by the aforementioned surrogate function in Equation (8). Here, ${\mathit{\beta }}^{\left(psol\right)}$ is introduced as a vector of split factors for ($psol$, C12an, C10en, C11al, cat) to capture the different separation behaviors of the mixture with the respective solvent:

$$\begin{array}{cc}\hfill {\mathit{\beta }}^{\left(\mathrm{DMF}\right)}& =({\varsigma }_{psol},0,0,0,0),\hfill \end{array}$$

(11)

$$\begin{array}{cc}\hfill {\mathit{\beta }}^{\left(\mathrm{DSUC}\right)}& =({\varsigma }_{psol},0,1,0,0),\hfill \end{array}$$

(12)

$$\begin{array}{cc}\hfill {\mathit{\beta }}^{\left(\mathrm{THPO}\right)}& =({\varsigma }_{psol},0,0,0,0).\hfill \end{array}$$

(13)

Here, ${\varsigma }_{psol}$ is used to make sure that only the circulated extraction solvent stream $F^{\left(psol\right)}$ leaves the distillation column through the distillate.

Note that the fourth and the third decanter are only connected through the less polar flow of the fourth decanter. By doing this, the extraction solvent stream of the cascade bypasses the fourth decanter, allowing for an additional degree of freedom. The molar flow entering the column, ${\mathit{F}}_{\mathrm{col}}^{\mathrm{in}}$, is composed of the more polar phase of the third decanter, ${\mathit{F}}_{\mathrm{col}},2^{\mathrm{in}}={\mathit{L}}_3$, and a fraction of the more polar phase of the fourth decanter, ${\mathit{F}}_{\mathrm{col}},1^{\mathrm{in}}={\mathit{L}}_4$,

$$\begin{array}{c}\hfill {\mathit{F}}_{\mathrm{col}}^{\mathrm{in}}={\mathit{L}}_4(1-\varsigma )+{\mathit{L}}_3.\end{array}$$

(14)

That fraction is determined by the split factor $\varsigma \in [0,1]$. A large split factor leads to a smaller flow into the column. However, the composition of the mixture at the bottom of the column, and hence the temperature, also depends on the split factor. The literature suggests that the aldehyde will react with itself to the unwanted side product aldol at temperatures above $403.15$ K [57]. Therefore, the temperature at the bottom of the column is restricted to $388.15$ K, the same temperature as inside the reactor, which should be safe. The split factor $\varsigma$ can be used by the optimizer to achieve this temperature.

The catalyst inside the decanter cascade is described by the catalyst mass balance

$$\begin{array}{c}\hfill F_{\mathrm{reac},\mathrm{cat}}^{\mathrm{out}}=L_{n_{\max }}^{\mathrm{cat}}+L_{n_{\max }-1}^{\mathrm{cat}}+V_1^{\mathrm{cat}},\end{array}$$

(15)

where ${\mathit{L}}^{\mathrm{cat}}$ is the catalyst in the more polar phase and ${\mathit{V}}^{\mathrm{cat}}$ is the catalyst in the less polar phase. The partitioning of the catalyst between the two phases is calculated with

$$\begin{array}{c}\hfill {log}_{10}{P}_{\mathrm{yx},n}={log}_{10}\left({\displaystyle \frac{V_n^{\mathrm{cat}}{L}_n}{V_n{L}_n^{\mathrm{cat}}}}\right),\end{array}$$

(16)

where ${log}_{10}{P}_{\mathrm{yx},n}$ is the logarithmic partition coefficient for stage n.

The activity coefficients $\mathit{\gamma }$ for the LLE are described by a modified UNIFAC (Dortmund) [58] implementation. The corresponding LLE conditions,

$$\begin{array}{c}\hfill x_i{\gamma }_i^{\left(I\right)}-y_i{\gamma }_i^{\left(II\right)}=0,\end{array}$$

(17)

are solved with a residuum less than 1 × 10⁻⁹. For each of the LLE points, the partition coefficients ${\mathit{P}}_{yx}$ for the catalyst/ligand is calculated with COSMOtherm [47], because the size of the catalyst/ligand complex prohibits the use of group contribution methods such as UNIFAC [59]. To determine the composition of the two liquid phases and the catalyst distribution, a surrogate model $\widehat{f}$ is employed,

$$\begin{array}{c}\hfill [x_{i,n},y_{i,n},{log}_{10}{P}_{\mathrm{yx},n}]=\widehat{\mathit{f}}(t_{1,n},t_{2,n}),\end{array}$$

(18)

where $(t_{1,n},t_{2,n})\in {[0,1]}^2$ is a parameterization of the binodal curve. To generate the input/output data for fitting the surrogate model, a parameter continuation [60] in the orthogonal direction of the tie lines of the LLE is used.

The surrogate has a two-dimensional input space and a nine-dimensional output space. It consists of a second-degree polynomial and an ANN with a single hidden layer fit to the residuals of the polynomial. The ANN is generated with the “train” command of MATLABs deep learning toolbox. A more thorough description of the algorithm used to generate the surrogate is not the scope of this article and can be found in another publication of the authors in the present special issue of Processes [61].

A previous optimization study yielded that the catalyst recovery plays an important role for the overall process cost [56]. Therefore, the accuracy of the employed surrogate with regard to the partition coefficient ${\mathit{P}}_{\mathrm{yx}}$ is of high importance. The maximum errors of the three surrogate models for ${log}_{10}{\mathit{P}}_{\mathrm{yx}}=f_{psol}^{{\mathrm{P}}_{\mathrm{yx}}}(t_1,t_2)$ at the sample points are:

$$\begin{array}{l}\phantom{a}{\overline{\delta }}_{psol}=max\left(f_{psol}^{{\mathrm{P}}_{\mathrm{yx}}}(t_1,t_2)-{\widehat{f}}_{psol}^{{\mathrm{P}}_{\mathrm{yx}}}(t_1,t_2)\right),\\ \phantom{a}{\underset{¯}{\delta }}_{psol}=min\left(f_{psol}^{{\mathrm{P}}_{\mathrm{yx}}}(t_1,t_2)-{\widehat{f}}_{psol}^{{\mathrm{P}}_{\mathrm{yx}}}(t_1,t_2)\right),\\ {\overline{\delta }}_{\mathrm{DMF}}\phantom{\mathrm{a}}=0.001087,\hspace{1em}{\underset{¯}{\delta }}_{\mathrm{DMF}}\phantom{\mathrm{a}}=-0.000684,\\ {\overline{\delta }}_{\mathrm{DSUC}}=0.000524,\hspace{1em}{\underset{¯}{\delta }}_{\mathrm{DSUC}}=-0.001201,\\ {\overline{\delta }}_{\mathrm{THPO}}=0.002302,\hspace{1em}{\underset{¯}{\delta }}_{\mathrm{THPO}}=-0.001656.\end{array}$$

(19)

2.2.3. Black Box

The process considered in this work is the first step in an attempt to optimize the hydroformylation simultaneously with a reductive amination (RA) step. In the RA the product of the hydroformylation (C11al) and diethylamine are used to generate the corresponding amine, which is a valuable product used as a building block for agrochemical and pharmaceutical chemicals as well as surfactants, coatings, and lubricants. Although we currently do not have a rigorous process model for the reductive amination, the underlying conditions influencing the hydroformylation are known from the design specifications for the reductive amination process. Therefore, it is modeled as a black box with some simplifying assumptions. We assume that the aldehyde is separated from the other species in the cause of the RA and that the other species are recycled to the hydroformylation process. The catalyst present in the product stream is assumed to be lost during this separation process.

2.2.4. Overall Process

In order to connect the process parts described above efficiently, some further considerations had to be taken into account.

Eleven species are necessary to describe the reactions and interactions inside the reactor, see Table 1. Not every species important for the reaction plays a role in the extraction cascade. The isomerized C10en and C11al can be seen to behave identically to their respective non-isomerized counterparts during the extraction. Therefore, the isomerized and non-isomerized species are treated as if they were a single species in the extraction, thereby reducing the number of species by two. Furthermore, C10an and C12an can both be used as a non-polar solvent for the extraction, and are therefore combined into one non-polar solvent stream assumed to behave like C12an since the amount of C12an is clearly higher. Furthermore, it is also assumed that only the liquid phase of the reactor goes into the extraction cascade. With these assumptions, the number of species in the extraction cascade can be reduced to five. The reactor outlet ${\mathit{F}}_{\mathrm{reac}}^{\mathrm{out}}$ then consists of the following species,

$$\begin{array}{c}\hfill {\mathit{F}}_{\mathrm{reac}}^{\mathrm{out}}=\left(\begin{array}{c}psol\\ \mathrm{C}12\mathrm{an}\\ \mathrm{C}10\mathrm{en}\\ \mathrm{C}11\mathrm{al}\\ \mathrm{cat}\end{array}\right),\end{array}$$

(20)

where $psol\in [\mathrm{DMF},\mathrm{DSUC},\mathrm{THPO}]$ is one of the investigated solvents and $\mathrm{cat}$ is the catalyst/ligand complex.

To keep the whole model consistent, the ratios of the combined species in the extraction cascade are assumed to be identical for the reactor outlet and the recycle streams, i.e., n-/iso-ratios and C10an/C12an ratio.

2.3. Process Optimization

The process optimizations in this work are mixed integer nonlinear programs (MINLPs) of the following form:

$$\begin{array}{c}\hfill \begin{array}{cc}\underset{\tilde{\mathit{x}}}{min}\hfill & J\left(\tilde{\mathit{x}}\right),\hfill \\ \mathrm{s}.\mathrm{t}.\hfill & \mathit{h}\left(\tilde{\mathit{x}}\right)=0,\hfill \\ \hfill & \mathit{g}\left(\tilde{\mathit{x}}\right)\le 0,\hfill \\ \hfill & \tilde{\mathit{x}}\in G,G\subseteq {\mathbb{R}}^n,\hfill \\ \hfill & {\tilde{x}}_i\in Z,Z\subseteq \mathbb{Z},\mathrm{for}\mathrm{all}i\in I,\hfill \end{array}\end{array}$$

(21)

where $J\left(\tilde{\mathit{x}}\right)$ is the objective function, $\mathit{h}\left(\tilde{\mathit{x}}\right)$ are the equality constraints, such as mass balances, and $\mathit{g}\left(\tilde{\mathit{x}}\right)$ are the inequality constraints, such as purity requirements. The variables $\tilde{\mathit{x}}$ include the degrees of freedom, of these variables ${\tilde{x}}_i$ with $i\in I$ have to fulfill integrality restrictions.

The degrees of freedom are the feed flows $\mathit{Z}$ of the reactant, the catalyst, and the solvents, the hold-up time in the reactor, characterized by the reactor volume, $V_{\mathrm{reac}}\in [2.1853,4000]$ m${}^3$, the reactor temperature, $T_{\mathrm{r}}\in [368.15,388.15]$ K, and pressure, $P_{\mathrm{r}}\in [10,20]$ bar, the amount of extraction solvent used in the extraction cascade, $F^{\left(psol\right)}\in [7,90]$ mol/s, the split factor $\varsigma \in [0,1]$, and the purge $\varpi \in [0.0001,0.01]$. Note, that the model would yield an NLP without the piecewise linear approximation of the activation function of the ANN describing the gas solubility in the reactor. This model reformulation yields an MINLP but is beneficial for reducing the computational burden [34].

Besides the mass balances, there is the overall process requirement of a product yield of 2 mol/s C11al, which roughly translates into 6500 metric tons of product per year.

The objective function is defined as

$$\begin{array}{cc}\hfill J& =TA{C}_{\mathrm{dec}}+TA{C}_{\mathrm{col}}+TA{C}_{\mathrm{cat}}+TA{C}_{\mathrm{reac}},\hfill \end{array}$$

(22)

$$\begin{array}{cc}\hfill TA{C}_{\mathrm{dec}}& ={\theta }_1\sum _n\sum _i{\left({\eta }_i{F}_{\mathrm{dec},n,i}^{\mathrm{in}}\right)}^{{\theta }_2},\hfill \end{array}$$

(23)

$$\begin{array}{cc}\hfill TA{C}_{\mathrm{col}}& ={\lambda }_{psol,1}+{\lambda }_{psol,2}·F^{\left(psol\right)}+{\lambda }_{psol,3}·rec+{\lambda }_{psol,4}·{\left(F^{\left(psol\right)}\right)}^2+{\lambda }_{psol,5}·re{c}^2,\hfill \end{array}$$

(24)

$$\begin{array}{cc}\hfill TA{C}_{\mathrm{cat}}& ={\theta }_3{V}_1^{\mathrm{cat}},\hfill \end{array}$$

(25)

$$\begin{array}{cc}\hfill TA{C}_{\mathrm{reac}}& ={\theta }_4{m}_{\mathrm{cat}}+{\theta }_5{V}_{\mathrm{reac}}^{{\theta }_6}+{\theta }_7\left(\sum _i{\theta }_{8,i}{Z}_i-2{\theta }_{8,\mathrm{C}10\mathrm{en}}\right),\hfill \end{array}$$

(26)

where Equation (23) describes the costs for the decanters [35], Equation (24) describes the costs for the extraction solvent column, Equation (25) describes the costs for the catalyst make-up, and Equation (26) describes the costs for the reactor [35], including the reactant stream and make-up streams for the solvents. The parameters $\mathit{\theta }$, $\mathit{\eta }$, and $\mathit{\lambda }$ can be found in Table A2.

The process optimizations are implemented as MINLPs and solved using the GAMS 26.1.0 framework with the multi-start heuristic of the deterministic global optimization software BARON 18.11.12., Cplex 12.8.0 is used as an LP/MIP subsolver, and CONOPT 4.09 is utilized as an NLP subsolver. The calculations are carried out on a Linux PC with 3.40 GHz Intel Core i7-6700 CPU and 16 GB memory.

3. Results

3.1. Process Optimization for Each Candidate Solvent

This section will present the process optimization results for each of the solvent candidates. The results were obtained by an optimization using BARONs multi-start heuristic with 140,000 starting points.

Some important molar flow rates can be found in Figure 2. For each of the three solvent candidates, a simplified version of the flowsheet is depicted. The highlighted boxes containing the molar flow rates point to the respective points in the flowsheet where these values are obtained. Additionally, the degrees of freedom missing from Figure 2 can be found in

Table 2. Optimization results.

Polar Solvent $\mathit{psol}$	Volume ${\mathit{V}}_{\mathbf{reac}}$ (m${}^3$)	Reactor Pressure ${\mathit{P}}_{\mathit{r}}$ (bar)	Temperature ${\mathit{T}}_{\mathit{r}}$ (K)	Recycle Purge $\mathit{\varpi }$ (-)	Extraction Cascade Split Factor $\mathit{\varsigma }$ (-)	Overall Process Costs $\mathit{TAC}$ ($/a)
DMF	273.3775	20	388.15	0.0013	0	1,934,522 (100 %)
DSUC	3458.1775	20	388.15	0.0018	0	6,154,228 (+218 %)
THPO	203.8447	20	388.15	0.0015	0	1,748,021 (−9.64 %)

.

3.1.1. Dimethylformamide (DMF)

The first of the three solvent candidates is the benchmark, DMF. It is well known to be ideally suitable for TMS in conjunction with an alkane. The catalyst recovery of 99.99992% in the extraction cascade is very good, therefore the catalyst concentration in the reactor can be high. A high catalyst concentration leads to faster reaction rates and thus the residence time can be low. The main parameter for the residence time is the reactor volume, $V_{\mathrm{reac}}$, which can be relatively small and with that the investment costs for the reactor are reduced. For the same reason, the reactor pressure, $P_{\mathrm{r}}$, and temperature, $T_{\mathrm{r}}$, are set to their maximum values to achieve fast reaction rates.

The extraction solvent stream in the extraction cascade can be relatively small because DMF has good selectivity for the catalyst. The split factor, $\varsigma$, is set to zero, meaning that all of the stream that could be bypassed directly into the recycle, ${\mathit{F}}_{\mathrm{col}}^{\mathrm{in}},1$, is let into the column. The first reason is that $F_{\mathrm{col}}^{\mathrm{in}},2,\mathrm{DMF}$ is close to the value of the circulated extraction solvent stream $F^{\left(\mathrm{DMF}\right)}$. Therefore, there would be only a little of the polar solvent left at the bottom of the distillation column, leading to a temperature that would be too high, and the distillation would be quite expensive because of the high recovery rate. The second reason is that ${\mathit{F}}_{\mathrm{col}}^{\mathrm{in}},1$ is composed of 92.91% DMF and the combined flow is composed of 92.97% DMF, which is quite close to the purity requirement at the top of the column. Therefore, the recovery rate of the polar solvent at the top of the column becomes smaller and the column becomes cheaper, though its diameter has to be a little larger because of the increased input stream.

A small fraction of the recycled product stream, $\varpi$, needs to be purged. The purge is the only sink for the side product C10an. Therefore, the purge is necessary, although it is quite expensive because of the purged reactant and solvents. The purged solvents need to be reintroduced through the feed stream, together with the amount of catalyst lost in the product stream and the necessary reactant.

The total annualized costs are 1,934,522 $/a.

We assume that the real value for the partition coefficient ${\mathit{P}}_{\mathrm{yx}},\mathrm{DMF}$ lies within the interval $[\widehat{f}(t_1,t_2)+{\underset{¯}{\delta }}_{\mathrm{DMF}},\widehat{f}(t_1,t_2)+{\overline{\delta }}_{\mathrm{DMF}}]$. To account for the surrogate model errors listed in Equation (19), an additional optimization with ${\mathit{P}}_{\mathrm{yx}}=\widehat{\mathit{f}}(t_1,t_2)+\overline{\mathit{\delta }}$ is conducted for each solvent candidate. By adding the maximum error at the sampling points of the surrogate, the catalyst distribution changes and more catalyst will be lost. The worst case annualized costs with regard to the approximation quality of the surrogate are 1,939,897 $/a.

3.1.2. Dimethyl Succinate (DSUC)

The catalyst recovery of 99.9972% in the extraction cascade achieved with DSUC is smaller than in the previous benchmark case. Therefore, the catalyst concentration has to be smaller in the reactor and consequently, the residence time has to be larger. The reactor volume has to be 12.6 times the size of the benchmark case, and the investment costs for the reactor are also higher. The reactor pressure, $P_{\mathrm{r}}$, and temperature, $T_{\mathrm{r}}$, are set to their maximum values to achieve fast reaction rates.

The extraction solvent stream has to be 1.75 times that of the benchmark case since DSUC is not as selective with regard to the catalyst as DMF is. Therefore, the distillation column and the decanters have to be larger, increasing the investment costs. Furthermore, the operating costs for the distillation column are increased because more solvent needs to be vaporized. In contrast to the benchmark case, ${\mathit{F}}_{\mathrm{col}}^{\mathrm{in}},1$ only contains 80.54% DSUC. However, if both streams, ${\mathit{F}}_{\mathrm{col}}^{\mathrm{in}},1$ and ${\mathit{F}}_{\mathrm{col}}^{\mathrm{in}},2$, are summed up the percentage of DSUC remains constant and the overall flux is larger, which reduces the costs due to the reduced recovery rate. Note, that the distillate is now composed of the extraction solvent and the C10en entering the column.

Again, a small fraction of the recycled product stream needs to be purged to remove C10an. In this test case, the overall molar flow rates are higher than in the previous benchmark case. DSUC drags more of the product, C11al, into the recycle. Therefore, to achieve the required amount of C11al in the product flow, the inlet flow of C11al into the extraction cascade needs to be doubled compared to the benchmark case. Furthermore, more solvent is necessary due to the fixed ratio of solvents to the other species in the reactor. Additionally, as there is more of the reactant in the reactor, more of the side products are generated. Consequently, the purge needs to be larger.

The total annualized costs are 6,154,228 $/a, 3.18 times the costs of the benchmark case.

The worst case annualized costs with regard to the approximation quality of the surrogate are 6,177,063 $/a.

3.1.3. Tetrahydropyranone (THPO)

The catalyst recovery in the extraction cascade is 99.99993% and thus even better than in the benchmark case. For that reason, the catalyst concentration in the reactor can be higher and therefore the reactor volume can be smaller. Again, the reactor pressure, $P_{\mathrm{r}}$, and temperature, $T_{\mathrm{r}}$, are set to their maximum values to achieve fast reaction rates.

The extraction solvent stream is a bit smaller than that of the benchmark case. Because the composition of ${\mathit{F}}_{\mathrm{col}},1^{\mathrm{in}}$ is also similar to that of the benchmark case, the same choice of the split factor is used by the optimizer and the entire stream is fed to the column to make the distillation cheaper.

The purge that removes the side-product C10an is a little larger than in the benchmark case.

Interestingly, the total annualized costs for the process with THPO as polar solvent are 1,748,021 $/a, 0.9036 times the costs of the benchmark. This is mainly due to the increased catalyst recovery. The catalyst make-up is the most expensive part of the whole process.

The worst case annualized costs with regard to the approximation quality of the surrogate are 1,758,083 $/a.

4. Conclusions

In this work, an approach for an integrated process and solvent design was presented. In our approach, unsuitable solvent candidates from known databases are eliminated during a screening procedure and a rigorous process optimization is conducted for the remaining solvent candidates.

DSUC and THPO have been identified by the screening procedure as ecologically benign solvent candidates for a hydroformylation process employing a TMS. They have been compared with the well known and developmental toxic benchmark solvent DMF.

The optimization results show that the selectivity of the solvent with respect to the catalyst and the product is of high importance. DSUC has a low selectivity for the catalyst and thus its use leads to an increased catalyst loss. This catalyst loss is compensated by the optimizer with a lower catalyst concentration in the reactor and thus a higher residence time. Furthermore, DSUC dissolves much of the product, leading to a high product recycle and consequently to higher overall fluxes and larger facilities. This renders DSUC unsuitable because of the greatly increased costs.

On the other hand, THPO has shown that it can recover the catalyst even better than DMF under the assumed operating conditions. Though it is slightly harder to separate from the other species by distillation than DMF, the use of THPO reduces the total annualized costs by 9.64%. This leads to the conclusion that the ability of the polar solvent to recover the catalyst while leaving the product in the less polar phase is of utmost importance.

Though the small savings lie within the economic uncertainties of the objective function, we have shown that a process with a well-established toxic solvent can be replaced by a process with an ecologically benign solvent with similar costs. A quantitative evaluation also depends on model accuracy, in particular the accuracy of the partition coefficient of the catalyst in the different solvent systems. In any case, experimental validation is necessary and will be part of future work.

Our present work was limited to already known solvent candidates from databases. To overcome this limitation, further research will be concerned with computer-aided molecular design approaches in Conjunction With Process Optimization.