Gibbs Reactor

The Gibbs Reactor is a chemical equilibrium reactor that computes outlet compositions by minimizing the total Gibbs free energy of the system. It is used for modeling chemical reactions at thermodynamic equilibrium.

Overview

The GibbsReactor class performs chemical equilibrium calculations using Gibbs free energy minimization with Lagrange multipliers. The reactor automatically determines the equilibrium composition based on:

Key Features

Mathematical Background

The reactor minimizes the objective function:

\[G = \sum_i n_i \left( \mu_i^0 + RT \ln(\phi_i y_i P) \right) - \sum_j \lambda_j \left( \sum_i a_{ij} n_i - b_j \right)\]

Where:

The Newton-Raphson method iteratively solves for compositions and Lagrange multipliers until convergence.

Basic Usage

import neqsim.process.equipment.reactor.GibbsReactor;
import neqsim.process.equipment.stream.Stream;
import neqsim.thermo.system.SystemSrkEos;

// Create inlet stream
SystemSrkEos system = new SystemSrkEos(298.15, 10.0);
system.addComponent("methane", 1.0, "mol/sec");
system.addComponent("oxygen", 2.0, "mol/sec");
system.addComponent("CO2", 0.0, "mol/sec");
system.addComponent("water", 0.0, "mol/sec");
system.setMixingRule(2);

Stream inlet = new Stream("inlet", system);
inlet.run();

// Create and configure reactor
GibbsReactor reactor = new GibbsReactor("combustion reactor", inlet);
reactor.setEnergyMode(GibbsReactor.EnergyMode.ISOTHERMAL);
reactor.setMaxIterations(5000);
reactor.setConvergenceTolerance(1e-6);
reactor.setDampingComposition(0.01);
reactor.run();

// Get results
Stream outlet = (Stream) reactor.getOutletStream();
System.out.println("Outlet temperature: " + outlet.getTemperature("C") + " °C");
System.out.println("Conversion completed: " + reactor.hasConverged());

Configuration Options

Energy Mode

// Isothermal: temperature remains constant
reactor.setEnergyMode(GibbsReactor.EnergyMode.ISOTHERMAL);

// Adiabatic: temperature changes based on reaction enthalpy
reactor.setEnergyMode(GibbsReactor.EnergyMode.ADIABATIC);

// Using string (case-insensitive)
reactor.setEnergyMode("adiabatic");

Solver Parameters

Parameter Method Default Description
Max Iterations setMaxIterations(int) 5000 Maximum Newton-Raphson iterations
Convergence Tolerance setConvergenceTolerance(double) 1e-3 Convergence criterion for delta norm
Damping Factor setDampingComposition(double) 0.05 Step size for composition updates 
reactor.setMaxIterations(10000);
reactor.setConvergenceTolerance(1e-8);
reactor.setDampingComposition(0.001);  // Smaller = more stable, slower

Inert Components

Mark components that should not participate in reactions:

// By name
reactor.setComponentAsInert("nitrogen");
reactor.setComponentAsInert("argon");

// By index
reactor.setComponentAsInert(0);

Database Species

// Use only components present in inlet stream (default)
reactor.setUseAllDatabaseSpecies(false);

// Add all species from Gibbs database (for product prediction)
reactor.setUseAllDatabaseSpecies(true);

Results and Diagnostics

Convergence Status

if (reactor.hasConverged()) {
    System.out.println("Solution converged in " + reactor.getActualIterations() + " iterations");
} else {
    System.out.println("Failed to converge. Final error: " + reactor.getFinalConvergenceError());
}

Thermodynamic Results

// Enthalpy of reaction (kJ)
double deltaH = reactor.getEnthalpyOfReactions();

// Temperature change in adiabatic mode (K)
double deltaT = reactor.getTemperatureChange();

// Reactor power (W, kW, or MW)
double powerW = reactor.getPower("W");
double powerKW = reactor.getPower("kW");

Mass Balance Verification

// Check mass balance closure
double massError = reactor.getMassBalanceError();  // Percentage error
boolean balanced = reactor.getMassBalanceConverged();  // True if error < 0.1%

// Element-wise balance
double[] elementIn = reactor.getElementMoleBalanceIn();
double[] elementOut = reactor.getElementMoleBalanceOut();
double[] elementDiff = reactor.getElementMoleBalanceDiff();
String[] elementNames = reactor.getElementNames();  // ["O", "N", "C", "H", "S", "Ar", "Z"]

Molar Flows

List<Double> inletMoles = reactor.getInletMoles();
List<Double> outletMoles = reactor.getOutletMoles();

Specialized Reactor: GibbsReactorCO2

For CO2/acid gas systems, use GibbsReactorCO2 which provides pre-configured reaction pathways:

import neqsim.process.equipment.reactor.GibbsReactorCO2;

GibbsReactorCO2 acidGasReactor = new GibbsReactorCO2("acid gas reactor", inlet);
acidGasReactor.run();

Important Limitations of GibbsReactorCO2:

Troubleshooting

Convergence Issues

  1. Reduce damping factor: Try setDampingComposition(0.001) or smaller
  2. Increase iterations: Use setMaxIterations(20000)
  3. Check initial compositions: Ensure products have small non-zero initial amounts
  4. Mark inerts: Components that don’t react should be marked as inert

Mass Balance Errors

If mass balance doesn’t close:

Numerical Instabilities

For stiff systems:

reactor.setDampingComposition(0.0001);  // Very small steps
reactor.setMaxIterations(50000);        // Allow more iterations
reactor.setConvergenceTolerance(1e-4);  // Relax tolerance slightly

Gibbs Database

The reactor uses thermodynamic data from CSV files in src/main/resources/data/GibbsReactDatabase/:

Supported Elements

The reactor tracks mass balance for: O, N, C, H, S, Ar, Z (charge)

Adding Custom Components

Custom components can be added to the database files following the existing format. Each component requires:

See Also