Multi-Phase Batch Reactor

Background

Many batch processes use different catalysts or operating regimes for successive reactions. Here, feed \(A\) is first converted into intermediate \(B\); after a catalyst change, \(B\) is converted into product \(P\). Switching too early leaves too little intermediate, while switching too late leaves too little time to form product. Higher temperature intensity accelerates both reactions but carries an operating cost.

This makes the process a genuine multi-phase problem: the dynamics change at an optimized event time, while all concentrations remain continuous.

Problem formulation

The normalized concentration vector in each phase is \(x_j=(A_j,B_j,P_j)\). During phase 1, \(0\le t\le t_s\),

\[ \begin{aligned} k_1(u_1)&=0.20+1.80u_1, \\ \dot A_1&=-k_1A_1, \\ \dot B_1&= k_1A_1, \\ \dot P_1&=0. \end{aligned} \]

During phase 2, \(t_s\le t\le2\ \mathrm h\),

\[ \begin{aligned} k_2(u_2)&=0.10+2.40u_2, \\ \dot A_2&=0, \\ \dot B_2&=-k_2B_2, \\ \dot P_2&=\eta k_2B_2, \end{aligned} \qquad \eta=0.92. \]

The initial condition, event connection, and bounds are

\[ \begin{aligned} x_1(0)&=(1,0,0), \\ x_1(t_s)&=x_2(t_s), \\ 0.20&\le t_s\le1.80, \\ 0&\le u_1(t),u_2(t)\le1. \end{aligned} \]

The terminal concentrations are free. The objective balances product yield against temperature effort:

\[ \min J=-P_2(2)+0.15\left( \int_0^{t_s}u_1(t)^2\,\mathrm dt+ \int_{t_s}^{2}u_2(t)^2\,\mathrm dt \right). \]

Variables and units

Symbol	Meaning	Unit
\(t\)	Batch time	h
\(A,B,P\)	Normalized feed, intermediate, and product concentrations	-
\(u_1,u_2\)	Normalized temperature intensities	-
\(k_1,k_2\)	Effective first-order rate constants	h\(^{-1}\)
\(t_s\)	Catalyst-switch time	h
\(\eta\)	Product yield from consumed \(B\)	-

Modeling choices

Each catalyst regime is represented by its own Lobatto phase with 24 mesh intervals and two interpolation points per interval. This keeps the concentration and temperature reconstructions piecewise linear and prevents between-node bound overshoot without increasing the number of decision points relative to the original higher-order mesh. Shared static parameters hold the three concentrations on both sides of the event, so using the same parameters in the terminal boundary of phase 1 and the initial boundary of phase 2 enforces continuity exactly. The switch time is another shared parameter and is optimized together with both control profiles.

The yield \(\eta<1\) represents material lost to unmodeled by-products. The effort weight produces an interior optimum in both phases instead of the uninformative solution that holds temperature at its upper bound throughout.

Run the example

python -m examples.multiphase_batch_reactor

Save the figure without opening a window:

python -m examples.multiphase_batch_reactor --save batch-reactor.png --no-show

Key implementation

system = System(["a_switch", "b_switch", "p_switch",
                 "a_final", "b_final", "p_final", "t_switch"])
a_s, b_s, p_s, a_f, b_f, p_f, t_s = system.s

phase_1.set_boundary_condition(
    [1.0, 0.0, 0.0], [a_s, b_s, p_s], 0.0, t_s
)
phase_2.set_boundary_condition(
    [a_s, b_s, p_s], [a_f, b_f, p_f], t_s, BATCH_TIME
)
system.set_phase([phase_1, phase_2])
system.set_objective(
    -p_f + EFFORT_WEIGHT * (phase_1.I[0] + phase_2.I[0])
)

Verified result

Ipopt terminates successfully with objective \(-0.45376749\). The optimized catalyst switch occurs at \(t_s=1.092691\ \mathrm h\), and the final product fraction is \(0.658482\). The time-average temperature intensities, computed as \((t_f-t_0)^{-1}\int u\,\mathrm dt\), are \(0.796591\) in phase 1 and \(0.860164\) in phase 2. The maximum concentration mismatch across the event is below \(2\times10^{-6}\), and a 2,001-point check in each phase shows no path-bound violation at the reported precision.

Source code

See the complete runnable example: examples/multiphase_batch_reactor.py.