2016-05-01 6 views
1

私はfzero関数を使って1つのパラメータ に応じて非線形方程式を解いていますが、私の方法には満足できません。私はこれらの問題があります:MatlabやOctaveでfzeroを使うとループや複雑な解を避ける

1)パラメータのforループは避けられますか?

2)複雑な解を避けるために、まずfzeroの有効区間を事前に計算する必要があります。 ここにはより良い解決策がありますか?

パラメータのステップサイズを小さくすると、実行時間が遅くなります。 の間隔をあらかじめ計算していないと、「区間のエンドポイントの関数値は有限で実数でなければならない」というエラーが表示されます。 Matlabでは 、Octaveでは "fzero:有効ではない初期ブラケット"となります。

ここコードMATLABに

% solve y = 90-asind(n*(sind(90-asind(sind(a0)/n)-y))) 

% set the equation paramaters 
n=1.48; a0=0:0.1:60; 

% loop over a0 
for i = 1:size(a0,2) 

    % for each a0 find where the argument of outer asind() 
    % will not give complex solutions, i.e. argument is between 1 and -1 

    fun1 = @(y) n*(sind(90-asind(sind(a0(i))/n)-y))-0.999; 
    y1 = fzero(fun1,[0 90]); 
    fun2 = @(y) n*(sind(90-asind(sind(a0(i))/n)-y))+0.999; 
    y2 = fzero(fun2,[0 90]); 


    % use y1, y2 as limits in fzero interval 

    fun3 = @(y) 90-asind(n*(sind(90-asind(sind(a0(i))/n)-y)))-y; 
    y(i) = fzero(fun3, [y1 y2]); 

end 

% plot the result 
figure; plot(y); grid minor; 
xlabel('Incident ray angle [Deg]'); 
ylabel('Lens surface tangent angle'); 

答えて

1

あり、Iは、以下の単純化されたループで以下のプロットを得ました。

for i = 1:size(a0,2) 
    fun3 = @(y) sind(90-y) - n*(sind(90-asind(sind(a0(i))/n)-y)); 
    y(i) = fzero(fun3, [0,90]); 
end 

違いは式の形です:私は90-y = asind(something)をsin(90-y)= somethingに置き換えました。 「何か」が絶対値で1より大きい場合、前者のバージョンはasindの複素数値のためにエラーを投げます。後者は、これが解ではないことを認識して、通常通り進行する(sin(90-y)は1より大きい何かに等しくない)。

ブラケットの事前計算は必要ありませんでした。[0,90]は単純に機能しました。私が作ったもう一つの変更は、正しい水平軸を得るために、plot(y)の代わりにplot(a0,y)のところにありました。

そしてここではforループを避けることはできません。また、心配する必要はありません。ベクトル化とは、コンテンツが低レベルの操作であるループを排除することを意味し、一部のC配列で操作することによってen masseを実行することができます。しかし、fzeroはまったくそうではありません。コードの実行に時間がかかる場合は、一連の方程式を解くのに時間がかかります。なぜなら、forというループがあるからです。

plot

+0

おかげで、コードは単純かつ明確に見えます。あなたのソリューションはOctave 4.0でも動作しますが、エラーや警告はありません。 – miquo

関連する問題