Foam drainage problem using pdepe
4 views (last 30 days)
Show older comments
The foam drainage equation ๐๐๐ ๐๐ก = โ ๐ ๐๐ง ๐1๐
2 ๐ ๐๐๐๐ 2 โ ๐2 ๐พ 2๐
๐๐ 0.5 ๐ ๐๐ง ๐๐ describes the spatio-temporal evolution of liquid in a foam due to gravitation and capillary forces. Herein, ๐๐(๐ฅ,๐ก) is the local liquid content, ๐ก is time, ๐ง is the spatial coordinate in the direction of the gravitational acceleration ๐, ๐ is the dynamic liquid viscosity, ๐พ is the interfacial tension, ๐ is the liquid density and ๐1 = 0.0066, ๐2 = 0.1610.5 are constants. In this task, the effect of bubble size ๐
on the drainage of liquid out of the foam shall be studied. You can assume the following conditions: โข The height ๐ป of the foam is constant โข The bubble size of constant in the whole foam โข At the foamโs bottom edge (๐ง = ๐ป), the liquid fraction is constant ๐๐ ๐ก, ๐ป = 0.36 โข At the foamโs top edge (๐ง = 0), the liquid flux is zero โข At ๐ก = 0, assume ๐๐(0, ๐ง) = 0.36 Multiscale Process Modelling and Process Analysis 2 ๐ง = 0 ๐ง = ๐ป ๐ Deliverable Task 02โFoam drainage โข Implement the foam drainage equation in MATLAB using the pdepe command with the respective initial and boundary conditions. Use the following parameter values โ ๐ = 0.001 ๐๐ โ
๐ , ๐พ = 0.04 ๐/๐, ๐ = 1000 ๐๐/๐3 , ๐ = 9.81 ๐ ๐ 2 โ ๐ป = 0.1 ๐, ฮโ = 0.0001 ๐, ฮ๐ก = 0.01 ๐ , ๐ = 100 ๐ โข Solve the model for ๐
= 0.0005 ๐, ๐
= 0.001 ๐, ๐
= 0.005 ๐ and generate a chart for each solution showing the spatial profile of the liquid fraction at times between ๐ก = 0 and ๐ก = ๐ in steps of ฮ๐ = 0.05๐ โข Compute the cumulated liquid flux out of the foam and compare the time-evolution with respect for each bubble size. โ Use the fact that the drained liquid ๐๐,๐ท can be computed at any time as the difference of the liquid being initially present in the foam and the liquid in the foam at time ๐ก, thus, ๐๐,๐ท(๐ก) = ๐๐,๐น ๐ก0 โ ๐๐,๐น(๐ก) โ Use the relation ๐๐,๐น = ๐ด๐ป (integral lower border z=0 upper border z=H) ๐๐ฟ ๐๐ง, where ๐ด๐น is the foamโs cross-section. The trapz command can be used to calculate the integral (compare to add-on of exercise 6). Assume ๐ด๐น = 0.002 m2
0 Comments
Answers (1)
SOUMNATH PAUL
on 22 Dec 2023
To my understanding of this problem, we need to define the partial differential equation, set the initial and the boundary conditions. In the second step I have solved the equation for the specified bubble sizes and plotted the spatial profile of the liquid fraction over time. Finally, I have computed the cumulated liquid flux out of the foam.
function foam_drainage
% Given parameters
eta = 0.001; % Pa.s
gamma_tension = 0.04; % N/m
rho = 1000; % kg/m^3
g = 9.81; % m/s^2
H = 0.1; % m
A_f = 0.002; % m^2 (cross-section of the foam)
k1 = 0.0066;
k2 = 0.161^0.5;
% Discretization parameters
dz = 0.0001; % m
dt = 0.01; % s
T = 100; % s
dTau = 0.05 * T; % s
% Bubble sizes to study
R_values = [0.0005, 0.001, 0.005]; % m
% Spatial mesh
zmesh = 0:dz:H;
% Time vector
tspan = 0:dt:T;
% Loop over the bubble sizes
for R = R_values
% Solve PDE
sol = pdepe(0, @(z,t,u,dudz) pdefun(z,t,u,dudz,k1,k2,R,eta,rho,g,gamma_tension), ...
@icfun, ...
@bcfun, ...
zmesh, tspan);
% Extract the solution for phi_l
phi_l = sol(:,:,1);
% Plot the spatial profile of the liquid fraction at specified times
figure;
hold on;
for t = 0:dTau:T
[~, tIdx] = min(abs(tspan - t)); % Find the closest time index
plot(zmesh, phi_l(tIdx,:), 'DisplayName', sprintf('t = %.2f s', t));
end
hold off;
title(sprintf('Spatial profile of liquid fraction over time (R = %.4f m)', R));
xlabel('Height z (m)');
ylabel('Liquid fraction \phi_l');
legend('show');
% Compute the cumulated liquid flux out of the foam
V_l_F_initial = A_f * H * trapz(zmesh, phi_l(1,:)); % Initial liquid volume in the foam
V_l_F = A_f * H * trapz(zmesh, phi_l, 2); % Liquid volume in the foam over time
V_l_D = V_l_F_initial - V_l_F; % Drained liquid volume over time
% Plot the time-evolution of the cumulated liquid flux for each bubble size
figure;
plot(tspan, V_l_D, 'DisplayName', sprintf('R = %.4f m', R));
title('Cumulated liquid flux out of the foam over time');
xlabel('Time t (s)');
ylabel('Drained liquid volume V_l_D (m^3)');
legend('show');
end
end
% Define the PDE function
function [c,f,s] = pdefun(z,t,u,dudz,k1,k2,R,eta,rho,g,gamma_tension)
c = 1;
f = -k1*R^2/(eta*rho*g) * u^2 - k2*gamma_tension^2/(2*R) * sqrt(u) * dudz;
s = 0;
end
% Define the initial condition function
function u0 = icfun(z)
u0 = 0.36; % Initial liquid fraction
end
% Define the boundary condition function
function [pl,ql,pr,qr] = bcfun(zl,ul,zr,ur,t)
pl = ul - 0.36; % At the bottom edge (z = H), the liquid fraction is constant
ql = 0;
pr = 0; % At the top edge (z = 0), the liquid flux is zero
qr = 1;
end
You can also see the examples section in below link for further reference:
Hope it helps!
Regards,
Soumnath
0 Comments
See Also
Categories
Find more on Medical Physics in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!