Asking for thermal configuration and reaction configuration in Haber-Bosch process

[KevinChengKaiLai]

I am designing a process system of Haber-Bosch process making ammonia with the recycle system and a flash for ammonia separation. However, I am confusing on how to defining the physical property of my reactor, especially in vapor liquid equilibrium form, I am using FTPX as my state definition while I reduce the degree of freedom to zero my feet initialization doesn't work therefore, I believe I set something wrong inside my thermal configuration.

Below is the thermal configuration dict:

thermal_config = {
    # Specifying components
    "components": {
        "nitrogen": {
            "type": Component,
            "valid_phase_types": [PhaseType.vaporPhase],
            "elemental_composition": {"N": 2},
            "dens_mol_liq_comp": Perrys,
            "enth_mol_liq_comp": Perrys,
            "enth_mol_ig_comp": RPP4,
            "pressure_sat_comp": RPP4,
            "phase_equilibrium_form": {("Vap", "Liq"): fugacity},
            # "visc_d_phase_comp": {"Vap": ChapmanEnskogLennardJones},
            # "therm_cond_phase_comp": {"Vap": Eucken},
            "parameter_data": {
                "mw": (28.0134e-3, pyunits.kg / pyunits.mol),
                "pressure_crit": (3.39e6, pyunits.Pa),
                "temperature_crit": (126.2, pyunits.K),
                "omega": 0.039,  # Added from gas_property.py
                "dens_mol_liq_comp_coeff": {
                    "eqn_type": 1,
                    "1": (3.2091, pyunits.kmol * pyunits.m**-3),
                    "2": (0.2861, None),
                    "3": (126.2, pyunits.K),
                    "4": (0.2966, None),
                },
                "cp_mol_ig_comp_coeff": {
                    "A": (3.115e1, pyunits.J / pyunits.mol / pyunits.K),
                    "B": (-1.357e-2, pyunits.J / pyunits.mol / pyunits.K**2),
                    "C": (2.680e-5, pyunits.J / pyunits.mol / pyunits.K**3),
                    "D": (-1.168e-8, pyunits.J / pyunits.mol / pyunits.K**4),
                },
                "cp_mol_liq_comp_coeff": {
                    "1": (2.8197e5, pyunits.J / pyunits.kmol / pyunits.K),
                    "2": (-1.2281e4, pyunits.J / pyunits.kmol / pyunits.K**2),
                    "3": (2.4800e2, pyunits.J / pyunits.kmol / pyunits.K**3),
                    "4": (-2.2182, pyunits.J / pyunits.kmol / pyunits.K**4),
                    "5": (7.4902e-3, pyunits.J / pyunits.kmol / pyunits.K**5),
                },
                "enth_mol_form_liq_comp_ref": (
                    8.669e3,
                    pyunits.J / pyunits.mol,
                ),
                "enth_mol_form_vap_comp_ref": (
                    0,
                    pyunits.J / pyunits.mol,
                ),
                "pressure_sat_comp_coeff": {
                    "A": (-6.09676, None),
                    "B": (1.13670, None),
                    "C": (-1.04072, None),
                    "D": (-1.93306, None),
                },
                # Added transport properties from gas_property.py
                # "lennard_jones_sigma": (3.798, pyunits.angstrom),
                # "lennard_jones_epsilon_reduced": (71.4, pyunits.K),
                "f_int_eucken": 1,
            },
        },
        "hydrogen": {
            "type": Component,
            "valid_phase_types": [PhaseType.vaporPhase],
            "elemental_composition": {"H": 2},
            "dens_mol_liq_comp": Perrys,
            "enth_mol_liq_comp": Perrys,
            "enth_mol_ig_comp": RPP4,
            "pressure_sat_comp": RPP4,
            "phase_equilibrium_form": {("Vap", "Liq"): fugacity},
            # "visc_d_phase_comp": {"Vap": ChapmanEnskogLennardJones},
            # "therm_cond_phase_comp": {"Vap": Eucken},
            "parameter_data": {
                "mw": (2.016e-3, pyunits.kg / pyunits.mol),
                "pressure_crit": (13e6, pyunits.Pa),
                "temperature_crit": (32.97, pyunits.K),
                "omega": -0.218,  # Added from gas_property.py
                "dens_mol_liq_comp_coeff": {
                    "eqn_type": 1,
                    "1": (5.414, pyunits.kmol * pyunits.m**-3),
                    "2": (0.34893, None),
                    "3": (33.19, pyunits.K),
                    "4": (0.2706, None),    
                },
                "cp_mol_ig_comp_coeff": {
                    "A": (2.714e1, pyunits.J / pyunits.mol / pyunits.K),
                    "B": (9.274e-3, pyunits.J / pyunits.mol / pyunits.K**2),
                    "C": (-1.381e-5, pyunits.J / pyunits.mol / pyunits.K**3),
                    "D": (7.645e-9, pyunits.J / pyunits.mol / pyunits.K**4),
                },
                "cp_mol_liq_comp_coeff": {
                    "1": (6.6653e1, pyunits.J / pyunits.kmol / pyunits.K),
                    "2": (6.7659e3, pyunits.J / pyunits.kmol / pyunits.K**2),
                    "3": (-1.2363e2, pyunits.J / pyunits.kmol / pyunits.K**3),
                    "4": (4.7827e2, pyunits.J / pyunits.kmol / pyunits.K**4),
                    "5": (0, pyunits.J / pyunits.kmol / pyunits.K**5),
                },
                "enth_mol_form_liq_comp_ref": (
                    7.9252e3,
                    pyunits.J / pyunits.mol,
                ),
                "enth_mol_form_vap_comp_ref": (
                    0,
                    pyunits.J / pyunits.mol,
                ),
                "pressure_sat_comp_coeff": {
                    "A": (-5.57929, None),
                    "B": (2.60012, None),
                    "C": (-0.85506, None),
                    "D": (1.70503, None),
                },
                # Added transport properties from gas_property.py
                # "lennard_jones_sigma": (2.826, pyunits.angstrom),
                # "lennard_jones_epsilon_reduced": (59.7, pyunits.K),
                "f_int_eucken": 1,
            },
        },
        "ammonia": {
            "type": Component,
            "valid_phase_types": [PhaseType.vaporPhase, PhaseType.liquidPhase],
            "elemental_composition": {"H": 3, "N": 1},
            "dens_mol_liq_comp": Perrys,
            "enth_mol_liq_comp": Perrys,
            "enth_mol_ig_comp": NIST,
            "pressure_sat_comp": RPP4,
            "phase_equilibrium_form": {("Vap", "Liq"): fugacity},
            # "visc_d_phase_comp": {"Vap": ChapmanEnskogLennardJones},
            # "therm_cond_phase_comp": {"Vap": Eucken},
            "parameter_data": {
                "mw": (17.031e-3, pyunits.kg / pyunits.mol),
                "pressure_crit": (11.3e6, pyunits.Pa),
                "temperature_crit": (405.56, pyunits.K),
                "omega": 0.253,  # Added typical value for ammonia
                "dens_mol_liq_comp_coeff": {
                    "eqn_type": 1,
                    "1": (3.5383, pyunits.kmol * pyunits.m**-3),
                    "2": (0.25443, None),
                    "3": (405.65, pyunits.K),
                    "4": (0.2888, None),    
                },
                "cp_mol_ig_comp_coeff": {
                    "A": (19.99563, pyunits.J / pyunits.mol / pyunits.K),
                    "B": (49.77119, pyunits.J / pyunits.mol / pyunits.K / pyunits.kK),
                    "C": (-15.37599,
                        pyunits.J / pyunits.mol / pyunits.K / pyunits.kK**2,
                    ),
                    "D": (
                        1.921168,
                        pyunits.J / pyunits.mol / pyunits.K / pyunits.kK**3,
                    ),
                    "E": (
                        0.189174,
                        pyunits.J / pyunits.mol / pyunits.K / pyunits.kK**-2,
                    ),
                    "F": (-53.30667, pyunits.kJ / pyunits.mol),
                    "G": (203.8591, pyunits.J / pyunits.mol / pyunits.K),
                    "H": (-45.89806, pyunits.kJ / pyunits.mol),
                },
                "cp_mol_liq_comp_coeff": {
                    "1": (6.1289e1, pyunits.J / pyunits.kmol / pyunits.K),
                    "2": (8.0925e4, pyunits.J / pyunits.kmol / pyunits.K**2),
                    "3": (7.994e2, pyunits.J / pyunits.kmol / pyunits.K**3),
                    "4": (-2.6510e3, pyunits.J / pyunits.kmol / pyunits.K**4),
                    "5": (0, pyunits.J / pyunits.kmol / pyunits.K**5),
                },
                "enth_mol_form_liq_comp_ref": (
                    -80.882e3,
                    pyunits.J / pyunits.mol,
                ),
                "enth_mol_form_vap_comp_ref": (
                    -45.9e3,
                    pyunits.J / pyunits.mol,
                ),
                "pressure_sat_comp_coeff": {
                    "A": (-18.122, None),
                    "B": (10.596, None),
                    "C": (-18.908, None),
                    "D": (-32.728, None),
                },
                # Added transport properties (typical values for ammonia)
                # "lennard_jones_sigma": (3.15, pyunits.angstrom),
                # "lennard_jones_epsilon_reduced": (320.0, pyunits.K),
                "f_int_eucken": 1,
            },
        },
    },
    
    # Specifying phases with both ideal and cubic EOS options
    "phases": {
        "Vap": {
            "type": VaporPhase,
            "equation_of_state": Cubic,  # Can be switched to Ideal if needed
            "equation_of_state_options": {"type": CubicType.PR},
            # "visc_d_phase": ViscosityWilke,
            "transport_property_options": {
                "viscosity_phi_ij_callback": wilke_phi_ij_callback,
            },
            "therm_cond_phase": ThermalConductivityWMS,
        },
        "Liq": {
            "type": LiquidPhase,
            "equation_of_state": Cubic,
            "equation_of_state_options": {"type": CubicType.PR},
            # "visc_d_phase": NoMethod,
            # "therm_cond_phase": NoMethod,
        }
    },
    
    # Vapor-liquid equilibrium configuration
    "phases_in_equilibrium": [("Vap", "Liq")],
    "phase_equilibrium_state": {("Vap", "Liq"): SmoothVLE},
    
    # Set base units of measurement
    "base_units": {
        "time": pyunits.s,
        "length": pyunits.m,
        "mass": pyunits.kg,
        "amount": pyunits.mol,
        "temperature": pyunits.K,
    },
    
    # Specifying state definition (keeping from original but making hybrid)
    "state_definition": FTPx,  # Can be switched to FTPx if needed
    "state_bounds": {
        "flow_mol": (0, 100, 1000, pyunits.mol / pyunits.s),
        "temperature": (273.15, 298.15, 2000, pyunits.K),
        "pressure": (5e4, 1e5, 1e8, pyunits.Pa),
    },
    "pressure_ref": (1e5, pyunits.Pa),
    "temperature_ref": (298.15, pyunits.K),
    
    # Binary interaction parameters
    "parameter_data": {
        "PR_kappa": {
            ("nitrogen", "nitrogen"): 0,
            ("nitrogen", "hydrogen"): 0,
            ("nitrogen", "ammonia"): -0.2090,

            ("hydrogen", "hydrogen"): 0,
            ("hydrogen", "nitrogen"): 0,
            ("hydrogen", "ammonia"): -0.13400,

            ("ammonia", "ammonia"): 0,
            ("ammonia", "nitrogen"): 0,
            ("ammonia", "hydrogen"): 0,
        }
    }
}

And the Reaction Configuration :

reaction_config = {
    "base_units": {
        "time": pyunits.s,
        "length": pyunits.m,
        "mass": pyunits.kg,
        "amount": pyunits.mol,
        "temperature": pyunits.K,
    },
    "rate_reactions": {
        "NH3_synthesis": {
            "stoichiometry": {
                ("Vap", "nitrogen"): -1,
                ("Vap", "hydrogen"): -3,
                ("Vap", "ammonia"): 2
            },
            "heat_of_reaction": constant_dh_rxn,
            "rate_constant": arrhenius,
            "rate_form": power_law_rate,
            "concentration_form": ConcentrationForm.moleFraction,
            "parameter_data":{

                "reaction_order": {("Vap", "hydrogen"): 1,
                                ("Vap", "nitrogen"): 1}, 
                "dh_rxn_ref": (-92.28e3, pyunits.J / pyunits.mol), # DONE
                "arrhenius_const": (1.5638e-9, 
                                     pyunits.mol *pyunits.m**-3 * pyunits.s**-1), # TODO
                "energy_activation": (40,765, pyunits.J / pyunits.mol), # DONE
            }
        }    
    }
}

Please let me know if I miss any property or define something wrong, thanks!!

[KevinChengKaiLai]

for the parameter blocks

# Add property  parameter blocks.
m.fs.thermo_params = GenericParameterBlock(**thermal_config)

# Add reaction parameter blocks.
m.fs.reaction_params = GenericReactionParameterBlock(
    property_package=m.fs.thermo_params, **reaction_config)

m.fs.Feed = Feed(property_package=m.fs.thermo_params)
m.fs.s01 = Arc(source=m.fs.Feed.outlet, destination=m.fs.C101.inlet)
...

there's no error here, however when I try to initialize:

m.fs.Feed.initialize()
propagate_state(arc=m.fs.s01)

The error jumps off and showing with the error log below:

• [INFO] idaes.init.fs.Feed.properties: Bubble, dew, and critical point initialization: infeasible - Converged to a locally infeasible point. Problem may be infeasible..
• [INFO] idaes.init.fs.Feed.properties: Equilibrium temperature initialization completed.
• [INFO] idaes.init.fs.Feed.properties: Phase equilibrium initialization: infeasible - Converged to a locally infeasible point. Problem may be infeasible..
• [INFO] idaes.init.fs.Feed.properties: Property initialization: maxIterations - Maximum Number of Iterations Exceeded..

InitializationError: fs.Feed.properties failed to initialize successfully. Please check the output logs for more information.

So I think the problem is that I need to do something with Dew point and Critical point and Phase equilibrium?