我尝试遵循本教程:离散泊松方程的数值解并添加边界条件,但我无法正确设置轴。然后我修改了一些变量,它似乎已经开始工作了,但我不是 100% 确定的。
这是我要解决的问题:
u_xx + u_yy = f(x,y), 0 <= x <= 1, 0 <= y <= 1.
u(x,y)=lower(x), if y=0
u(x,y)=upper(x), if y=1
u(x,y)=left(y), if x=0
u(x,y)=right(y), if x=1
Solver 函数如下:(让我担心的地方标有 %WARNING%):
function [ hor_padded ] = poisson_solver( f, upper, lower, left, right ,m, n)
%m, n - Number of points (same in x and y direction)
x_a=0;
x_b=1;
y_a=0;
y_b=1;
dx = (x_b-x_a)/m; %dx=dy
x=linspace(x_a,x_b,m);
y=linspace(y_a,y_b,n);
g=zeros(m,n);
for i=1:n,
for j=1:m,
%WARNING 1%
g(j,i)=-f(x(i),y(j))*dx*dx;
if j==2
g(j,i)=g(j,i)+lower(x(i));
end
if j==n-1
g(j,i)=g(j,i)+upper(x(i));
end
if i==2
g(j,i)=g(j,i)+left(y(j));
end
if i==n-1
g(j,i)=g(j,i)+right(y(j));
end
end
end
b = reshape(g(2:m-1, 2:n-1),(m-2)*(n-2),1);
A = gallery('poisson',m-2);
U = A\b;
u = reshape(U,m-2,n-2);
upper_part=upper(x(2:m-1));
lower_part=lower(x(2:m-1));
left_part = left(y');
right_part = right(y');
%WARNING 2%
ver_padded = vertcat(lower_part, u, upper_part);
hor_padded = horzcat(left_part, ver_padded, right_part);
end
- 警告 1) Matlab 似乎在矩阵内混合了 x 和 y 轴。(如此处所述:Axis labeling question)。所以我在for循环中改变了X和Y方向。它是否正确?
- 警告 2)我一开始就有了 upper_part 和 lower_part,但结果图中的边界没有相遇。我检查了变量,尝试了相反的方法,它开始工作了。它是否正确?
我的测试程序在这里:
clear;
m = 100;
n = 100;
x_a=0;
x_b=1;
y_a=0;
y_b=1;
e=exp(1);
f=@(x,y)(e*x+e*y);
%Boundary conditions
upper = @(x) (x);
lower = @(x) (3-x);
left = @(y) (3-3*y);
right = @(y) (2-y);
u = poisson_solver(f, upper, lower, left, right, m, n);
x=linspace(x_a,x_b,m);
y=linspace(y_a,y_b,n);
surf(x,y,u);
xlabel('x');
ylabel('y');
还有另一个测试用例:(其他一切都和上面的例子一样):
upper = @(x) (exp(x));
lower = @(x) (cos(x*2*pi));
left = @(y) (cos(y*2*pi));
right = @(y) (exp(y));
名称“上”、“下”、“左”和“右”指的是结果图,而不是求解器函数中矩阵 u 中的相应位置。