0

このチュートリアルに従ってみました:離散ポアソン方程式の数値解法と境界条件の追加を試みましたが、軸を正しく取得できませんでした。次に、いくつかの変数をハックして、機能し始めたようですが、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

ソルバー関数は次のとおりです: (私が心配する場所は %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. 警告 1) Matlab は行列内で x 軸と y 軸を混在させているようです。(ここで説明されているように:軸のラベル付けの質問)。そのため、for サイクルで X 方向と Y 方向を変更しました。これは正しいです?
  2. 警告 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));

"upper"、"lower"、"left"、"right" という名前は、ソルバー関数の行列 u 内の対応する場所ではなく、結果のグラフを指します。

4

1 に答える 1

1

警告 1) Matlab は行列内で x 軸と y 軸を混在させているようです。(ここで説明されているように:軸のラベル付けの質問)。そのため、for サイクルで X 方向と Y 方向を変更しました。これは正しいです?

はい、行列の最初のインデックスは行インデックス、2 番目のインデックスは列インデックスです。

[ (1,1)  (1,2)  . . .   (1,n);
  (2,1)  (2,2)  . . .   (2,n);
    :      :    .         :
    :      :      .       :
    :      :        .     :
  (m,1)  (m,2)  . . .   (m,n)];

警告 2) 最初に upper_part と lower_part を逆にしましたが、結果のグラフで境界が一致しませんでした。変数を調べ、逆の方法で試したところ、機能し始めました。これは正しいです?

行列を xy 平面として定義する場合は、インデックスの定義に固執する必要があります。インデックスの例でわかるように、「原点」( (1,1)) は左上にあり、「y 軸」は「下」を指しています。それを踏まえると、上と下の変更は正しいです。

余談ですが、もちろん、すべてを別の方法で定義することもできます。それは、定義した方法で配置および格納された数字の集まりです!

于 2012-12-10T18:21:36.617 に答える