5

我有一些关于 的不等式{x,y}，它们满足以下等式：

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

请注意，xandy必须是整数。

图形上可以表示如下，蓝色区域是满足上述不等式的区域：

替代文字

现在的问题是，Matlab 中是否有任何函数可以找到每个可接受的对{x,y}？如果有一种算法可以做这种事情，我也会很高兴听到它。

当然，我们总是可以使用的一种方法是蛮力方法，我们测试每种可能的组合{x,y}以查看是否满足不等式。但这是最后的手段，因为它很耗时。我正在寻找一个聪明的算法来做到这一点，或者在最好的情况下，一个我可以立即使用的现有库。

这些x^2+y^2>=100 and x^2+y^2<=200只是示例；实际上f，并且g可以是任何次数的任何多项式函数。

编辑：也欢迎 C# 代码。

4

1 回答 1

4

对于一般的多项式不等式，这肯定是不可能通过除枚举搜索之外的任何方法来完成的，即使存在有限数量的解。（也许我应该说不是微不足道的，因为它是可能的。枚举搜索会起作用，但会受到浮点问题的影响。）请注意，对于高阶不等式，不需要简单地连接感兴趣的域。

编辑： 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

我在 64 位机器上运行这个，有 8 gig 的 ram。但即便如此，测试本身还是占用了 CPU。

Elapsed time is 50.182385 seconds.

请注意，浮点考虑有时会导致找到不同数量的点，具体取决于计算的完成方式。

最后，如果您的约束方程更复杂，您可能需要在表达式中使用 y 边界的根，以帮助确定满足约束的位置。这里的好处是它仍然适用于更复杂的多项式边界。

于 2010-09-25T13:51:40.580 回答