Genetic Algorithms for Chemical Engineering Optimization Problems

Abstract

Chemical engineering processes are frequently composed of multiple complex phenomena. These systems can be represented by a set of several equations, which are referred to as mathematical model of the process. Optimization in chemical engineering utilizes specialized techniques to determine the values of the decision variables at which the performance of the process, measured as the objective function(s), is minimum or maximum. The profitability of the process improves remarkably as a result of this selection. This benefit has encouraged the broad application of optimization for important industrial challenges. However, many problems in chemical engineering processes are hard to find the optimum using gradient-based algorithms. For example, the cases when the objective functions of the processes are multimodal, discontinuous, or implicit. Genetic algorithms (GAs) are a kind of metaheuristic searching optimization methods, which are inspired by nature, the mechanics of natural evolution and genetics. Genetic algorithms have received significant attention due to their remarkable advantages over classical algorithms. Compared with traditional optimization approaches, GAs are straightforward, robust, capable of handling the non-differentiable, discontinuous, or multimodal problems. The purpose of this paper is to give several case studies using genetic algorithms in chemical engineering optimization problems.

1. Introduction

Chemical engineering processes are frequently composed of multiple complex phenomena. These systems can be represented by a set of several equations, such as z = f(d; x; p), where z is the vector state of the system. The system of equations is referred to as the mathematical model of the process. The performance of the process is predicted by the model from the assigned data of several “input” variables, d and x, and a group of parameters, p. Among the input variables, some, referred as x, can be changed and are known as design variables, while others, referred as d, are predetermined. The performance of the process can be evaluated through a set of output variables, y = g (z), referred as the function of state of the system. Optimization in chemical engineering utilizes specialized techniques to determine the values of the decision variables, x, at which the performance of the process, measured as the objective function(s), I(x), is minimum or maximum. The profitability of the process improves remarkably as a result of this selection of input/operating/decision variables. This benefit has encouraged the broad application of optimization for important industrial challenges. However, many problems in chemical engineering processes are hard to find the optimum using gradient-based algorithms. For example, the cases when the objective functions of the processes are multimodal, discontinuous, or implicit. Genetic algorithms (GAs) are a kind of metaheuristic searching optimization methods, which are inspired by nature, the mechanics of natural evolution, and genetics [1]. Genetic algorithms have received significant attention due to their remarkable advantages over classical algorithms. Compared with traditional optimization approaches, GAs are straightforward, robust, capable of handling the non-differentiable, discontinuous, or multimodal problems. GAs have been effectively employed in a wide range of various engineering, manufacturing, and management applications [2]. The purpose of this chapter is to give several case studies using genetic algorithms in chemical engineering optimization problems. The case studies include the optimization of an autothermal ammonia synthesis reactor, a separation module using membrane technology, and data-driven modeling optimizations of solid oxide fuel cells.

2. Optimization of an autothermal ammonia synthesis reactor

In the chemical process industries, ammonia is one of the most widely manufactured inorganic compounds [3]. The majority of ammonia produced commercially is consumed in fertilizers, with the rest going into plastics, synthetic fibers and resins, pharmaceuticals, explosives, papers, and refrigeration [4]. As a result, modeling and optimization of ammonia synthesis process have received a significant attention from both the academia and industry. Ammonia is produced predominantly from the combination of elements such as nitrogen and hydrogen in a catalytic process using a promoted iron catalyst firstly established by Haber and Bosch as the reaction [4]:

$$N_2+3H_2\rightleftharpoons 2{\mathit{NH}}_3$$E1

The reaction is reversible and exothermic, releasing a significant amount of heat. In order to achieve a high conversion, the heat of the reaction should be removed. Therefore, the process is typically carried out in an autothermal synthesis reactor, in which the heat of reaction is utilized to preheat the feed gas and ensure the suitable temperature inside. The production of ammonia depends on several factors such as the reactor length, the operating pressure, temperature of the feed and reacted gas, the flow rate, and composition of the gas mixture. The optimization problem of the process is to maximize the economic return. Many studies discussing the modeling, simulation, and optimization of an autothermal ammonia synthesis reactor can be found in literature. Some of them can be mentioned here as in Babu et al. [5], Babu and Angira [6], Carvalho et al. [7], Edgar et al. [8], Ksasy et al. [9], Murase et al. [10], Upreti and Deb [11], Yusup et al. [12]. However, the model discussed in the studies of Edgar et al. [8], Murase et al. [10] has some minor errors and has been corrected in Upreti and Deb [11]. Moreover, the studies primarily focus on optimizing reactor length for a specific reactor top temperature, usually 694 K [6, 7, 12], or for a limited set of temperatures [9, 11]. However, as reported in some studies [11, 12], the economic return is determined by the top temperature and also the reactor length (the temperature of feed gas entering to the reaction zone). As a result, rather than a single variable problem of reactor length, the optimization problem should be viewed as a multivariable problem.

In the study [13], both the reactor length and the reactor top temperature are considered in the design variables for maximizing the profit return of the process. In order to solve the multivariate optimization problem, the cyclic coordinate search technique was employed. This method alters the value of one decision variable at a time, and for each coordinate direction, the golden section search was utilized to solve the single variable optimum problem. However, this traditional searching approach is prone to get caught in local optima. Therefore, the genetic algorithm has higher chance to obtain the global optimum profit of the process.

2.1 Problem formulation of ammonia synthesis reactor

The system discussed here is an autothermal synthesis reactor, which is described in [10] and contains the correction of the objective function reported in [6, 11]. The feed gas contains 21.75 mole% nitrogen, 65.25 mole% hydrogen, 5.0 mole% ammonia, 4.0 mole% methane, and 4.0 mole% argon. In an autothermal reactor, the feed gas mixture enters from the bottom of the reactor, flows upward, enters the catalyst zone from the top, and moves downward. In the catalyst zone, the reaction takes place at around 500°C and 200 atm of pressure. The heat generated by the reaction is utilized to preheat the feed gas mixture in counter current flow. Figure 1 shows the schematic diagram of an autothermal ammonia synthesis reactor. The considered factors affecting the synthesis process are the temperature of feed gas at the entrance of the reaction zone (top temperature) and the reactor length. The goal of the optimal design is to determine the conditions that will give the highest economic return from the reactor operation.

Figure 1.

Schematic diagram of an autothermal ammonia reactor [10].

2.1.1 Objective function

The return of the process, which is calculated from the value of the product gas (heating value and ammonia value), subtract the cost of feed gas (as a source of heat only) and minus the amortization of reactor capital expenses, is the objective function for maximization (F). Other operating costs are not considered [11].

$$\begin{array}{l}F=1.3356\times {10}^7-1.708\times {10}^4{N}_{N_2}+704.09\left(T_g-T_0\right)-699.27\left(T_f-T_0\right)\\ -{\left(3.4566\times {10}^7+1.9837\times {10}^{9}L\right)}^{\frac{1}{2}}\end{array}$$E2

in which,

T_f is the temperature of the feed gas.

T_g is the temperature of the reacting gas (the gas in the catalyst zone).

T₀ is the top temperature or temperature at the inlet of the catalyst zone.

$N_{N_2}$ is the flow rate of nitrogen.

L is the length of the reactor.

2.1.2 Equality constraints

The heat balance for the feed gas and the reacting gas and the mass balance for the nitrogen flow along the catalyst zone, respectively, give the mathematical model for the system:

$$\frac{{\mathit{dT}}_f}{\mathit{dx}}=-\frac{{\mathit{US}}_1}{{\mathit{WC}}_{\mathit{pf}}}\left(T_g-T_f\right)$$E3

$$\frac{{\mathit{dT}}_g}{\mathit{dx}}=-\frac{{\mathit{US}}_1}{{\mathit{WC}}_{\mathit{pg}}}\left(T_g-T_f\right)+\frac{\left(-\mathrm{\Delta }H\right)S_2}{{\mathit{WC}}_{\mathit{pg}}}\left(-\frac{{\mathit{dN}}_2}{\mathit{dx}}\right)$$E4

$$\frac{{\mathit{dN}}_{N_2}}{\mathit{dx}}=-f\left(k_1\frac{p_{N_2}{p}_{H_2}^{1.5}}{p_{{\mathit{NH}}_3}}-k_2\frac{p_{{\mathit{NH}}_3}}{p_{H_2}^{1.5}}\right)$$E5

in which

$$k_1=1.78954\times {10}^{4}exp\left(\frac{-20800}{{\mathit{RT}}_g}\right)$$E6

$$k_2=2.5714\times {10}^{16}exp\left(\frac{-47400}{{\mathit{RT}}_g}\right)$$E7

$$p_{H_2}=3p_{N_2}$$E8

$$p_{N_2}=\frac{286N_{N_2}}{2.598N_{N_2}^0+2N_N}$$E9

$$p_{{\mathit{NH}}_3}=\frac{286\left(2.23N_N^0-2N_{N_2}\right)}{2.598N_N^0+2N_N}$$E10

The differential equations are valid in the interval [0, L], in which L is the reactor length (or the length of the catalyst zone).

The notations $T_f^0$, $T_g^0$, and $N_{N_2}^0$ denote the initial value at x = 0 (at the inlet of the catalyst zone) for T_f, T_g, and $N_{N_2}$, respectively, and are given by

$$T_f^0=T_g^0=T_0,N_{N_2}^0=701.2\text{kmol}/{\mathrm{hm}}^2$$E11

Other notations of the system are summarized in Table 1.

Notation
Cpf	Heat capacity of the feed gas
Cpg	Heat capacity of the reacting gas
f	Catalyst activity
ΔH	Heat of reaction
N	Mass flow of component designed by subscript
k	Reaction rate constant
p	Partial pressure of component designated by subscript
R	Universal gas constant
S1	Surface area of catalyst tubes per unit length of reactor
S2	Cross-sectional area of catalyst zone
U	Overall heat transfer coefficient
W	Total mass transfer flow rate

Table 1.

Notation of the synthesis system.

2.1.3 Box constraints

The variables are subjected to the following physical constraints, as is typical in industries [10]:

$$0<L\le 10,400\le T_f\le 800,0\le N_{N_2}\le 3220$$E12

The length of the reactor and the top temperature are chosen as the design variables. The remaining variables (T_f, T_g, and $N_{N_2}$) can be calculated from the model by three differential Eqs. (3), (4), and (5). Then, the objective function is determined by Eq. (2). Due to the constraints of temperatures, the variable T₀ is also set to be within 400 and 800.

The optimal design problem is summarized as follows:

$$\begin{array}{l}maximizeF=F\left(x,T_0\right)\\ s.t.\left\{\begin{array}{l}\frac{{\mathit{dT}}_f}{\mathit{dx}}=-\frac{{\mathit{US}}_1}{{\mathit{WC}}_{\mathit{pf}}}\left(T_g-T_f\right)\\ \frac{{\mathit{dT}}_f}{\mathit{dx}}=-\frac{{\mathit{US}}_1}{{\mathit{WC}}_{\mathit{pg}}}\left(T_g-T_f\right)+\frac{\left(-\mathrm{\Delta }H\right)S_2}{{\mathit{WC}}_{\mathit{pg}}}\left(-\frac{{\mathit{dN}}_2}{\mathit{dx}}\right)\\ \frac{{\mathit{dN}}_2}{\mathit{dx}}=-f\left(k_1\frac{p_{N_2}{p}_{H_2}^{1.5}}{p_{{\mathit{NH}}_3}}-k_2\frac{p_{{\mathit{NH}}_3}}{p_{H_2}^{1.5}}\right)\\ T_f^0=T_g^0=T_0,N_{N_2}^0=701.2\\ 400\le T_f\le 800,0\le N_{N_2}\le 3220\\ 0<L\le 10,400\le T_0\le 800\end{array}\right.\end{array}$$E13

2.2 Optimization strategy

The system of ordinary differential Eqs. (3), (4), and (5) with initial conditions (11) was solved by Runge–Kutta fourth-order method. The system is well defined when the top temperature (T₀) and the interval of integration (the reactor length L) are assigned. After that, the economic return is clearly evaluated, and it can be considered as two-variable function of the top temperature and reactor length. The genetic algorithm is employed to find the optimal solution. In order to handle the constraints, a penalty or barrier is defined to the objective function whenever any of variable limits is violated. The proposed strategy is detailed as follows.

2.2.1 Genetic algorithm for optimization problem

The range of the design variables is $0<x\le 10,400\le T_0\le 800$. The parameters of GA such as population size, crossover probability, mutation probability values were set to be 50, 1.0, and 0.30, respectively. The selection was chosen as roulette wheel selection with elitism. The number of generations was set to be 500.

2.2.2 Barrier method for constrained optimization

In barrier or penalty methods, the objective function will receive an undesired value when one of the constraints is violated. Therefore, the solution will be kept in the feasible region. The objective function has been modified as

$$F=\left\{\begin{array}{l}F\text{if}400\le T_f\le 800,0\le N_{N_2}\le 3220\\ 0\text{otherwise}\end{array}\right.$$E14

2.2.3 Parameters

The parameters were obtained from the literature [7, 10] and summarized in Table 2.

Parameter	Value	Unit
Cpf	0.707	kcal/kg K
Cpg	0.719	kcal/kg K
f	1.0
ΔH	−26,000	kcal/kmol
R	1.987	kcal/kmol K
S1	10	M
S2	0.78	m^2
U	500	kcal/h m^2 K
W	26,400	kg/h

Table 2.

Model parameters [7].

2.3 Optimization results

Figure 2 shows the fitness values as a function of generation. As can be observed, the fitness function value achieved the highest after roughly 20 generations and then stayed unchanged. After 100 generations, it was obtained that the reactor length should be 6.772 m, and the top temperature should be 707.09 K. The process produces a profit of 5.018× 10⁶ $ per year. The other parameters of the process are summarized in Table 3 and compared with the findings of a cyclic coordinate search [13]. The profit value is slightly higher than those reported in the literature, which focused solely on reactor length optimization. From the results, the temperature at the entrance of the catalyst zone should be slightly higher, and that the reactor length should also be slightly longer than previously reported.

Figure 2.

Fitness value versus generations.

Variables	Interval	Cyclic coordinate [13]	Genetic algorithm
x (m)	[0,10]	6.724	6.772
T0 (K)	[600,800]	700.27	707.09
Tf (K)	[400,800]	400.00	401.09
Tg (K)		629.94	631.12
NN2	[0,3220]	490.68	490.68
F (10^6$/year)		5.018	5.018

Table 3.

Maximization results.

3. Economic optimization of membrane module for ultrafiltration of protein solution

The behavior of permeate flux has a significant impact on the performance of cross-flow ultrafiltration. Many factors cause flux declination, such as solution properties, membrane properties, and operation conditions. The majority of current research has centered on increasing membrane performance in terms of permeability and selectivity [14, 15]. Just a few studies have paid attention to the configuration and operation of the membrane module [16].

Various factors determine the decision of membrane module geometry for a given application, including fabrication method, power consumption, and fouling potential [17]. Manufacturers frequently recommend the membrane module design from the fabrication standpoint [17]. There is virtually no evidence that their approach prioritizes the energy efficiency. Currently, with a growing in energy concern and a falling in membrane cost, the membrane module design should place a higher attention on energy efficiency. As a result, it is necessary to propose a module design methodology that takes into account the energy factor.

Furthermore, membrane operating conditions are usually decided by user experience, a handbook, or a manual from the membrane supplier. However, the permeate flux equation governing the performance of the membrane system varies greatly between different situations. In this aspect, for any specific application, a general methodology for the design and operation conditions should be studied.

In cross-flow ultrafiltration of protein solution, Nguyen et al. [18] proposed a simple combined model, which simultaneously considers pore blockage and cake filtration, to describe the flux declination. Then, in the study [19], the correlation between the steady-state permeate flux and operation parameters was reported. From the steady-state operation equation, optimal design and operation conditions for each particular application could be established.

However, just a few reports on the optimization of membrane processes and cost estimation have been published, or the cost estimation is too general. For example, Wiley et al. [17] optimized the membrane module configurations for brackish water desalination. However, the operation mode is single-pass and only the membrane cost and energy cost were taken into account. Sethi and Wiesner [20] developed the cost model for the removal of natural organic matter, but the study has not conducted the optimization. In membrane technology, the feed and bleed operation mode, which combines the batch and the single-pass configurations, is commonly utilized for continuous full-scale filtration [21, 22]. Therefore, the optimization of a membrane module operated in feed and bleed mode for protein ultrafiltration is considered. The membrane geometry dimensions and operating conditions are design variables in the problem. The system is represented by a set of ordinary differential equations. The objective function is the annual cost, which consists of various types of capital investments and an operating expense. The capital investments are classified into several categories, which are individually correlated to plant scale, particularly the membrane area. The operating expense is the power consumption.

3.1 Process configuration and model calculation

3.1.1 Membrane plant configuration

The configuration of filtration system is continuous feed and bleed, which is shown schematically in Figure 3. The notations are summarized in Table 4. There are two main pumps in this operation: the feed pump provides the necessary trans-membrane pressure, while the recirculation pump maintains the cross-flow rate through the modules. The concentrate is continually withdrawn from the system at a flow rate (R).

Figure 3.

Schematic configuration of feed-and-bleed mode membrane system.

Notation	Name and units
F	Feed flow rate [m^3/hr]
R	Retentate (concentrate) flow rate [m^3/hr]
Q	Recirculation flow rate [m^3/hr]
Flow	Flow rate in membrane module [m^3/hr]
P	Permeation flow rate [m^3/hr]
PF	Pressure at outlet of feed pump [kPa]
Pi	Pressure at the inlet of membrane module [kPa]
Po	Pressure at the inlet of membrane module [kPa]
EP	Energy consumed by the feed pump [kW]
EQ	Energy consumed by the recirculation pump [kW]
ϕ0	Initial concentration of protein solution [m^3/m^3]
ϕi	Inlet concentration of protein solution [m^3/m^3]
ϕf	Final concentration of protein solution [m^3/m^3]
ϕP	Concentration of protein in permeate flux [m^3/m^3]
u	Fluid flow velocity [m/s]
ρ	Fluid density [kg/m^3]
μ	Fluid viscosity [kg/(m·s)]
w, h, L, dh	Width, height, length, hydraulic diameter of the membrane module

Table 4.

Summary of system configuration notations.

3.1.2 Modeling of membrane modules

The material balance for total mass and protein give:

$$\begin{array}{l}F=R+P\\ F{\varphi }_0=R{\varphi }_{\mathrm{f}}+P{\varphi }_{\mathrm{p}}\\ F{\varphi }_0+Q{\varphi }_{\mathrm{f}}=\left(F+Q\right){\varphi }_{\mathrm{i}}\end{array}$$E15

The viscosity and density of protein solution correlate to its concentration [23]:

$$\begin{array}{l}{\mu }_{\mathrm{m}}=8.94\times {10}^{-4}exp\left(13.5482{\varphi }_{\mathrm{m}}\right)\\ {\rho }_{\mathrm{m}}=1000\left(1-{\varphi }_{\mathrm{m}}\right)+1360{\varphi }_{\mathrm{m}}\end{array}$$E16

in which, μ_m, ϕ_m are the viscosity and concentration of the mixture (solution), respectively. ρ_m is the density of protein solution, 1000 kg/m³ is the density of water, and 1360 kg/m³ is the density of dry protein powder.

The permeate flux through the membrane is [19].

$$\frac{J}{u}=3.66\times {10}^{-7}{\left(\frac{P}{{\rho }_{\mathrm{m}}{u}^2}\right)}^{0.27}{\left(\frac{{\rho }_{\mathrm{m}}{\mathit{ud}}_{\mathrm{h}}}{{\mu }_{\mathrm{m}}}\right)}^{0.52}$$E17

in which P is the trans-membrane pressure, u is the cross-flow velocity, d_h is the hydraulic diameter of the flow channel.

The equation for permeate flux can be rewritten as:

$$J=5.124\times {10}^{-9}{\left(\frac{P}{{\rho }_{\mathrm{m}}{u}^2}\right)}^{0.27}{\left(\frac{{\rho }_{\mathrm{m}}{\mathit{ud}}_{\mathrm{h}}}{{\mu }_{\mathrm{m}}}\right)}^{0.52}\left(\frac{u}{d_{\mathrm{h}}}\right)$$E18

or in terms of shear rate $\dot{\gamma }=\frac{6u}{h}=\frac{12u}{d_{\mathrm{h}}}$, where h is the channel height:

$$J=4.27\times {10}^{-10}{\left(\frac{P}{{\rho }_{\mathrm{m}}{u}^2}\right)}^{0.27}{\left(\frac{{\rho }_{\mathrm{m}}{\mathit{ud}}_{\mathrm{h}}}{{\mu }_{\mathrm{m}}}\right)}^{0.52}\dot{\gamma }$$E19

The flow rate/velocity drop and channel length change are calculated from the total mass balance and component balance within the control volume dz across the channel. Protein is assumed to be fully rejected by the membrane:

$$d\left(\mathit{\text{Flow}}\times {\varphi }_{\mathrm{m}}\right)=0$$E20

$$d\left(\mathit{\text{Flow}}\right)=-J\times d\left(\mathit{wz}\right)=-Jw\mathit{dz}$$E21

The pressure loss is estimated by the Darcy-Weisbach Equation [24, 25].

$$\mathit{dP}=-4f{\rho }_{\mathrm{m}}\frac{\mathit{dz}}{d_{\mathrm{h}}}\frac{u^2}{2}$$E22

in which f is the friction factor.

The set of ordinary equations that describes the membrane module system was established as follows [26].

$$\left\{\begin{array}{l}\frac{d\left(\mathit{\text{Flow}}\right)}{d{\varphi }_{\mathrm{m}}}=-\frac{\mathit{\text{Flow}}}{{\varphi }_{\mathrm{m}}}\\ \frac{\mathit{dz}}{d{\varphi }_{\mathrm{m}}}=\frac{\mathit{\text{Flow}}}{{\varphi }_{\mathrm{m}}wJ}\\ \frac{\mathit{dP}}{d{\varphi }_{\mathrm{m}}}=-4f{\rho }_{\mathrm{m}}\frac{1}{d_{\mathrm{h}}}\frac{u^2}{2}\frac{\mathit{\text{Flow}}}{{\varphi }_{\mathrm{m}}}\frac{1}{wJ}\\ {\mu }_{\mathrm{m}}=8.94\times {10}^{-4}exp\left(13.5482{\varphi }_{\mathrm{m}}\right)\\ {\rho }_{\mathrm{m}}=1000\left(1-{\varphi }_{\mathrm{m}}\right)+1360{\varphi }_{\mathrm{m}}\\ J=5.124\times {10}^{-9}{\left(\frac{P}{{\rho }_{\mathrm{m}}{u}^2}\right)}^{0.27}{\left(\frac{{\rho }_{\mathrm{m}}{\mathit{ud}}_{\mathrm{h}}}{{\mu }_{\mathrm{m}}}\right)}^{0.52}\left(\frac{u}{d_{\mathrm{h}}}\right)\\ u=\frac{\mathit{\text{Flow}}}{h\times w}\\ f=\frac{24}{Re}\text{if}Re<2000\\ f=\frac{0.079}{{Re}^{0.25}}\text{if}Re>2000\end{array}\right.$$E23

in the range of concentration [ϕ_i, ϕ_f] and the initial condition

$$\left\{\begin{array}{l}\mathit{\text{Flow}}\left({\varphi }_{\mathrm{i}}\right)=F+Q\\ z\left({\varphi }_{\mathrm{i}}\right)=0\\ P\left({\varphi }_{\mathrm{i}}\right)=P_{\mathrm{i}}\end{array}\right.$$E24

The system of the ordinary equations can be solved numerically by Runge–Kutta fourth-order method [27] to obtain the flow rate, the length, and the pressure. From that, the two important factors determining the total cost, membrane area, and total energy were calculated:

$$A_{\text{membrane}}=w\times L$$E25

$$\begin{array}{l}{E}_{\text{pumps}}=E_{\mathrm{P}}+E_{\mathrm{Q}}\\ =F\times P_i+Q\times \mathrm{\Delta }{P}_{\text{drop}}\end{array}$$E26

In this equation,

E_P, E_Q are the power supplied by the feed pump and recirculation pump, respectively.

P_F, P_i, P_o are the pressure at the outlet of the feed pump, at the inlet and outlet of the membrane module, respectively.

ΔP_drop is the pressure drop in the membrane module.

3.1.3 Cost estimation

3.1.3.1 Operating cost

The operating cost consists of power consumption of the pumps and membrane replacement. The annual energy expense of the pumps is calculated as

$$C_{\text{energy}}=\frac{E_{\text{pumps}}}{\eta }\times 8000\times \frac{3600}{1000}\times \text{electricity price}\left[\text{\$}/\text{year}\right]$$E27

E_pumps is from Eq. (26). The system is assumed to work 8000 hr./year (24 hr./day and 333 days/year). η is the efficiency of the pumps, which is set to be 0.7. The electricity price is supposed to be 0.08 $/kWh, which is the price for the industrial sector in the United States [28].

The membrane replacement cost is calculated as

$$C_{\text{membrane}}=A_{\text{membrane}}\times c_{\text{membrane}}\times \left(\frac{A}{P}\right)\left[\text{\$}/\text{year}\right]$$E28

where C_membrane [$/m²/year] is the membrane replacement cost calculated per year.

c_membrane [$] is the membrane price per unit area,

$\left(\frac{A}{P}\right)$ is the amortization factor, which presents the time value of money [29] and is calculated as a function of interest rate i and the membrane life.

The membrane price is usually about 200 $/m² [9], and membrane life is 12–18 months. Therefore, the membrane replacement cost per year is roughly estimated as 200 $/m²/year for the interest of i = 8%.

3.1.3.2 Capital cost

It is widely observed that capital costs are correlated to the size in the power-law form [30]:

$$\text{cost}=k{\left(\text{size}\right)}^n$$E29

In order to achieve higher accuracy, rather than simply predicting the whole capital cost of the membrane plant to capacity, Sethi and Wiesner [20] divided the capital investment into several major categories, which was correlated to the size independently. The major categories include pumps and other manufactured equipment.

1. Pump capital cost

   The pumps capital cost can be estimated as (Perry et al. [31]):

   $${C^{\ast }}_{\text{pump}}=I\times f_1\times f_2\times C_L\times 81.27\times {\left(Q\times P\right)}^{0.4}$$E30

   in which.

   I: a cost index ratio for updating the cost to the recent year.

   f₁: an adjust factor for pump construction material.

   f₂: an adjust factor for suction pressure range.

   C_L: a factor used to incorporate labor costs.

   Q: flow capacity of the pump [m³/h].

   P: pressure outlet of the pump [kPa].

   The cost index, I in Eq. (30), can be referred to as the chemical engineering (CE) index to update the cost. It can be obtained from [32], and the value of 2.4 is used. C_L = 1.4 with the assumption that 40% of the cost is required to install the equipment [33]. The factors f₁ and f₂ can be found in [31]. f₁ = 1.5 when the material is stainless steel. f₂ = 1.0 when the pump pressure is below 10 bar (1 MPa).

   The pump size (Q × P in Eq. (30)) of the two pumps (feed and recirculation pump) is

   $$\text{pump size}=\left(F+Q\right)\times P_{\mathrm{i}}$$E31

2. Capital cost of other equipment

   In membrane application, the membrane area is the key parameter, which determines the plant capacity [34]. Thus, the membrane area is chosen as the basic for the estimation of various components in the capital costs.

   Non-membrane equipment and facilities, excluding the pumps, were grouped into four main categories: (1) pipes and valves; (2) instruments and controls; (3) tanks and frames; and (4) miscellaneous. The capital cost of each is correlated to the membrane area as follows (Sethi and Wiesner [20])

   1. Pipes and valves

      $${C^{\ast }}_{\mathrm{PV}}=6000{\left(A_{\text{membrane}}\right)}^{0.42}\left[\text{\$}\right]$$E32

   2. Instruments and controls

      $${C^{\ast }}_{\mathrm{IC}}=1500{\left(A_{\text{membrane}}\right)}^{0.66}\left[\text{\$}\right]$$E33

   3. Tanks and frames

      $${C^{\ast }}_{\mathrm{TF}}=3100{\left(A_{\text{membrane}}\right)}^{0.53}\left[\text{\$}\right]$$E34

   4. Miscellaneous

      $${C^{\ast }}_{\mathrm{MI}}=8000{\left(A_{\text{membrane}}\right)}^{0.57}\left[\text{\$}\right]$$E35

3. Annual capital cost

   The capital cost can be annualized using the amortization factor as

   $$C_{\text{capital cost}}={C^{\ast }}_{\text{capital}}\times \left(\frac{A}{P}\right)\left[\text{\$}/\text{year}\right]$$E36

   For the plant design year of 20 years and the interest rate 8%, the amortization factor will be about 0.1.

3.2 Formulizations of the problem

3.2.1 Fix parameters and design variables

In the problem, some variables, called input variables, are fixed due to the requirement of the design. In membrane design, these are feed flow F, inlet concentration ϕ₀, and outlet concentration ϕ_f. The protein is assumed to be entirely rejected by the membrane (ϕ_p = 0).

The design variables were: channel geometry (width × length × height), the inlet pressure (P_i), and recirculation flow rate (Q).

3.2.2 Objective function

The objective function is the sum of capital cost and operating cost, which were annualized:

$$\begin{array}{l}\text{minimum}f\left(\mathbf{x}\right)=\text{annual total cost}=\text{annual capital cost}+\text{annual operating cost}\\ =C_{\mathit{\text{energy}}}+C_{\mathit{\text{membrane}}}+C_{\mathit{\text{capital}}\left(\mathit{\text{pumps}}+\mathit{\text{other}}\right)}\end{array}$$E37

3.2.3 Constraints

The pressure at the outlet point should be positive. This constraint is satisfied by assigning a high value to the objective function if the outlet pressure is negative.

The decision variables are frequently limited on a finite range

$${\mathbf{x}}_l\le \mathbf{x}\le {\mathbf{x}}_u$$E38

in which x is the vector of decision variables x = (P_i,Q,w,h), subscripts l and u indicate the lower and upper bound.

The system parameters and variables are summarized in Table 5.

Parameters	Value
Feed flow rate (m^3/hr)	0.02–200
Inlet pressure (kPa)	200–1000
Recirculation flow rate (m^3/hr)	0–50
Initial solid fraction (m^3/m^3)	0.1
Final solid fraction (m^3/m^3)	0.4
Plant design year (year)	20
Interest rate (%)	8
Energy price ($/kWh)	0.08
Efficiency of pumps (%)	70
Operating temperature (°C)	25
Module height (mm)	0–100
Module width (m)	0–30

Table 5.

System parameters and variables.

3.2.4 Optimization by genetic algorithm

The parameters of GA such as population size, crossover probability, mutation probability values were set to be, 100, 1.0, and 0.30, respectively. The selection was based on roulette wheel with elitism, which means the most fit individual is guaranteed a place in the next generation. The number of generations was assigned to be 500. Because the problem is to minimize the cost, the fitness function was defined as:

$$\text{fitness}=\left\{\begin{array}{l}0\text{if}f\left(x\right)>5\times {10}^5\\ 5\times {10}^5-f\left(x\right)\text{otherwise}\end{array}\right.$$E39

3.3 Optimum design

For the demonstration of this method, optimum designs of several feed flow rates have been carried out. The lower limit of the membrane width is 0.1 m, the lower limit for the module height is 0.5 mm. The designs are shown in Table 6.

Feed [m^3/hr]	ϕ0 [−]	ϕf [−]	Pressure [kPa]	Recirculation [m^3/hr]	width [m]	height [mm]	total cost [$/yr]
0.02	0.1	0.4	523	2.8	0.1	5.0	1.29 × 10^3
0.2	0.1	0.4	1000	4.9	0.1	8.9	4.30 × 10^3
2	0.1	0.4	1000	0.2	0.1	6.9	1.18 × 10^4
20	0.1	0.4	987	0.2	1.3	5.0	5.60 × 10^4
200	0.1	0.4	1000	0.8	11.1	5.0	3.65 × 10^5

Table 6.

Optimum designs of membrane module.

Figure 4 presents the optimum total cost per unit of feed flow. The cost per unit of feed flow decreases with an increase in plant capacity. It reflects the economies of scale.

Figure 4.

The behavior of cost per unit flow rate design in optimum condition with plant capacity.

The results also suggest that the membrane module dimensions and operation condition will change greatly depending on the process requirements, such as the required feed capacity. It is challenging to predict the direction. It might be concluded that the permeate flux also greatly affects the geometric design and operation strategy in membrane separation processes. It is difficult to find a general rule for the design, for each specific system, the correlation between the permeate flux and operating conditions and membrane geometry should be investigated.

4. Modeling and optimization of the BSCF-based single-chamber solid oxide fuel cell by artificial neural network and genetic algorithm

Fuel cells that are highly effective and green technology for converting chemical energy stored in fuel to useable power are currently regarded as one of the most promising approaches for future energy requirements [35]. The solid oxide fuel cell (SOFC) has demonstrated an exceptional integration of advantages, such as high efficiency, fuel flexibility, wide contamination acceptance, and low pollution [36, 37]. Modeling and simulation are valuable tools for determining the impact of various design factors and operating conditions on cell performance, as well as for improving fuel cells [38, 39, 40]. Plenty of models have been reported to add to the understanding of fuel cells. Modeling approaches can be categorized into two types: theoretical and empirical one [39, 41, 42]. In the theoretical approach, the spatial dimensions of the models range from simple 0 (0-D) [43, 44] and 1 (1-D) [42, 45, 46, 47], to more complicated 2 (2-D) [48, 49, 50, 51] and 3 (3-D) [52, 53, 54], all with various characteristics and directed at different objectives. The mathematical models, which are based on conservation principles, require a lot of data on parameters and properties of fuel cell, as well as complicated equations and time-consuming calculation.

Empirical or data-driven approach may be more feasible for fuel cell users since the behavior can be quickly and simply deduced without a comprehensive understanding of the internal components, just based on the experimental data [39, 41]. Least squares support vector machine (LS-SVM) [55], Hammerstein model [56, 57] are examples of these approaches. In this approach, artificial neural network (ANN) shows several advantages, including high nonlinearity, rapid computation, a low degree of error in matching experimental data. Using ANNs to model SOFCs appears to be a very promising method.

In this section, an ANN was used to model the performance of the BSCF/GDC-based cathode SOFC. The cell voltage was predicted from cathode sintering temperature, cell operating temperature, and cell current. Several network architectures were examined to find the best structure, and the network was trained using back-propagation methods. The data for training, validation, and testing were taken from our study [58]. The genetic algorithm and the developed ANN were then used to find the best conditions for achieving maximum power.

4.1 Artificial neural network models

Artificial neural networks (ANNs), which were analogous to biological nervous systems, consist of interconnected nodes known as neurons to receive and transfer data [59]. The most basic form, feed-forward architecture, is made up of an input layer, one hidden layers, and an output layer. The input and output layers have the same number of neurons as the number of inputs and outputs in the system to be modeled. Weighted connections connect each neuron to every other neuron in the next layer. In any layer except the input, the weighted sum of data from the previous layer is the input of a neuron. The neuron then activates the data using a function and transfers the response to all neurons in the next layer. The size of the hidden layers is a significant factor that affects the estimation precision because it can make the network become insufficient or overfitting [60]. The number of neurons in hidden layer is generally determined through trials. Figure 5 illustrates a 3–5-1 feed forward artificial neural network with operating temperature, sintering temperature, and current as inputs.

Figure 5.

Artificial neural network (3–5-1) structure.

The activation function employed in this model is the logistic sigmoid$f\left(x\right)=1/\left(1+e^{-x}\right)$.

The input data (x_i) are scaled to normalized value (x_norm) to enhance the performance of the network:

$$x_{\mathrm{norm}}=0.8\frac{\left(x-x_{\min }\right)}{x_{\max }-x_{\min }}+0.1$$E40

in which x_max and x_min are the bounded interval of the experimental data.

To assess the performance of ANN, the mean squared error (MSE) and coefficient of determination (R²) are usually used [61].

Various factors affect the performance of fuel cells such as cathode and anode structure, electrolyte material and thickness, cell temperature, inlet and outlet gas compositions. Two important factors, cathode sintered temperature and cell operating temperature, were considered in this model. The sintered temperature is from 1000–1050°C, whereas the operating temperature ranges from 625–700°C. The sintered temperature affects the structure of the obtained cathode as reported in [58]. The explanation of the range for the investigated parameters can be found in [62].

An ANN with one input layer, one hidden layer, and one single output layer was proposed. Current density, sintered temperature of the cathode, and cell operating temperature are the inputs. Back-propagation algorithm [63] was used to train the network. The maximum number of iteration and minimum performance gradient were set to 400 and 10⁻⁵, respectively, to stop the training. The proper network structure is determined through a series of trial tests. The data were split into three subsets at random: training, validation, and test, each containing 70, 20, 10% of the total samples, respectively. The validation and test sets are necessary for evaluating the validation and power of the networks.

The parameters of the neural network were saved and utilized in the next stage to optimize the power density using genetic algorithms.

4.2 Optimization by genetic algorithm

The objective function is the power density of the fuel cell

$$P=I\times V$$E41

where P is the power density (mW.cm⁻²), I is the current density (mA.cm⁻²) and is the input to the ANN model, and V is the cell voltage (V), which is calculated from the model.

The design variables and their corresponding ranges are summarized as follows:

• sintered temperature of the cathode, [1000–1050] (°C).

• operating temperature of the cell, [625–700] (°C).

• electric current of the cell, [0–1500] (mA.cm⁻²)

The parameters of GA as population size, mutation probability values were set to be 100 and 0.10, respectively. The survival of the individuals was decided by roulette wheel with elitism. The number of generations was 500.

4.3 Optimization results

Figure 6 depicts the fitness values (maximum and mean) of the population versus generation. As indicated in the figure, after about 20 generations, the value of fitness function attained to a maximum value and then remained unchanged. After 100 generations, the maximum fuel cell power density of 451.64 mW/cm² could be achieved at the sintered temperature of 1005°C, operating temperature of 668°C, and current density of 777 mA/cm².

Figure 6.

The fitness values versus generation.

5. Conclusions

The application of genetic algorithm in chemical engineering processes has been illustrated by three case studies. The results suggest that the optimum conditions of complex chemical problems can be easily obtained using genetic algorithm. The successes of genetic algorithm for the challenging problems reported herein, the development of many faster and flexible versions of GA, the improvement of computing ability all suggest the continually increasing impact of metaheuristic methods in chemical engineering systems.