用差分法求出数值解之后,求与精确解误差出问题了。

1 vue (au cours des 30 derniers jours)
raecsha
raecsha le 11 Oct 2022
求误差得到的结果是原函数的图不知道那里的语句出问题了,图片左为数值解图像右边误差除了问题,求解
程序如下:
clear all; close all;
a=4; h=0.05; x=[0:h:1];
tau=0.0125; t=[0:tau:1];
s=a*tau/h;
N=length(x)-1; M=length(t)-1;
[T X]=meshgrid(t,x);
%构造矩阵
e=s^2*ones(N-1,1);
A=spdiags([e 2*(1-e) e],[-1 0 1],N-1,N-1);
%设置初始条件和边界条件
u=zeros(N+1,M+1);
u(:,1)=sin(pi*x); u(:,2)=0;
u(1,:)=0; u(end,:)=0;
%有限差分
for n=2:M
u(2:N,n+1)=A*u(2:N,n)-u(2:N,n-1);
u(2,n+1)=u(2,n+1)+s^2*u(1,n);
u(N,n+1)=u(N,n+1)+s^2*u(end,n);
end
%画图
subplot(2,2,1)
mesh(t,x,u), view(10,10), xlabel t, ylabel x, zlabel u
subplot(2,2,2)
%与解析解的误差
m=sin(pi*X).*cos(4*pi*T)
mesh(t,x,u-m), view(10,8), axis([0 1 0 1 -10 10])
xlabel t, ylabel x, zlabel Error

Réponse acceptée

华纳娱乐游戏网址【359663.tv】
直接写了个对的,也没写成复杂的矩阵运算格式,你去看看我推荐的讲偏微分方程的书籍的对应章节再对照代码吧。
clear; clc; close all;
c = 4;
X_min = 0; X_Step = 1/10; X_max = 1;
T_min = 0; T_Step = 0.025; T_max = 1;
F = @( x ) sin( pi * x );
G = @( x ) 0 * x;
X_Vector = X_min : X_Step : X_max;
T_Vector = T_min : T_Step : T_max;
lambda = c * T_Step / X_Step;
U_Matrix = NaN * ones( numel( X_Vector ), numel( T_Vector ) );
U_Matrix( 1, : ) = 0; U_Matrix( end, : ) = 0;
U_Matrix( :, 1 ) = F( X_Vector ).';
U_Matrix( 2 : end - 1, 2 ) = ( lambda^2 * ( F( X_Vector( 1 : end - 2 ) ) + F( X_Vector( 3 : end ) ) ) / 2 + ...
    ( 1 - lambda^2 ) * F( X_Vector( 2 : end - 1 ) )  + T_Step * G( X_Vector( 2 : end - 1 ) ) ).';
for jj = 3 : 1 : numel( T_Vector )
    U_Matrix( 2 : end - 1, jj ) = lambda^2 * ( U_Matrix( 1 : end - 2, jj - 1 ) + U_Matrix( 3 : end, jj - 1 ) ) + ...
        ( 1 - lambda^2 ) * 2 * U_Matrix( 2 : end - 1, jj - 1 )  - U_Matrix( 2 : end - 1, jj - 2 );
end
[ X_Matrix, T_Matrix ] = meshgrid( X_Vector, T_Vector );
U_Analytic = @( x, t ) sin( pi * x ) .* cos( c * pi * t );
U_AnaMatrix = U_Analytic( X_Matrix, T_Matrix );
%%
figure;
subplot( 1, 2, 1 )
surf( X_Matrix, T_Matrix, U_Matrix.',...
    'FaceColor', 'r', 'FaceAlpha', 0.7, 'MeshStyle', 'default', 'EdgeColor', 'k', 'EdgeAlpha', 0.6, ...
    'AlignVertexCenters', 'on', 'LineStyle', ':', 'LineWidth', 0.8 );
xlabel( 'x', 'FontSize', 18 ); ylabel( 't', 'FontSize', 18 ); zlabel( 'u_{numeric}', 'FontSize', 18 );
axis equal; grid on;
subplot( 1, 2 ,2 )
surf( X_Matrix, T_Matrix, U_AnaMatrix,...
    'FaceColor', 'b', 'FaceAlpha', 0.7, 'MeshStyle', 'default', 'EdgeColor', 'k', 'EdgeAlpha', 0.6, ...
    'AlignVertexCenters', 'on', 'LineStyle', ':', 'LineWidth', 0.8 );
xlabel( 'x', 'FontSize', 18 ); ylabel( 't', 'FontSize', 18 ); zlabel( 'u_{analytic}', 'FontSize', 18 );
axis equal; grid on;

Plus de réponses (0)

Catégories

En savoir plus sur Logical dans Help Center et File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!