0 
k = 0.019;
Pstar = 100;
H = 33;
h = 0.1;
X = 36;
N = round(X/h);
t = zeros(1,N+1);
P = zeros(1,N+1);
P(1) = 84;
t(1) = 0;
yHeun = zeros(1,N+1);
yHeun(1)=84;
a = 1; b = 100;
while b-a >0.5
    c = (a+b)/2;
    for n = 1:N
        t(n+1) = t(n) + h;
        Inside = nthroot(sin(2*pi*t/12),15);
        Harvest = c*0.5*(Inside+1);
        P(n+1) = P(n) + h*(k*P(n)*(Pstar-P(n))-Harvest(n));
        if P < 0
            P = 0;
        end
        yHeun(n+1) = yHeun(n) + h*0.5*((k*P(n)*(Pstar-P(n))-Harvest(n))+(k*P(n+1)*(Pstar-P(n+1))-Harvest(n+1)));
    end
    if sign(yHeun(c)) == sign(yHeun(a))
        c = a;
    else
        c = b;
    end
end
disp(['The root is between ' num2str(a) ' and ' num2str(b) '.'])

これは私が実行しようとしているコードであり、おそらくひどいことはわかっていますが、コーディングがひどいので、コードを実行しようとするたびに、次のように表示されます。

yHeun(50.5) にアクセスしようとしました。index は正の整数または論理値でなければなりません。

sign(yHeun(c)) == sign(yHeun(a)) の場合、Matlab3Q4 (30 行目) のエラー

yHeun(c または a など) が整数になるものを返すようにする方法がわかりません。while+for ループも正しく実行したとは思いません。

質問: 「H の上限を 100 (高い値を指定すると、36 か月後に母集団が 0 になる) から始めます。下限を 1 にします。上記の問題 #3 のソルバーを while ループの途中に置き、上限と下限の差が 0.5 未満になるまで、H の上限と下限を二等分し続けます。」

4

1 に答える 1

1

30行目（エラーあり）はこれだと思います：

if sign(yHeun(c)) == sign(yHeun(a))

ここで、上記の結果として、は と等しいと思い ます (ところで、デバッグによって推測が正しかったかどうかを確認できます - 30 行目の前に追加してみてください)。c50.5c = (a+b)/2disp(c)

数値を強制的に整数にするには、次を使用しますfloor。

c = floor((a+b)/2);

ある種の分割統治アルゴリズムを使用しようとしているようです。b - aが 1 に等しいときに停止するだけで十分です。

于 2015-06-01T01:39:54.110 に答える