Implementing numerical method for PDE

3 次查看（过去 30 天）
schlang
schlang 2022-10-13
编辑： Davide Masiello 2022-10-13
Hello
I am trying to solve the following PDE with intital and boundary conditions such that . I used the second order centered finite difference discrtization for and then want solve the ode system using ode15s
Here is my attempt. When I plot the solution obtained from ode15s and compare it to the exact solution they are different. I am not if I made a mistake somewhere. Help is really appreciated
clc,clear,close all
% parameters
t0 = 0;
T = 1.0;
tspan = [t0 T];
xl = 0;
xr = 1;
m = 20;
x = linspace(xl,xr,m + 1);
dx = 1/m;
Uexact = @(t,x) exp(1i*(x-t));
% initial conditions
U0 = Uexact(0,x)';
U0 = U0(2:end-1);
% solve
fn = @(t,U) ODE(t,U,m,dx);
opts = odeset('RelTol',1e-13, 'AbsTol',1e-15);
[t,U] = ode15s(fn, tspan, U0, opts);
%compare with exact solution
plot(x(2:end-1),U(end,:))
hold on
plot(x(2:end-1),Uexact(T,x(2:end-1)))
function dUdt = ODE(t,U,m,dx)
A = eye(m-1);
A = A * (-2);
A = A + diag(ones(m-2,1),1);
A = A + diag(ones(m-2,1),-1);
A = (1/dx^2) * A;
g = zeros(m-1,1);
g(1) = g(1) + (1/dx^2) * exp(1i*(-1*t));
g(end) = g(end) + (1/dx^2) * exp(1i*(1-t));
dUdt = (1i) * (A*U) + g;
end
Thanks

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

采纳的回答

Davide Masiello
Davide Masiello 2022-10-13
编辑：Davide Masiello 2022-10-13
Ran in:
clear,clc
tspan = [0 1];
N = 100;
x = linspace(0,1,N);
dx = 1/(N-1);
Uexact = @(t,x) exp(1i*(x-t));
U0 = Uexact(0,x);
M = eye(N);
M(1,1) = 0;
M(N,N) = 0;
opts = odeset('Mass',M,'MassSingular','yes');
[t,U] = ode15s(@(t,U)yourPDE(t,U,N,dx), tspan, U0, opts);
plot(x,real(U(end,:)),'k',x(1:4:end),real(Uexact(1,x(1:4:end))),'r.')
xlabel('x')
ylabel('U')
title('At final time')
legend('Numerical','Exact','Location','Best')
plot(x,real(U(1:3:end,:)),'k',x(1:3:end),real(Uexact(t(1:3:end),x(1:3:end))),'r.')
xlabel('x')
ylabel('U')
title('At several times')
function dUdt = yourPDE(t,U,N,dx)
dUdt(1,1) = U(1)-exp(-1i*t);
dUdt(2:N-1,1) = 1i*(U(1:end-2)-2*U(2:end-1)+U(3:end))/dx^2;
dUdt(N,1) = U(end)-exp(1i*(1-t));
end

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

产品


版本

R2022b

Community Treasure Hunt

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

Start Hunting!

Translated by