6

我有一个问题可以归结为找到一种将三角矩阵映射到跳过对角线的向量的方法。

基本上我需要使用 Gecode 库翻译这个 C++ 代码

// implied constraints
for (int k=0, i=0; i<n-1; i++)
  for (int j=i+1; j<n; j++, k++)
    rel(*this, d[k], IRT_GQ, (j-i)*(j-i+1)/2);

进入这个 MiniZinc(功能)代码

constraint
   forall ( i in 1..m-1 , j in i+1..m )
      ( (differences[?]) >=  (floor(int2float(( j-i )*( j-i+1 )) / int2float(2)) ));

我需要找出differences[?].

MiniZinc 是一种没有适当 for 循环的函数/数学语言。所以我必须映射那些索引 i 和 j ,它们接触上三角矩阵的所有且仅接触上三角矩阵的单元格,跳过其对角线,将这些单元格从 0 编号到任何值。

如果这是一个规则的三角矩阵(不是),这样的解决方案可以了

index = x + (y+1)*y/2

我正在处理的矩阵是一个n*n索引从 0 到 n-1 的方阵,但是为矩阵提供更通用的解决方案会很好n*m

这是完整的 Minizinc 代码

% modified version of the file found at https://github.com/MiniZinc/minizinc-benchmarks/blob/master/golomb/golomb.mzn

include "alldifferent.mzn";

int: m;
int: n = m*m;
array[1..m] of var 0..n: mark;
array[int] of var 0..n: differences = [mark[j] - mark[i] | i in 1..m, j in i+1..m];

constraint mark[1] = 0;

constraint forall ( i in 1..m-1 ) ( mark[i] < mark[i+1] );

% this version of the constraint works
constraint forall ( i in 1..m-1 , j in i+1..m )
    ( (mark[j] - mark[i]) >= (floor(int2float(( j-i )*( j-i+1 )) / int2float(2))) );

%this version does not
%constraint forall ( i in 1..m-1, j in i+1..m )
%    ( (differences[(i-1) + ((j-2)*(j-1)) div 2]) >= (floor(int2float(( j-i )*( j-i+1 )) / int2float(2))) );

constraint alldifferent(differences);

constraint differences[1] < differences[(m*(m-1)) div 2];

solve :: int_search(mark, input_order, indomain, complete) minimize mark[m];

output ["golomb ", show(mark), "\n"];

谢谢。

4

1 回答 1

5

当心。您从该链接中找到的公式index = x + (y+1)*y/2包括对角线条目,并且适用于下三角矩阵,我认为这不是您想要的。您正在寻找的确切公式实际上是index = x + ((y-1)y)/2 (参见:https ://math.stackexchange.com/questions/646117/how-to-find-a-function-mapping-matrix-indices )。

再次注意,我给你的这个公式假设你的索引:x,y,是从零开始的您的 MiniZinc 代码使用从 1 (1 <= i <= m), 1 <= j <= m))开始的索引 i,j 。对于从 1 开始的索引,公式为T(i,j) = i + ((j-2)(j-1))/2。所以你的代码应该是这样的:

constraint
   forall ( i in 1..m-1 , j in i+1..m )
      ((distances[(i + ((j-2)*(j-1)) div 2]) >= ...

请注意,它(j-2)(j-1)始终是 2 的倍数,因此我们可以只使用除数为 2 的整数除法(无需担心转换为浮点数/从浮点数转换)。


以上假设您使用的是m*m方阵。
为了推广到M*N矩形矩阵,一个公式可以是:

通式

其中 0 <= i < M, 0<= j < N [如果您再次需要索引从 1 开始,将 i 替换为 i-1 并将 j 替换为 j-1 在上面的公式中]。这会触及上三角矩阵的所有单元格以及当 N > M 时出现的正方形的“边上的额外块”。也就是说,它触及所有单元格 (i,j),使得 i < j 为 0 <= i < M, 0 <= j < N。

侧面的额外块


完整代码:

% original: https://github.com/MiniZinc/minizinc-benchmarks/blob/master/golomb/golomb.mzn

include "alldifferent.mzn";

int: m;
int: n = m*m;
array[1..m] of var 0..n: mark;
array[1..(m*(m-1)) div 2] of var 0..n: differences;

constraint mark[1] = 0;
constraint forall ( i in 1..m-1 ) ( mark[i] < mark[i+1] );
constraint alldifferent(differences);
constraint forall (i,j in 1..m where j > i)
    (differences[i + ((j-1)*(j-2)) div 2] = mark[j] - mark[i]);
constraint forall (i,j in 1..m where j > i)
    (differences[i + ((j-1)*(j-2)) div 2] >= (floor(int2float(( j-i )*( j-i+1 )) / int2float(2))));
constraint differences[1] < differences[(m*(m-1)) div 2];

solve :: int_search(mark, input_order, indomain, complete)
    minimize mark[m];

output ["golomb ", show(mark), "\n"];

下三角版本(取之前的代码并在必要时交换 i 和 j):

% original: https://github.com/MiniZinc/minizinc-benchmarks/blob/master/golomb/golomb.mzn

include "alldifferent.mzn";

int: m;
int: n = m*m;
array[1..m] of var 0..n: mark;
array[1..(m*(m-1)) div 2] of var 0..n: differences;

constraint mark[1] = 0;
constraint forall ( i in 1..m-1 ) ( mark[i] < mark[i+1] );
constraint alldifferent(differences);
constraint forall (i,j in 1..m where i > j)
    (differences[j + ((i-1)*(i-2)) div 2] = mark[i] - mark[j]);
constraint forall (i,j in 1..m where i > j)
    (differences[j + ((i-1)*(i-2)) div 2] >= (floor(int2float(( i-j )*( i-j+1 )) / int2float(2))));
constraint differences[1] < differences[(m*(m-1)) div 2];

solve :: int_search(mark, input_order, indomain, complete)
    minimize mark[m];

output ["golomb ", show(mark), "\n"];
于 2014-12-05T23:04:49.180 回答