Handling inert species in Gibbs Reactor

[JavalVyas2000]

The inert_species_balance function in gibbs reactor does not create the balance for all the inert species, but only for the ones that are linearly independent with element balance as per the logic in the function. This increases the degrees of freedom as the variables are created for all the species including inerts but there are no constraints addressing the linearly independent with element balance.

The previous version of gibbs reactor use to create the constraints for all the inert species commit id - 119656f

To give an example of this, look at the below example

from idaes.models.unit_models import GibbsReactor
import pyomo.environ as pyo
from idaes.models.properties.modular_properties.base.generic_property import (
    GenericParameterBlock,
)
from idaes.models.properties.modular_properties.base.generic_reaction import (
    GenericReactionParameterBlock,
)
from idaes.models_extra.power_generation.properties.natural_gas_PR import (
    get_prop,
    get_rxn,
)

m=pyo.ConcreteModel()
m.fs=FlowsheetBlock(dynamic=False)
NG_config = get_prop(
        components=[
            "H2",
            "CO",
            "H2O",
            "CO2",
            "CH4",
            "C2H6",
            "C3H8",
            "C4H10",
            "N2",
            "O2",
            "Ar",
        ]
    )
m.fs.NG_props = GenericParameterBlock(**NG_config)
m.fs.gr=GibbsReactor(has_heat_transfer=True,
        has_pressure_change=True,
        inert_species=["N2", "Ar", "O2"],
        property_package=m.fs.NG_props)

### This is the inert_species_balance function
inert_species=["N2", "Ar", "O2"]
for j in inert_species:
    cv = m.fs.gr.control_volume
    e_comp = cv.properties_out[0].config.parameters.element_comp

    # Check for linear dependence with element balances
    # If an inert species is the only source of element e,
    # the inert species balance would be linearly dependent on the
    # element balance for e.
    dependent = True

    if len(m.fs.gr.control_volume.properties_in.phase_list) > 1:
        # Multiple phases avoid linear dependency
        dependent = False
    else:
        for e in m.fs.gr.config.property_package.element_list:
            if e_comp[j][e] == 0:
                # Element e not in component j, no effect
                continue
            else:
                for i in m.fs.gr.control_volume.properties_in.component_list:
                    if i == j:
                        continue
                    else:
                        # If comp j shares element e with comp i
                        # cannot be linearly dependent
                        if e_comp[i][e] != 0:
                            print(j,e)
                            dependent = False
    print(j,dependent)
    if (
        not dependent
        and ('Vap', j) in m.fs.gr.control_volume.properties_in.phase_component_set
    ):
        print(j)
        # return 0 == (
        #     cv.properties_in[t].get_material_flow_terms(p, j)
        #     - cv.properties_out[t].get_material_flow_terms(p, j)
        # )
    else:
        print(j)
        # return Constraint.Skip

This gives shows that the model does not create constraints for N2 and Ar as they are linearly independent, but creates constraints for O2 as it has O. This makes the dof to increase by 2 of the flowsheet for each Gibbs Reactor.

[andrewlee94]

So, my first reaction here is that this is user error. How can O2 be inert if you also have CO2, CO, H2O, and hydrocarbons? What happens if O2 is removed from the inert species list?

[JavalVyas2000]

I checked the DOF after removing O2 and it was the same.

[JavalVyas2000]

Another possible issue I found with this is that the Gibbs Reactor makes the lagrange_mult variable for all the elements in the property package but does not consider if the element is linearly independent. As the linearly independent elements have the lagrange_mult as a free variable. This made the DOF to be 0 in this case.