I am trying to write the MATLAB code of the question ofresearch paper Drive link of paper-http​s://drive.​google.com​/file/d/12​N5ugCymHKB​ubfqOlM2LU​3DVBxdrCoU​R/view?usp​=drive_lin​k

1 次查看(过去 30 天)
function batch_distillation()
% Binary Data
ki = [1.51; 1.13; 0.95; 0.91];
alphaij = [5.91, 2.51, 1.29, 1.00;
1/5.91, 1, 1/2.51, 1/1.29;
1/2.51, 1/2.51, 1, 1/1.29;
1/1.29, 1/1.29, 1/1.29, 1];
% Given data
F_feed = 400; % kmol/day
F_i = [100; 100; 100; 100]; % kmol of each component in the feed
x_i = [0.97; 0.97; 0.97; 0.97]; % molar fractions in the feed
K = ki .* x_i; % Equilibrium constants
alpha = alphaij .* (x_i ./ x_i'); % Relative volatilities
% Time parameters
t_span = [0 10]; % Set the time span for batch distillation (adjust as needed)
t_steps = 100; % Number of time steps
% Initial conditions
L0 = 0; % Initial liquid holdup
V0 = 0; % Initial vapor holdup
z_L0 = x_i; % Initial liquid composition (same as feed composition)
z_V0 = x_i; % Initial vapor composition (same as feed composition)
% Solve the batch distillation using ODE solver
options = odeset('RelTol', 1e-6, 'AbsTol', 1e-9);
[t, y] = ode45(@(t,y) batch_distillation_ode(t, y, F_feed, ki, alphaij, x_i, t_span), t_span, [L0; V0; z_L0; z_V0], options);
% Extract results
L = y(:, 1);
V = y(:, 2);
z_L = y(:, 3:6);
z_V = y(:, 7:10);
% Plot results
figure;
plot(t, L, 'b', 'LineWidth', 2);
hold on;
plot(t, V, 'r', 'LineWidth', 2);
xlabel('Time (hours)');
ylabel('Moles (kmol)');
legend('Liquid holdup (L)', 'Vapor holdup (V)');
title('Moles vs Time');
figure;
plot(t, z_L, 'LineWidth', 2);
hold on;
plot(t, z_V, 'LineWidth', 2);
xlabel('Time (hours)');
ylabel('Molar Fraction');
legend('Benzene (L)', 'Toluene (L)', 'Ethylbenzene (L)', 'o-xylene (L)',...
'Benzene (V)', 'Toluene (V)', 'Ethylbenzene (V)', 'o-xylene (V)');
title('Molar Fractions vs Time');
end
function dydt = batch_distillation_ode(t, y, F_feed, ki, alphaij, x_i, t_span)
% Unpack variables
L = y(1);
V = y(2);
z_L = y(3:6);
z_V = y(7:10);
% Calculate liquid and vapor flow rates
F_L = L / t_span(end);
F_V = V / t_span(end);
% Calculate liquid and vapor compositions
x_L = z_L ./ sum(z_L);
x_V = z_V ./ sum(z_V);
% Calculate liquid and vapor phase rates
r_L = F_L .* x_L;
r_V = F_V .* x_V;
% Calculate material balances
dLdt = F_feed - F_L;
dVdt = F_L - F_V;
dz_Ldt = (r_L - z_L) ./ L;
dz_Vdt = (r_V - z_V) ./ V;
% Pack the derivatives
dydt = [dLdt; dVdt; dz_Ldt; dz_Vdt];
end
  3 个评论
Aditya
Aditya 2023-11-17
It seems that due to your initialization all the values are coming to be INF for z_L and z_V, please do check your initialization. if the issue still persists please do share the paper that you are refering to.

请先登录,再进行评论。

回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Thermal Analysis 的更多信息

产品


版本

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by