How can I modify my code from explicit to implicit method ?
3 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Hello all,
I have written my code explicitly. I would like to see how my code works implicitly. Can someone help me share their thoughts about it ?
Thank you
clear; clc; close all;
W =5; L=0.25;
Nx=10; Ny=10;
x = linspace(0,L,Nx);
y = linspace(0,W,Ny);
[X, Y] =meshgrid(x,y);
dx = L/Nx; dy = W/Ny;
k1 = dx/dy;
k2 = dx^2/dy^2;
rhoo = 0.8; mu =0.8E-5;
per = 1e-11; %1e-09
u = zeros(Ny, Nx);
v = zeros(Ny, Nx);
P = 101325 * ones(Ny, Nx);
% Boundary conditions
u(1,:) = u(2,:); %0 % left
v(1,:) = v(2,:); % left
P(1,:) = 101325;
u(:,1) = u(:,2); %0 % top - open to atmosphere
v(:,1) = v(:,2); % top - open to atmosphere
P(:,1) = 101325;
u(Ny,:) = u(Ny-1,:); % right - symmetry
v(Ny,:) = v(Ny-1,:); % right - symmetry
P(Ny,:) = P(Ny-1,:); % right - symmetry
u(:,Nx) = 0; % bottom
v(:,Nx) = 0; % bottom
P(:,Nx) = P(:,Nx-1);
S = [
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
0.1458 0.1929 0.2123 0 0 0 0 0 0 0
];
uold = u; vold = v; Pold = P;
err_uvelocity = 1; err_vvelocity = 1; err_pressure = 1;
while err_uvelocity > 1e-10 || err_vvelocity > 1e-10 || err_pressure > 1e-10
%for test = 1:1
for i=2:Ny-1
for j=2:Nx-1
%%%%%%%%%%%%%%%%%%%%%
u(1,:) = u(2,:); %0 % left
v(1,:) = v(2,:); % left
u(:,1) = u(:,2); %0 % top - open to atmosphere
v(:,1) = v(:,2); % top - open to atmosphere
u(Ny,:) = u(Ny-1,:); % right - symmetry
v(Ny,:) = v(Ny-1,:); % right - symmetry
P(Ny,:) = P(Ny-1,:); % right - symmetry
P(:,Nx) = P(:,Nx-1); % bottom
% P(:,1) = 101325;
%%%%%%%%%%%%%%%%%%%%%
%%Calculating pressure
XX = (S(i,j)*mu*(dx^2))/ (per*rhoo) ;
YY = -k2*(P(i,j-1) + P(i,j+1)) ;
ZZ = -P(i-1,j) - P(i+1,j);
denom = -2 -2*k2;
P(i,j) = (XX + YY + ZZ) / denom;
%%%%%%%%%%%%%%%%%%%%%
%%Calculating velocities
%AlongX = (P(i+1,j) - P(i-1,j))/(2*dx);
%AlongY = (P(i,j+1) - P(i,j-1))/(2*dy);
AlongX = (P(i+1,j) - P(i,j))/(dx);
AlongY = (P(i,j+1) - P(i,j))/(dy);
Const = -per/mu;
u(i,j) = Const*AlongX;
v(i,j) = Const*AlongY;
%%%%%%%%%%%%%%%%%%%%%
err_uvelocity = norm(uold(i,j) - u(i,j));
err_vvelocity = norm(vold(i,j) - v(i,j));
err_pressure = norm(Pold(i,j) - P(i,j));
end
end
uold=u;
vold=v;
Pold = P;
%PPrime = P.^2;
end
% contourf(x, y, v);
% colorbar;
## quiver(x,y,u,v);
2 comentarios
Torsten
el 3 de Mzo. de 2023
Editada: Torsten
el 3 de Mzo. de 2023
I don't see anything "explicit" in the code above that could be made "implicit".
Seems this is an iteration between pressure and velocity, maybe a pressure-velocity coupling method like SIMPLE or something similar. It might be true that there are explicit and implicit methods to handle this coupling, but I think it's unrealistic that you will get help for this in a MATLAB forum. Maybe you can ask CFD experts.
Respuestas (0)
Ver también
Categorías
Más información sobre Assumptions en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!