Can you help me write this question in finite volume method? This is my first time writing FVM. I wrote using FDM below :(

Question

Shengwei Bie 2022-3-8

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1666174-can-you-help-me-write-this-question-in-finite-volume-method-this-is-my-first-time-writing-fvm-i-wr

回答： Brahmadev 2024-2-9

%Upwind scheme
n=320;
X=n+1;
Y=n+1;
u=0.5;
v=0.5;
Lx=4;
Ly=4;
dX=Lx/n;
dY=Ly/n;
c(:,:,1)=zeros(Y,X);
c(:,:,2)=zeros(Y,X);
for a=(2+0.9/dX):1:(1.1/dX)
    for b=(2+0.9/dY):1:(1.1/dY)
        x=(a-1)*dX;
        y=(b-1)*dY;
        c(a,b,1)=sin(5*pi*(x-0.9))*sin(5*pi*(y-0.9));
    end
end
dT=0.001;
for T=1:1:4
    for n=2:1:(round(T/dT)+1)
        for j=2:1:(Lx/dX)
            for i=2:1:(Ly/dY)
                c(i,j,2)=c(i,j,1)+(v*dT/dX)*(c(i-1,j,1)-c(i,j,1))+(u*dT/dY)*(c(i,j-1,1)-c(i,j,1));
                
            end
        end
        c(:,:,1)=c(:,:,2);
    end
    dT=dT+0.001;
    fig = figure();
    mesh(0:dX:Lx,0:dY:Ly,c(:,:,2));
    xlabel('X');
    ylabel('Y');
    zlabel('c');
    grid on;
    axis([0 4 0 4 -1 1]);
end

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Brahmadev 2024-2-9

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1666174-can-you-help-me-write-this-question-in-finite-volume-method-this-is-my-first-time-writing-fvm-i-wr#answer_1405556

在 MATLAB Online 中打开

Hi @Shengwei Bie,

The MATLAB code for First order upwind seems to be on the right track, here are some more things to consider:

The time loop has an increment dT=dT+0.001; at the end, which doesn't make sense in this context since dT is your time step and should remain constant as defined at the beginning.
The boundary conditions are not enforced in your code. You should set the edges of the computational domain to zero at each time step.

Here is a revised version of the code that takes care of these considerations:

% Define the domain and grid size
n = 320;
Lx = 4;
Ly = 4;
dx = Lx / n;
dy = Ly / n;
x = linspace(0, Lx, n+1);
y = linspace(0, Ly, n+1);
[X, Y] = meshgrid(x, y);
% Define the time domain
dt = 0.001;
t_final = 4;
Nt = round(t_final / dt);
% Define the initial condition
c = zeros(n+1, n+1);
for i = 2:n
    for j = 2:n
        if (x(i) > 0.9 && x(i) < 1.1) && (y(j) > 0.9 && y(j) < 1.1)
            c(j, i) = sin(5 * (x(i) - 0.9)) * sin(5 * (y(j) - 0.9));
        end
    end
end
% Define velocities
u = 0.5;
v = 0.5;
% FVM implementation
for t = 1:Nt
    c_new = c;
    
    % Compute fluxes and update the scalar field c
    for i = 2:n
        for j = 2:n
            % Compute fluxes using upwind scheme
            flux_x = u * (c(j, i) - c(j, i-1)) / dx;
            flux_y = v * (c(j, i) - c(j-1, i)) / dy;
            
            % Update c based on net fluxes
            c_new(j, i) = c(j, i) - dt * (flux_x + flux_y);
        end
    end
    
    % Apply boundary conditions
    c_new(1, :) = 0;   % Bottom boundary
    c_new(end, :) = 0; % Top boundary
    c_new(:, 1) = 0;   % Left boundary
    c_new(:, end) = 0; % Right boundary
    
    c = c_new;
    
    % Plot at specific time steps
    if mod(t, Nt / 4) == 0 || t == 1
        figure;
        surf(X, Y, c, 'EdgeColor', 'none');
        title(sprintf('FVM with Upwind Scheme at t = %.3f', t * dt));
        xlabel('x');
        ylabel('y');
        zlabel('c');
        view(3);
        colorbar;
        pause(0.1); % To visualize the evolution
    end
end