1

我很尴尬地问这个问题，因为我相信我可能遗漏了一些明显的东西，但我就是看不出我哪里出错了。作为一个更大项目的一部分，我正在研究应用离散化方法来近似方波上的对流方程。但是，我注意到在某些情况下，我的方波的边界被错误地应用了。当 10.25<=X<=10.5 时，它的初始条件应该是 1，其他地方应该是 0。这是该问题的一个示例：

L=20;                   %domain length
dx=0.01;                %space step
nx=(L/dx);              %number of steps in space
x=0:dx:(nx-1)*dx;       %steps along x
sqr=zeros(1,nx);        %pre-allocate array space

for j=1:nx
    if ((10.25<=x(j))&&(x(j)<=10.5))  
        sqr(j)=1;
    else
        sqr(j)=0;
end
d=plot(x,sqr,'r','LineWidth',2);    axis([10.1 10.6 0 1.1]);
drawnow;

在这种情况下，波形显示不正确，x=10.5 取值为 0 而不是 1，如下所示：

https://dl.dropboxusercontent.com/u/8037738/project/square_wave.png

奇怪的是，如果我将域长度更改为不同的值，它有时会正确显示。这是当域长度设置为 30 并且正确显示时：

https://dl.dropboxusercontent.com/u/8037738/project/square_wave_correct.png

我真的不明白，因为我的 x 数组总是以 0.01 的间隔离散，所以它在循环时永远不会“错过”10.5。我希望我已经充分解释了这个问题，如果这是我的愚蠢错误，我提前道歉。

4

2 回答 2

1

你忘记了end当你完成if

L=20;                   %domain length
dx=0.01;                %space step
nx=(L/dx);              %number of steps in space
x=0:dx:(nx-1)*dx;       %steps along x
sqr=zeros(1,nx);        %pre-allocate array space

eps = 1.e-8;    

for j=1:nx
    if ((10.25-eps<=x(j))&&(x(j)<=10.5+eps))  
        sqr(j)=1;
    else
        sqr(j)=0;
    end
end
d=plot(x,sqr,'r','LineWidth',2);    axis([10.1 10.6 0 1.1]);
drawnow;

我还建议您在使用与浮点值进行比较时使用较小的容差，例如x(j)<=10.5

如果您在那一点上看到差异，matlab 会告诉您

x(1051)-10.5

ans =

     1.776356839400250e-15

因此，由于四舍五入，x(j) 的值略大于 10.5，与您预期的不相等

于 2014-03-06T22:14:17.343 回答
1

您的问题是 x 值不完全落在空间步长上。例如，如果您查看以下内容：

j=1051;
format long;
x(j)

ans =

10.500000000000002

sebas 的答案通过为每个边界添加容差来解决此问题。在执行 if 语句中的逻辑之前，您还可以尝试将 x(j) 舍入到最接近的第 100 位：

roundn(x(j), -2)
于 2014-03-06T22:24:21.093 回答