@Torsten thanks for your answer in my previous question ,your comments lmake me notice that its not necessarut the value 0.5 but its sufficient to have ;u(x_i,y,t_i)=u(x_j,y,t_j) , i also forget to tell you that im comparing my result with the resylts in the follox=wing article :Numerical Solution of the 1 + 2 Dimensional Fisher’s Equation by Finite Elements and the Galerkin Method by J. ROESSLER AND H. H~~SSNER
issue in approximation of wave speed of propagation
16 次查看(过去 30 天)
显示 更早的评论
0
im exploring Fisher KPP equation and diffrent numerical methods to study its behaviour specifically travelling waves and speed of propagation
here is the code in Matlab where iuse finite difference scheme:
C = 0.05; % Courant number
D = 1; % Diffusion coefficient
r = 1; % Reaction rate
Lx = 10; % Width of the domain
Ly = 10; % Height of the domain
nx = 51; % Number of grid points in x-direction
ny = 51; % Number of grid points in y-direction
nt = 5000; % Number of time steps
dx = Lx / (nx-1); % Spacing of grid in x-direction
dy = Ly / (ny-1); % Spacing of grid in y-direction
dt = C * dx; % Time step size
% Initialize concentration profile
un = zeros(ny, nx);
x = linspace(-Lx/2, Lx/2, nx);
y = linspace(-Ly/2, Ly/2, ny);
[X, Y] = meshgrid(x, y);
% Set initial condition
initial_condition = (X.^2 + 2*Y.^2 < 0.25);
un = 0.5 * initial_condition;
% Initialize matrix to save solutions over time
U = zeros(nt, ny, nx);
U(1,:,:) = un;
% Loop over time steps to solve the PDE
for n = 1:nt
uc = un;
% Update concentration profile using finite difference method
for i = 2:nx-1
for j = 2:ny-1
un(j, i) = uc(j, i) + ...
dt * (uc(j-1, i) - 2 * uc(j, i) + uc(j+1, i))/dy^2 + ...
dt * (uc(j, i-1) - 2 * uc(j, i) + uc(j, i+1))/dx^2 + ...
dt * uc(j, i) * (1 - uc(j,i));
end
end
% Apply Dirichlet boundary conditions
un(1,:) = 0;
un(end,:) = 0;
un(:,1) = 0;
un(:,end) = 0;
% Store solution
U(n+1,:,:) = un;
end
% Initialize arrays to store pairs of (x, t)
x_values = [];
t_values = [];
% Loop over time steps
for ti = 1:nt
for tj = ti+1:nt % to avoid comparing the same time step twice
% Extract concentration profiles at times ti and tj
ui = squeeze(U(ti,:,:));
uj = squeeze(U(tj,:,:));
% Find pairs (xi, xj) where ui(xi, y) = uj(xj, y) for all y
[xi, xj] = find(all(abs(ui - uj) < 1e-2, 1)); % tolerance for numerical solution
% Store pairs (xi, ti) and (xj, tj)
for k = 1:length(xi)
x_values = [x_values; x(xi(k)), x(xj(k))];
t_values = [t_values; ti, tj];
end
end
end
% Plot pairs of (ti, tj) against (xi, xj)
figure;
plot(t_values, x_values, 'o');
xlabel('Temps');
ylabel('Position x');
title('Paires (ti, tj) où u(xi, y, ti) = u(xj, y, tj) pour tout y');
grid on;
i'm expecting that the slope graph equal or large than 2 but with this code i dont get that ? what is the issue? the slope graph her is the approximate speed whixh u calculte by :c=(x_j-x_i)/(t_j-t_i) such that u(x_i,y,t_i)=u(x_j,y,t_j)can you please tell me the issue and if you have another suggestion please tell me
4 个评论
回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Fit Postprocessing 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!