Process Heat Integration

Calculates the minimum heating and cooling demand of an industrial process, as well as useful excess heat.

Introduction

The function HeatIntegration.m can be used to calculate the minimum heating and cooling demand of an industrial process, as well as the amount of useful excess heat the process produces. The model is based on the work of Papoulias and Grossmann (https://doi.org/10.1016/0098-1354(83)85023-6) and I extended it to determine the amount of useful excess heat from the process.

Short Description of Functions

• HeatIntegration() calculates the minimum heating demand and cooling demand, as well as useful excess heat. Input arguments are the process streams and latent process heats, as well as the minimum temperature difference for heat exchange and a threshold temperature above which excess heat is considered useful.
• getEnthalpy() Calculates the enthalpy for a liquid or gaseous process stream. For liquid process streams, only constant heat capacity values are supported, while for gaseous process streams, constant heat capacity values and the DIPPR equation 107 are supported. Adding other heat capacity models is easy.

Examples

Steam Methane Reforming

This is a simplified version of a steam methane reforming process, assuming methane and liquid water enter the process at 300 K and 10 bar, are heated to 1,000 K, and then react in an isothermal reactor to reach chemical equilibrium. The product mixture is cooled to 700 K to undergo water gas shift reaction (not included).

% Each species is specified by [heat of formation at 298.15 K in kJ/mol,
% [Parameters A, B, C, D, and E for ideal gas heat capacity according to DIPPR equation 107 in kJ/(mol-K)]
% or a constant heat capacity value
properties = struct("Hydrogen", struct("h_0", 0, "Cp_model", "DIPPR107", "Cp_params", [0.027617, 0.00956, 2466, 0.00376, 567.6]), ...
                    "Methane", struct("h_0", -74.6, "Cp_model", "DIPPR107", "Cp_params", [0.033298, 0.079933, 2086.9, 0.041602, 991.96]), ...
                    "Carbon_Monoxide", struct("h_0", -110.53, "Cp_model", "DIPPR107", "Cp_params", [0.029108, 0.008773, 3085.1, 0.0084553, 1538.2]), ...
                    "Water_g", struct("h_0", -241.81, "Cp_model", "DIPPR107", "Cp_params", [0.033363, 0.02679, 2610.5, 0.008896, 1169]), ...
                    "Water_l", struct("h_0", -285.83, "Cp_model", "const", "Cp_params", 0.07538));

% ASSEMBLE STREAMS
% Inlet stream 1
components_1 = [struct("Name", 'Methane', "x", 1, "Props", properties.Methane)];
stream_1 = struct("n", 1, "T_in", 298.15, "T_out", 1000, "Components", components_1, "State", "gas");

% Inlet stream 2
% Assume 10 bar pressure (boiling point of water is 179.91°C) and 1000 K for the reactor
components_2_l = [struct("Name", 'Water', "x", 1, "Props", properties.Water_l)];
stream_2_l = struct("n", 1, "T_in", 298.15, "T_out", 453.06, "Components", components_2_l, "State", "liquid");
components_2_g = [struct("Name", 'Water', "x", 1, "Props", properties.Water_g)];
stream_2_g = struct("n", 1, "T_in", 453.06, "T_out", 1000, "Components", components_2_g, "State", "gas");

% Outlet stream
% ...the equilibrium conversion is 0.41 under these conditions
X = 0.41; % [-]
% ...calculate stream composition
n_CH4_3 = stream_1.n*(1-X);
n_H2O_3 = stream_2_g.n*(1-X);
n_CO_3 = stream_1.n*X;
n_H2_3 = 3*stream_1.n*X;
n_tot_3 = n_CH4_3 + n_H2O_3 + n_CO_3 + n_H2_3;
components_3 = [struct("Name", 'Methane', "x", n_CH4_3/n_tot_3, "Props", properties.Methane), struct("Name", 'Water', "x", n_H2O_3/n_tot_3, "Props", properties.Water_g), struct("Name", 'Carbon_Monoxide', "x", n_CO_3/n_tot_3, "Props", properties.Carbon_Monoxide), struct("Name", 'Hydrogen', "x", n_H2_3/n_tot_3, "Props", properties.Hydrogen)];
stream_3 = struct("n", n_tot_3, "T_in", 1000, "T_out", 700, "Components", components_3, "State", "gas");

streams = [stream_1, stream_2_l, stream_2_g, stream_3];

% LATENT HEAT
% Water evaporation at 179.91°C
Q_vap = getEnthalpy(struct("n", stream_2_g.n, "T", 453.06, "Components", components_2_g, "State", "gas")) - getEnthalpy(struct("n", stream_2_l.n, "T", 453.06, "Components", components_2_l, "State", "liquid"));
% Heat of reaction
Q_R = getEnthalpy(struct("n", stream_3.n, "T", 1000, "Components", components_3, "State", "gas")) - getEnthalpy(struct("n", stream_2_g.n, "T", 1000, "Components", components_2_g, "State", "gas")) - getEnthalpy(struct("n", stream_1.n, "T", 1000, "Components", components_1, "State", "gas"));

latentHeats = [struct("Q", Q_vap, "T", 453.06), struct("Q", -Q_R, "T", 1000)];

% DETERMINE MINIMUM HEAT REQUIREMENT AND USEFUL EXCESS HEAT
dT_min = 10; % [K]
T_use = 373.15; % [K]
[Q_H, Q_C, Q_use] = HeatIntegration(streams, latentHeats, dT_min, T_use);

Where to Find Data

• Standard enthalpy values for gaseous species can be obtained from cccbdb.nist.gov, DIPPR 107 parameters for the heat capacities from Aspen Plus or databases.
• Heat capacity and standard enthalpy values for liquids are reported in the NIST Chemistry WebBook.

Contact

Feel free to reach out to me under [e-mail].

Disclaimer

The functions are provided 'as is', without warranty of any kind, express or implied.