Documentation Index
Fetch the complete documentation index at: https://feasible-1447f9c5.mintlify.app/llms.txt
Use this file to discover all available pages before exploring further.
In chemical reactions, an equilibrium is reached when the reaction rate forming a compound is equal to the one reversing it. The concentrations of the reactants and the product are then stable over time.
For a process at constant pressure and temperature, this means that the Gibbs free energy is minimized. It can be calculated as
G(T,p,x)=i∑xi[ci+log∑xixi]
where xi is the number of moles of each compound, so xi/∑xi is the relative fraction of that compound.
As you can see from the equation above, this is a nonlinear problem (because of the logarithm).
The constraints for chemical equilibrium are the atom balances for each atom: Obviously we cannot produce or destroy basic atoms, they can only form molecules with other atoms or stay on their own.
So for our example below where we subject ½N₂H₄ + ½O₂ to high pressure and temperature, for each mole of inputs we have:
- ½N₂ = 1 mole of nitrogen
- ½H₄ = 2 moles of hydrogen
- ½O₂ = 1 mole of oxygen
Chemical equilibrium problem
The problem considered is the determination of the equilibrium composition
resulting from subjecting ammonium (NH₄) + oxygen (½O₂) to a temperature
of 3500°K and a pressure of 750 psi (51.7 bar). In the table below we show for each compound
*j* of 10 possible compounds (where the monotonic atoms are termed compounds)
the Gibbs free energy function (*F°ⱼ/RT*)ⱼ, the computed value of *c*ⱼ
for *P* = 750 psi, and the number of atoms of H, N, and O per molecule.
The number of atomic weights of H, N, and O in the mixture are assumed
to be *b*₁ = 2, *b*₂ = 1, and *b*₃ = 1.
Data:
Data on ½N₂H₄ + ½O₂ at 3500°K, 750 psi
| *j* | Compound | (*F°/RT*)ⱼ | *c*ⱼ | *i* = 1 (H) | *i* = 2 (N) | *i* = 3 (O) |
|-----|----------|------------|--------|:-----------:|:-----------:|:-----------:|
| 1 | H | −10.021 | −6.089 | 1 | | |
| 2 | H₂ | −21.096 | −17.164 | 2 | | |
| 3 | H₂O | −37.986 | −34.054 | 2 | | 1 |
| 4 | N | −9.846 | −5.914 | | 1 | |
| 5 | N₂ | −28.653 | −24.721 | | 2 | |
| 6 | NH | −18.918 | −14.986 | 1 | 1 | |
| 7 | NO | −28.032 | −24.100 | | 1 | 1 |
| 8 | O | −14.640 | −10.708 | | | 1 |
| 9 | O₂ | −30.594 | −26.662 | | | 2 |
| 10 | OH | −26.111 | −22.179 | 1 | | 1 |
Learnings
The main difference to previous examples is the selection of the solver. For nonlinear problems, Feasible will use Ipopt. As all other used solvers, it is open source.
Ipopt uses information from the first and second derivative of its objective function. For linear or quadratic programs, the derivates could easily be calculated by the solver from just the coefficients.
For example, the first derivative of 2x is 2, and the second derivative is 0.For nonlinear programs with arbitrary functions, calculation is more complex; specifically, the chain and product rules are required.
As another example, the first derivative of x⋅logx is logx+1, and the second derivative is then 1/x.Therefore, the derivates need to be passed to Ipopt as well. Feasible, through it’s JuMP interface,
calculates these derivatives using forward mode automatic differentiation.
Since Ipopt is based on information from the derivatives, it cannot distinguish between a local minimum and a global one. This is why its return code on success is LOCALLY_SOLVED.
For some problems, namely convex ones, one can prove that this local optimum is in fact a global one. The function in this example, x⋅logx, is convex for x>0. For x<0, it is not defined.
This can be seen from the second derivative 1/x>0 for all x>0. So we know that this local minimum is in fact a global one.