How to generate a random matrix with binary entries satisfying certain conditions

299 views
Skip to first unread message

ghjklmn

unread,
Jun 28, 2020, 7:54:51 AM6/28/20
to AMPL Modeling Language
I am trying to generate a random matrix A = {X, Y} such that each entry of A is either 0 or 1, and no two entries of A in the same column and two successive rows, with the first row is always odd-numbered (e.g. entries in rows 1 & 2, rows 3&4, rows 5&6, etc), are both equal to 1 (they can be both 0 though). Below is my attempt to do this:

param n;
set X = 1..n ;
set Y = 1....(2*n+1) ;

#generate a random matrix with entries = 0 or 1
param A {X, Y} := ceil(Uniform(0,2))-1

for (j in Y, i in X} {
   if (ord(i) mod 2 == 1) then 
            if (A[i, j] == A[i+1, j] and A[i, j] == 1) then let (A[i, j] = 0 or A[i+1, j]=0); 
   } ;

display A;

Question. Is this code correct to achieve my goal above (especially the nested IF lines within the `for` loop?) I tried putting this into a data file and set param n := 5 : but I keep getting the error message as in the attached file. Can someone please help me overcome this issue, and/or help verify the code above?




Screen Shot 2020-06-28 at 12.34.03 AM.png

AMPL Google Group

unread,
Jun 28, 2020, 12:44:05 PM6/28/20
to AMPL Modeling Language
I already commented on an earlier version of this question. In this version you can't write "let (A[i, j] = 0 or A[i+1, j]=0" since that does not tell AMPL which assignment to make. (Also, AMPL's "let" statements use the := operator rather than the = operator.) Instead you could randomly assign one or the other element to zero:

if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
else let A[i+1, j] := 0;

You need to correct a lot of simple syntax errors as well, such as writing { instead of ( to begin a set expression, adding a semicolon at the end of every statement, etc.

If you have any other questions, respond only to this answer.


--
Robert Fourer
am...@googlegroups.com
{#HS:1207695670-81605#}
--
You received this message because you are subscribed to the Google Groups "AMPL Modeling Language" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ampl+uns...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ampl/bc82a815-4acc-4cd8-a869-05485f61b5e8o%40googlegroups.com.

ghjklmn

unread,
Jun 28, 2020, 6:58:17 PM6/28/20
to AMPL Modeling Language
Thank you very much for your assistance on that last part. I have a further question: if I only want to check if the two successive rows in certain columns of a param A = {X,Y} are not both 1 (for example, no two 1s in rows 1 and 2 of column 1, no two 1s in rows 3 and 4 of column 2, etc), then is the logic of the following nested IF correct? Do I need parentheses like ( ) or { } around some of the outer IF...then?

param n;
set X := 1..n ;
set Y := 1....(2*n+1) ;

#generate a random matrix with entries = 0 or 1
param A {X, Y} := floor(Uniform(0,2)) ;

for {j in Y, i in X}{ 
   if (ord(i) mod 2 = 1 and ord(j) = ceil(ord(i)/2) ) then 
            if (A[i, j] = A[i+1, j] and A[i, j] = 1) then 
                if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
                else let A[i+1, j] := 0; 
   } ;

To unsubscribe from this group and stop receiving emails from it, send an email to am...@googlegroups.com.

AMPL Google Group

unread,
Jun 29, 2020, 8:43:45 AM6/29/20
to AMPL Modeling Language
You can write the loop like this,

for {j in Y, i in X}{
   if (i mod 2 = 1 and j = ceil(i/2) ) then
      if (A[i, j] = A[i+1, j] and A[i, j] = 1) then {
         if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
         else let A[i+1, j] := 0;
      };
} ;

or like this,

for {j in Y, i in X} {
   if (i mod 2 = 1 and j = ceil(i/2) ) then
      if (A[i, j] = A[i+1, j] and A[i, j] = 1) then
         if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
         else let A[i+1, j] := 0;
      else;
}


--
Robert Fourer
am...@googlegroups.com
{#HS:1207695670-81605#}
On Sun, Jun 28, 2020 at 10:58 PM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
Thank you very much for your assistance on that last part. I have a further question: if I only want to check if the two successive rows in certain columns of a param A = {X,Y} are not both 1 (for example, no two 1s in rows 1 and 2 of column 1, no two 1s in rows 3 and 4 of column 2, etc), then is the logic of the following nested IF correct? Do I need parentheses like ( ) or { } around some of the outer IF...then?

param n;
set X := 1..n ;
set Y := 1....(2*n+1) ;

#generate a random matrix with entries = 0 or 1
param A {X, Y} := floor(Uniform(0,2)) ;

for {j in Y, i in X}{
if (ord(i) mod 2 = 1 and ord(j) = ceil(ord(i)/2) ) then
if (A[i, j] = A[i+1, j] and A[i, j] = 1) then
if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
else let A[i+1, j] := 0;
} ;

On Sun, Jun 28, 2020 at 4:43 PM UTC, AMPL Google Group <am...@googlegroups.com> wrote:
I already commented on an earlier version of this question. In this version you can't write "let (A[i, j] = 0 or A[i+1, j]=0" since that does not tell AMPL which assignment to make. (Also, AMPL's "let" statements use the := operator rather than the = operator.) Instead you could randomly assign one or the other element to zero:

if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
else let A[i+1, j] := 0;

You need to correct a lot of simple syntax errors as well, such as writing { instead of ( to begin a set expression, adding a semicolon at the end of every statement, etc.

If you have any other questions, respond only to this answer.


--
Robert Fourer
am...@googlegroups.com
On Sun, Jun 28, 2020 at 11:54 AM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
I am trying to generate a random matrix A = {X, Y} such that each entry of A is either 0 or 1, and no two entries of A in the same column and two successive rows, with the first row is always odd-numbered (e.g. entries in rows 1 & 2, rows 3&4, rows 5&6, etc), are both equal to 1 (they can be both 0 though). Below is my attempt to do this:

param n;
set X = 1..n ;
set Y = 1....(2*n+1) ;

#generate a random matrix with entries = 0 or 1
param A {X, Y} := ceil(Uniform(0,2))-1

for (j in Y, i in X} {
if (ord(i) mod 2 == 1) then
if (A[i, j] == A[i+1, j] and A[i, j] == 1) then let (A[i, j] = 0 or A[i+1, j]=0);
} ;

display A;

Question. Is this code correct to achieve my goal above (especially the nested IF lines within the `for` loop?) I tried putting this into a data file and set param n := 5 : but I keep getting the error message as in the attached file. Can someone please help me overcome this issue, and/or help verify the code above?




--
You received this message because you are subscribed to the Google Groups "AMPL Modeling Language" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ampl+uns...@googlegroups.com.

ghjklmn

unread,
Jun 29, 2020, 9:09:08 PM6/29/20
to AMPL Modeling Language
Thank you very much for your help. I have a follow-up question though: keeping the same constraints, if I need to constrain on the number of 1s appearing in matrix A (for example, number of 1 to number of 0 is 1:9 for a given dim(A) = dim(Y) x dim(X) , how can I achieve that? The only possible way I can think of is to build a count variable to count the total number of 1 after a matrix A is randomly generated, and then evaluate the ratio of that to dim(A) = (2n+1)*n to see if the ratio is 0.1 (most likely not). If not, then change some 0s to 1 until the ratio is equal to 0.1. But this seems to be way too difficult to achieve.

For the first step (after A is generated), I can have like:

param rat ;
param fixed_rat := 1/9 ;  #pre-defined ratio of 1 to 0 

let count := 0 ;

for {i in X, j in Y}  {
    if(A[i, j] = 1) then count = count + 1;
}

rat = count/(card(X)*card(Y)) ;

if (rat <> fixed_rat) then {how do to this?} 
To unsubscribe from this group and stop receiving emails from it, send an email to am...@googlegroups.com.

ghjklmn

unread,
Jun 30, 2020, 12:02:49 AM6/30/20
to AMPL Modeling Language
So I ran into a problem that is so weird logically. First, part of my model file is as follows:

set X ordered ;

set DONORS ordered ;


set Z :=  X ordered ;


param ARCS {DONORS, X}; #denote an arc

param n; 


Now, I have the data file as:


param n := 100000;

set X := 1..n ; 

set DONORS := 1..2*n ; 


param ARCS := floor(Uniform(0,2)); 


for {j in X, i in DONORS} {

if (ord(i) mod 2 = 1 and ord(j) = ceil(ord(i)/2)) then

if (ARCS[i, j] = ARCS[i+1, j] and ARCS[i, j] = 1) then {

if Uniform(0,1) <= 0.5 then let ARCS[i, j] := 0;

                        else let ARCS[i+1, j] := 0;

} ;

} ;


i kept getting the error in the images, but I don't understand why they were errors. Like, how should we *correctly* declare a two-dimensional table array with random entries of 0 and 1?

To unsubscribe from this group and stop receiving emails from it, send an email to am...@googlegroups.com.
Screen Shot 2020-06-29 at 10.47.31 PM.png
Screen Shot 2020-06-29 at 10.44.07 PM.png

AMPL Google Group

unread,
Jun 30, 2020, 12:17:45 PM6/30/20
to AMPL Modeling Language
AMPL expressions are not accepted in data files. To define a set or parameter by an expression, give the expression in the set's definition in the model file instead:

param n;
set X = 1..n ordered;
set DONORS = 1..2*n ordered;
param ARCS {DONORS, X} = floor(Uniform(0,2));

Then you only need to give "param n := 100000;" in the data file. In your situation, however, this is not going to work for ARCS, because you want to change some its values; when you try to execute a statement like "let ARCS[i, j] := 0;" you will get an error like this:

ARCS has an = assignment in the model.

So instead your model should define only

param ARCS {DONORS, X};

and then later you should use "let" to initialize the ARCS values:

let {j in X, i in DONORS} ARCS[i,j] := floor(Uniform(0,2));
for {j in X, i in DONORS} { . . . change some ARCS values . . . }

I put all of this together in the attached file.


--
Robert Fourer
am...@googlegroups.com
{#HS:1207695670-81605#}
On Tue, Jun 30, 2020 at 4:03 AM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
So I ran into a problem that is so weird logically. First, part of my model file is as follows:

set X ordered ;
set DONORS ordered ;

set Z := X ordered ;

param ARCS {DONORS, X}; #denote an arc
param n;

Now, I have the data file as:

param n := 100000;
set X := 1..n ;
set DONORS := 1..2*n ;

param ARCS := floor(Uniform(0,2));

for {j in X, i in DONORS} {
if (ord(i) mod 2 = 1 and ord(j) = ceil(ord(i)/2)) then
if (ARCS[i, j] = ARCS[i+1, j] and ARCS[i, j] = 1) then {
if Uniform(0,1) <= 0.5 then let ARCS[i, j] := 0;
else let ARCS[i+1, j] := 0;
} ;
} ;

i kept getting the error in the images, but I don't understand why they were errors. Like, how should we *correctly* declare a two-dimensional table array with random entries of 0 and 1?

On Tue, Jun 30, 2020 at 1:09 AM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
Thank you very much for your help. I have a follow-up question though: keeping the same constraints, if I need to constrain on the number of 1s appearing in matrix A (for example, number of 1 to number of 0 is 1:9 for a given dim(A) = dim(Y) x dim(X) , how can I achieve that? The only possible way I can think of is to build a count variable to count the total number of 1 after a matrix A is randomly generated, and then evaluate the ratio of that to dim(A) = (2n+1)*n to see if the ratio is 0.1 (most likely not). If not, then change some 0s to 1 until the ratio is equal to 0.1. But this seems to be way too difficult to achieve.

For the first step (after A is generated), I can have like:

param rat ;
param fixed_rat := 1/9 ; #pre-defined ratio of 1 to 0

let count := 0 ;

for {i in X, j in Y} {
if(A[i, j] = 1) then count = count + 1;
}

rat = count/(card(X)*card(Y)) ;

if (rat <> fixed_rat) then {how do to this?}

On Mon, Jun 29, 2020 at 12:43 PM UTC, AMPL Google Group <am...@googlegroups.com> wrote:
You can write the loop like this,

for {j in Y, i in X}{
   if (i mod 2 = 1 and j = ceil(i/2) ) then
      if (A[i, j] = A[i+1, j] and A[i, j] = 1) then {
         if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
         else let A[i+1, j] := 0;
      };
} ;

or like this,

for {j in Y, i in X} {
   if (i mod 2 = 1 and j = ceil(i/2) ) then
      if (A[i, j] = A[i+1, j] and A[i, j] = 1) then
         if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
         else let A[i+1, j] := 0;
      else;
}


--
Robert Fourer
am...@googlegroups.com
On Sun, Jun 28, 2020 at 10:58 PM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
Thank you very much for your assistance on that last part. I have a further question: if I only want to check if the two successive rows in certain columns of a param A = {X,Y} are not both 1 (for example, no two 1s in rows 1 and 2 of column 1, no two 1s in rows 3 and 4 of column 2, etc), then is the logic of the following nested IF correct? Do I need parentheses like ( ) or { } around some of the outer IF...then?

param n;
set X := 1..n ;
set Y := 1....(2*n+1) ;

#generate a random matrix with entries = 0 or 1
param A {X, Y} := floor(Uniform(0,2)) ;

for {j in Y, i in X}{
if (ord(i) mod 2 = 1 and ord(j) = ceil(ord(i)/2) ) then
if (A[i, j] = A[i+1, j] and A[i, j] = 1) then

if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
else let A[i+1, j] := 0;
} ;

On Sun, Jun 28, 2020 at 4:43 PM UTC, AMPL Google Group <am...@googlegroups.com> wrote:
I already commented on an earlier version of this question. In this version you can't write "let (A[i, j] = 0 or A[i+1, j]=0" since that does not tell AMPL which assignment to make. (Also, AMPL's "let" statements use the := operator rather than the = operator.) Instead you could randomly assign one or the other element to zero:

if Uniform(0,1) <= 0.5 then let A[i, j] := 0;
else let A[i+1, j] := 0;

You need to correct a lot of simple syntax errors as well, such as writing { instead of ( to begin a set expression, adding a semicolon at the end of every statement, etc.

If you have any other questions, respond only to this answer.


--
Robert Fourer
am...@googlegroups.com
On Sun, Jun 28, 2020 at 11:54 AM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
I am trying to generate a random matrix A = {X, Y} such that each entry of A is either 0 or 1, and no two entries of A in the same column and two successive rows, with the first row is always odd-numbered (e.g. entries in rows 1 & 2, rows 3&4, rows 5&6, etc), are both equal to 1 (they can be both 0 though). Below is my attempt to do this:

param n;
set X = 1..n ;
set Y = 1....(2*n+1) ;

#generate a random matrix with entries = 0 or 1
param A {X, Y} := ceil(Uniform(0,2))-1

for (j in Y, i in X} {
if (ord(i) mod 2 == 1) then
if (A[i, j] == A[i+1, j] and A[i, j] == 1) then let (A[i, j] = 0 or A[i+1, j]=0);
} ;

display A;

Question. Is this code correct to achieve my goal above (especially the nested IF lines within the `for` loop?) I tried putting this into a data file and set param n := 5 : but I keep getting the error message as in the attached file. Can someone please help me overcome this issue, and/or help verify the code above?




--
You received this message because you are subscribed to the Google Groups "AMPL Modeling Language" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ampl+uns...@googlegroups.com.
test2.run
Message has been deleted

ghjklmn

unread,
Jul 1, 2020, 9:43:49 PM7/1/20
to AMPL Modeling Language
Thanks to your great explanation and example file, I was finally able to make it run and solve the optimization problem. 

I have two additional questions related to the param ARCS and n in the model file:

1. How to adjust the resulting randomly-generated param ARCS with one additional constraint: the # of 1 to # of 0, in the randomly-generated matrix, satisfies a pre-defined ratio (for example, (# of 1s): (# of 0s) = 1:10) ?? I posted the question and my best attempt already, but I think you were too busy that you overlooked it. Note that this would be placed at the bottom of your Test2.run above.

let fixed_rat := 1/9 ; #pre-defined ratio 
let count := 0 ;
for {j in X, i in DONORS} {
      if(ARCS[i, j] = 1) then count = count + 1;
}

let rat := count/(card(X)*card(Y)) ;
if (rat <> fixed_rat) then {how do I continue here??}

2. I would like to conduct a simulation on the objective values and cplex solver's time for different values of n (for example, run the current code for n from 500 to 4500 with an incremental step size = 500), how can we do this without changing the entire code in the model file? It seems to me I have to declare set n = 500..T by 500 in model file, but this does not work because n is used in the definition of set X and Y. Rewriting the entire definition of X and Y in terms of the indices of fails to work either.


To unsubscribe from this group and stop receiving emails from it, send an email to am...@googlegroups.com.

AMPL Google Group

unread,
Jul 2, 2020, 9:21:11 AM7/2/20
to AMPL Modeling Language
1. If there are too many 1s, pick an element of ARCS at random, and if it is a 1 that is eligible to be a 0, change it to 0. Repeat until the number of 1s is reduced to the desired number. If there are not enough ones, use the same approach but change 0s to 1s. An AMPL "while" loop can be useful for this purpose; the general outline for the "too many 1s" case would be

while (number of 1s in ARCS is > desired number)
choose random i and j
if ARCS[i,j] = 1 and it is allowed to be zero
then let ARCS[i,j] := 0;
}

2. Increasing n will automatically increase the size of your problem. So instead of setting n in a data file, you could consider a loop like this:

for {size in 500..4500 by 500} {
let n := size;
generate a new ARCS matrix
solve;
record results
}


--
Robert Fourer
am...@googlegroups.com
{#HS:1207695670-81605#}
On Thu, Jul 2, 2020 at 1:44 AM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
Thanks to your great explanation and example file, I was finally able to make it run and solve the optimization problem.

I have two additional questions related to the param ARCS and n in the model file:

1. How to adjust the resulting randomly-generated param ARCS with one additional constraint: the # of 1 to # of 0, in the randomly-generated matrix, satisfies a pre-defined ratio (for example, (# of 1s): (# of 0s) = 1:10) ?? I posted the question and my best attempt already, but I think you were too busy that you overlooked it. Note that this would be placed at the bottom of your Test2.run above.

let fixed_rat := 1/9 ; #pre-defined ratio
let count := 0 ;
for {j in X, i in DONORS} {
if(ARCS[i, j] = 1) then count = count + 1;
}

let rat := count/(card(X)*card(Y)) ;
if (rat <> fixed_rat) then {how do I continue here??}

2. I would like to conduct a simulation on the objective values and cplex solver's time for different values of n (for example, run the current code for n from 500 to 4500 with an incremental step size = 500), how can we do this without changing the entire code in the model file? It seems to me I have to declare set n = 500..T by 500 in model file, but this does not work because n is used in the definition of set X and Y. Rewriting the entire definition of X and Y in terms of the indices of fails to work either.

On Wed, Jul 1, 2020 at 5:06 AM UTC, AMPL Modeling Language <am...@googlegroups.com> wrote:
Thanks to your great explanation and example file, I was finally able to make it run and solve the IP!! I cannot believe choosing between the model file and the run files to put the correct syntax, all for the declaration and value initialization purposes, into matters that much in AMPL (?!!) The result obtained is quite surprising (and weird) to me though, but I sincerely appreciate your help on this, Dr. Fourer. Below is what I did in the run file based on your hint:

model chain.mod;

data sample2.dat;

option randseed 0;


let {j in X, i in DONORS} ARCS[i,j] := floor(Uniform(0,2));

for {j in X, i in DONORS} {

if (ord(i) mod 2 = 1 and ord(j) = ceil(ord(i)/2)) then

if (ARCS[i, j] = ARCS[i+1, j] and ARCS[i, j] = 1) then {

if Uniform(0,1) <= 0.5 then let ARCS[i, j] := 0;

else let ARCS[i+1, j] := 0;

} ;

} ;

option solver cplex, solution_round 6;

solve;

option omit_zero_rows 1; #supress all zero-rows

option omit_zero_col 1; #supress all zero-columns

display _solve_time;

Finally, I have two additional questions related to the param ARCS and n in the data file:

1. How to adjust the resulting randomly-generated param ARCS with one additional constraint: the ratio of 1 to 0, in the resulting random matrix, must follow certain pre-defined ratio (e.g. # of 1s/# of 0s = 1:10) ?? I posted the question and my best attempt already, but I think you were too busy that you overlooked it. Note that this would be placed at the bottom of your Test2.run above.


let fixed_rat := 1/9 ; #pre-defined ratio
let count := 0 ;
for {j in X, i in DONORS} {
if(ARCS[i, j] = 1) then count = count + 1;
}

let rat := count/(card(X)*card(Y)) ;
if (rat <> fixed_rat) then {how do I continue here??}

2. I would like to conduct a Monte Carlo simulation on the objective values and cplex solver's time for different values of n (for example, run the current code for n from 500 to 4500 with an incremental step = 500), how can we do this without changing the entire code in the model file? It seems to me I have to declare set n = 500..T by 500 in model file, but this will mess up with the cardinality of set X and set Y in the current working code (because set X:= 1..n, so n should be an integer, rather than a set. But rewriting the entire definition of set X and Y in terms of the indices of n is way too complicated for me, and I am not sure if it works).

3. My current AMPL code runs out of memory usage when n:= 5000 (my laptop is MacAir 2014, and the AMPL's Cplex took 2.5 - 3 hours to solve when n = 4000. The obtained result is straightforward, and the same as other cases with smaller value of n). But I have not tried using all 12 threads available in MacAir (I am not sure how cplex solver solves my Integer Programming problem, and how many threads it uses to solve primal, dual and barrier, if any). Do you think if I specify 'option cplex_options 'concurrentopt,' this would solve the memory overload issue?

For your convenience and concrete reference, attached is my fully working code. I apologize for having so many questions.

On Tue, Jun 30, 2020 at 4:17 PM UTC, AMPL Google Group <am...@googlegroups.com> wrote:
AMPL expressions are not accepted in data files. To define a set or parameter by an expression, give the expression in the set's definition in the model file instead:

param n;
set X = 1..n ordered;
set DONORS = 1..2*n ordered;
param ARCS {DONORS, X} = floor(Uniform(0,2));

Then you only need to give "param n := 100000;" in the data file. In your situation, however, this is not going to work for ARCS, because you want to change some its values; when you try to execute a statement like "let ARCS[i, j] := 0;" you will get an error like this:

ARCS has an = assignment in the model.

So instead your model should define only

param ARCS {DONORS, X};

and then later you should use "let" to initialize the ARCS values:

let {j in X, i in DONORS} ARCS[i,j] := floor(Uniform(0,2));
for {j in X, i in DONORS} { . . . change some ARCS values . . . }

I put all of this together in the attached file.


--
Robert Fourer
am...@googlegroups.com
To unsubscribe from this group and stop receiving emails from it, send an email to ampl+uns...@googlegroups.com.

ghjklmn

unread,
Jul 6, 2020, 9:05:33 PM7/6/20
to AMPL Modeling Language
Dear Dr. Fourer,

Thank you very much for your help on both parts (they are very helpful to me!). However, for 1, to choose random i and j, can we write this as:

let i:= floor(Uniform(1, card(X)) 
let j:= floor(Uniform(1, card(DONORS)))
if (2*ord(i) - 1 = ord(j) and ARCS[i,j] = 1) or (2*ord(i) = ord(j) and ARCS[i, j]=1)
then let ARCS[i,j] := 0;

2. Attached is the model , data and run files that I made. Everything executed perfectly, except the constraint E in the mod.file is not satisfied at all, based on the solutions ?! I tried to use both the solvers cplex or ilogcp, but none of them gives the correct answer. Below is the constraint E:

subject to E {j in PATIENTS, k in PATIENTS, i in DONORS, m in DONORS: j <> k and i <> m and (2*ord(j)-1 = ord(i)) and (2*ord(k)-1 = ord(m)) and ARCS[i, j] > 0 and ARCS[m, k] > 0}:
y[i, j] + y[nextw(i,DONORS, 1), j] + y[m, k] + y[nextw(m, DONORS,1),k] + y[i, k] + y[nextw(i,DONORS, 1), k] 
  + y[m, j] + y[nextw(m, DONORS,1),j] + 1 <= (2*x[j] + 2*x[k]);

Result (by executing "include chain.run")

x [*] :=

p01  1

p02  1

p03  1

p04  0

p05  1


y [*,*]

:   p01 p02 p03 p04 p05    :=

d01   0   1   0   0   0

d02   0   1   0   0   0

d03   1   0   0   0   0

d04   1   0   0   0   0

d05   0   0   1   0   0

d06   0   0   0   0   1

d07   0   0   0   0   0

d08   0   0   0   0   0

d09   0   0   1   0   0

d10   0   0   0   0   1


E :=

p03 p04 d05 d07   0

p04 p03 d07 d05   0


As you might see, constraint E is not satisfied with the triples (p_01, d_01, d_02) and (p_02, d_03, d_04) (yellow color in solution y) or between (p_03, d_05, d_06) and (p_05, d_09, d_10) (blue color in solution y),  as the displayed E shows. I don't understand why the solver gives a supposedly infeasible solutions. Basically, what constraint E prevents the cases: 
If x_{p_01} = x_{p_02} = 1, 
then 
y_{d_01, p_01} + y_{d_02, p_01} + y_{d_01, p_02} +  y_{d_02, p_02} +  y_{d_03, p_02} + y_{d_03, p_01} + y_{d_04, p_01} + y_{d_04, p_02} + 1 <= 2(x_{p_01} + x_{p_02})  ------ Adding the slack constant 1 is because linear programming cannot deal with strict inequality constraint ("<" in this case), nor the cplex or ilogcp solver. But then if E is satisfied, then the solver cplex or ilogcp should not output the above results.

Question. Can you please help me figure out why this "violation of constraint E" occurs? I'm getting stuck on this weird issue for quite some time.

 
To unsubscribe from this group and stop receiving emails from it, send an email to am...@googlegroups.com.
sol.run
example.mod
sample.dat

AMPL Google Group

unread,
Jul 7, 2020, 1:16:34 PM7/7/20
to AMPL Modeling Language
1. You are only considering a change to ARCS[i,j] when ord(j) equals 2*ord(i) - 1 or 2*ord(i). But you could also consider changing other entries in order to adjust the number of 1s in ARCS. Anyhow, for whatever you write, you can try some tests to see whether the ARCS matrix is coming out correctly.

2. When the current solution does not satisfy a constraint, AMPL shows the "slack" on the constraint to be negative. Thus you can check whether CPLEX's solution satisfies the E constraints by using a display statement like this:

ampl: include sol.run;
CPLEX 12.10.0.0: optimal integer solution; objective 4
0 MIP simplex iterations
0 branch-and-bound nodes

ampl: display E.slack;
E.slack :=
p03 p04 d05 d07   0
p04 p03 d07 d05   0
;

The output of display shows that there are two E constraints, and that both have a slack of zero in the solution that CPLEX returned -- meaning that the solution does satisfy these constraints.

Use "expand E;" to see exactly the constraints that AMPL generated. You can study these constraints to determine why they are allowing a wrong solution, and then based on what you find, you can fix the model formulation.


--
Robert Fourer
am...@googlegroups.com
{#HS:1207695670-81605#}
On Thu, Jul 2, 2020 at 1:20 PM UTC, AMPL Google Group <am...@googlegroups.com> wrote:
1. If there are too many 1s, pick an element of ARCS at random, and if it is a 1 that is eligible to be a 0, change it to 0. Repeat until the number of 1s is reduced to the desired number. If there are not enough ones, use the same approach but change 0s to 1s. An AMPL "while" loop can be useful for this purpose; the general outline for the "too many 1s" case would be

while (number of 1s in ARCS is > desired number)
choose random i and j
if ARCS[i,j] = 1 and it is allowed to be zero
then let ARCS[i,j] := 0;
}

2. Increasing n will automatically increase the size of your problem. So instead of setting n in a data file, you could consider a loop like this:

for {size in 500..4500 by 500} {
let n := size;
generate a new ARCS matrix
solve;
record results
}


--
Robert Fourer
am...@googlegroups.com
Reply all
Reply to author
Forward
0 new messages