不等式の拘束剤を満たす{x、y}の離散ペアを見つけます
質問
私はいくつかの不平等を持っています {x,y}
, 、それは次の方程式を満たします:
x>=0
y>=0
f(x,y)=x^2+y^2>=100
g(x,y)=x^2+y^2<=200
ご了承ください x
と y
整数でなければなりません。
グラフィカルには、次のように表現できます。青い領域は、上記の不平等を満たす領域です。
問題は、Matlabにすべての許容できるペアを見つける機能があるということです。 {x,y}
?この種のことをするアルゴリズムがある場合、私もそれについて聞いてうれしいです。
もちろん、私たちが常に使用できるアプローチの1つは、可能なすべての組み合わせをテストするブルートフォースアプローチです {x,y}
不平等が満たされているかどうかを確認します。しかし、これは最後の手段です。なぜなら、それは時間がかかるからです。これを行う巧妙なアルゴリズム、または最良の場合、すぐに使用できる既存のライブラリを探しています。
x^2+y^2>=100
and x^2+y^2<=200
単なる例です。現実に f
と g
任意の程度の多項式機能になることがあります。
編集:C#コードも歓迎されます。
解決
これは、有限数のソリューションがあっても、列挙検索以外の方法で、一般的に多項式の不平等の一般的なセットに対して行うことはできません。 (おそらく、些細なことではないと言うべきです。列挙検索は、浮動ポイントの問題を条件として機能します。)関心のある領域は、高次の不平等のために単純に接続する必要はないことに注意してください。
編集:OPは、どのように検索を行うかについて尋ねました。
問題を考慮してください
x^3 + y^3 >= 1e12
x^4 + y^4 <= 1e16
x >= 0, y >= 0
このシステムのすべての整数ソリューションを解決します。すべての整数ソリューションが要求されているため、ここではいかなる形式でも整数プログラミングでは十分ではないことに注意してください。
ここでメッシュグリッドを使用すると、ドメイン(0:10000)x(0:10000)のポイントを見ることができます。そのため、1E8ポイントのセットをサンプリングすることを強制し、すべてのポイントをテストして、制約を満たすかどうかを確認します。
単純なループは、それよりも効率的になる可能性がありますが、それでもある程度の努力が必要です。
% Note that I will store these in a cell array,
% since I cannot preallocate the results.
tic
xmax = 10000;
xy = cell(1,xmax);
for x = 0:xmax
% solve for y, given x. This requires us to
% solve for those values of y such that
% y^3 >= 1e12 - x.^3
% y^4 <= 1e16 - x.^4
% These are simple expressions to solve for.
y = ceil((1e12 - x.^3).^(1/3)):floor((1e16 - x.^4).^0.25);
n = numel(y);
if n > 0
xy{x+1} = [repmat(x,1,n);y];
end
end
% flatten the cell array
xy = cell2mat(xy);
toc
必要な時間は...
Elapsed time is 0.600419 seconds.
私たちがテストしたかもしれない100020001の組み合わせのうち、いくつのソリューションを見つけましたか?
size(xy)
ans =
2 4371264
確かに、徹底的な検索の書き込みが簡単です。
tic
[x,y] = meshgrid(0:10000);
k = (x.^3 + y.^3 >= 1e12) & (x.^4 + y.^4 <= 1e16);
xy = [x(k),y(k)];
toc
8ギグのRAMで、64ビットマシンでこれを実行しました。しかし、それでもテスト自体はCPU HOGでした。
Elapsed time is 50.182385 seconds.
フローティングポイントの考慮事項は、計算の完了方法に応じて、異なる数のポイントを見つけることがあることに注意してください。
最後に、制約方程式がより複雑な場合、yの境界の式に根を使用して、制約が満たされる場所を特定する必要があるかもしれません。ここでの良いことは、それがまだより複雑な多項式の境界のために機能することです。