Finden Sie die Discrete Paar {x, y}, dass Satisfy Ungleichungsrestriktionen

Question

Ich habe einige Ungleichheiten in Bezug auf {x,y}, dass die folgenden Gleichungen erfüllt:

x>=0
y>=0
f(x,y)=x^2+y^2>=100
g(x,y)=x^2+y^2<=200

Beachten Sie, dass x und y integer sein müssen.

Grafisch kann wie folgt dargestellt werden, die blaue Region ist die Region, die die obigen Ungleichungen erfüllt:

Die Frage ist nun, gibt es eine Funktion in Matlab, die jedes zulässiges Paar {x,y} findet? Wenn ein Algorithmus ist diese Art der Sache zu tun, würde ich mich freuen, um es so gut zu hören.

Natürlich wird ein Ansatz können wir immer den Einsatz Brute-Force-Ansatz, wo wir jede mögliche Kombination von {x,y} testen, um zu sehen, ob die Ungleichheiten erfüllt sind. Aber dies ist der letzte Ausweg, weil es zeitaufwendig. Ich bin auf der Suche nach einem cleveren Algorithmus, der dies tut, oder im besten Fall eine bestehende Bibliothek, dass ich gerade weg verwenden kann.

Die x^2+y^2>=100 and x^2+y^2<=200 sind nur Beispiele; in Wirklichkeit f und g können beliebige Polynomfunktionen von jedem Grad sein.

Edit:. C # -Code sind ebenfalls willkommen

Answer 1

Dies ist sicherlich nicht möglich, für einen allgemeinen Satz von Polynom Ungleichheiten im Allgemeinen zu tun, durch eine andere Methode als enumerative suchen, auch wenn es eine endliche Anzahl von Lösungen. (Vielleicht sollte ich nicht trivial sagen, wie es möglich ist. Enumerative Suche funktionieren wird, vorbehaltlich Gleitkomma Fragen.) Beachten Sie, dass die Domäne von Interesse muss nicht nur für eine höhere Ordnung Ungleichheiten verbunden werden.

Edit:. Die OP hat gefragt, wie man vorgehen könnte eine Suche tun

Betrachten Sie das Problem

x^3 + y^3 >= 1e12
x^4 + y^4 <= 1e16

x >= 0, y >= 0

Lösen Sie für alle ganzzahligen Lösungen dieses Systems. Beachten Sie, dass Integer-Programmierung in irgendeiner Form hier nicht ausreichen wird, da alle ganzzahligen Lösungen angefordert werden.

Die Verwendung von meshgrid hier würde uns zwingen, an den Punkten in der Domäne zu erhalten (0: 10000) X (0: 10000). So würde es uns zwingen, eine Reihe von 1e8 Punkten abzutasten, jeden Punkt zu testen, um zu sehen, ob sie die Bedingungen erfüllen.

Eine einfache Schleife kann möglicherweise effizienter sein als das, obwohl es noch einige Anstrengungen erfordern.

% 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

Die Zeit, die erforderlich war ...

Elapsed time is 0.600419 seconds.

Von den 100020001 Kombinationen dass wir getestet haben könnte, wie viele Lösungen haben wir finden?

size(xy)
ans =
           2     4371264

Zugegeben, die erschöpfende Suche ist einfacher zu schreiben.

tic
[x,y] = meshgrid(0:10000);
k = (x.^3 + y.^3 >= 1e12) & (x.^4 + y.^4 <= 1e16);
xy = [x(k),y(k)];
toc

Ich lief dies auf einer 64-Bit-Maschine, mit 8 GB RAM. Aber auch so der Test selbst war eine CPU hog.

Elapsed time is 50.182385 seconds.

Beachten Sie, dass Gleitkomma-Überlegungen werden manchmal eine unterschiedliche Anzahl von Punkten verursachen gefunden werden, je nachdem, wie die Berechnungen werden durchgeführt.

Schließlich, wenn Ihre Bedingungsgleichungen komplexer sind, müssen Sie möglicherweise Wurzeln für die Grenzen auf y im Ausdruck verwenden, um Hilfe zu identifizieren, in denen die Bedingungen erfüllt sind. Die nette Sache hier ist es immer noch mehr kompliziert Polynom Grenzen arbeitet.