Minimize a function f(x) over a given domain - CLP - OPTMODEL

[harmonic]

Hello, 

I would like to understand if there is a way to find the best configuration to minimize a function over a given domain of variables x, y, z:

For example using default values:
xij like and list of arrays {3 , 4 ,5}, y {6, 7 , 8}, z{1, 2 , 3}

Vi {100, 200, 300}

min f(x) = sum (xij*Vi)

The best configuration will be 3*100 + 6*200 + 1*300 

I am staring to check the CLP and the OPTMODEL, do you have real example for this kind of problem?

Thank you

@RobPratt 

[harmonic]

Hello, 

I would like to understand if there is a way to find the best configuration to minimize a function over a given domain of variables x, y, z:

For example using default values:
xij like and list of arrays {3 , 4 ,5}, y {6, 7 , 8}, z{1, 2 , 3}

Vi {100, 200, 300}

min f(x) = sum (xij*Vi)

The best configuration will be 3*100 + 6*200 + 1*300 

I am staring to check the CLP and the OPTMODEL, do you have real example for this kind of problem?

Thank you

[Ksharp]

You'd better post it at OR forum.

https://communities.sas.com/t5/Mathematical-Optimization/bd-p/operations_research

and calling @RobPratt 

Here is an example:

proc optmodel;
set idx={1..9};
num x{idx}=[3    4   5  6   7   8   1   2    3];
num v{idx}=[100 100 100 200 200 200 300 300 300];
var a{idx} binary;

min obj=sum{i in idx} a[i]*x[i]*v[i];

con con1:sum{i in 1..3} a[i]=1;
con con2:sum{j in 4..6} a[j]=1;
con con3:sum{k in 7..9} a[k]=1;

solve with clp / findallsolns;

print obj a x v;
quit;

[harmonic]

Thank you,

I will focus on my real problem.

1. I have two different types of data inputs and I want to determine which one is better for my objective function.

2. Next, I would like to create different arrays from the input in the procedure as follows (maybe using read data into):

- vol_M1: {100, 100}
- pcs_M1: {10, 10.5}
- vol_S1: {300, 300}
- pcs_S1: {11, 11.5}
- vol_S2: {400, 400}
- pcs_S2: {12}

3. After that, I need to determine the best configuration using all combinations to minimize the absolute value of 'f'.

The minimum absolute value of 'f' is calculated as follows:

Min |f| = vol_M1_i * pcs_M1_i – [(vol_s1_j * pcs_S1_j) + (vol_S2_k * pcs_S2_k)]

Where:
- i ranges from 1 to the dimension of pcs_M1
- j ranges from 1 to the dimension of pcs_S1
- k ranges from 1 to the dimension of pcs_S2

4. How can I implement this in Proc OPTMODEL?

This is a sample, but I have 7000 different "cod" items.

@RobPratt @Ksharp 

[RobPratt]

I don't quite understand.  Can you please do the following to help clarify?

First describe the business problem, without algebra or explicit data.  Maybe trying to match supplies and demands as closely as possible?  Then express the objective function explicitly using your small sample data.  Finally, provide the expected optimal solution for your sample data.

Thanks,

Rob

[harmonic]

I will repost my answer, I don't know why I can't find it.

My business solution is to find the best configuration to minimize the obj function with a clp or a optmodel because with a for cycle there will be a lot of iterations

For example using these arrays in a for cycle:

vol_M1: {100, 100}

pcs_M1: {10, 10.5}

vol_S1: {300, 300}

pcs_S1: {11, 11.5}

vol_S2: {400, 400}

pcs_S2: {12}

Min |f| = vol_M1_i * pcs_M1_i – [(vol_s1_j * pcs_S1_j) + (vol_S2_k * pcs_S2_k)]

1. 100*10 - [(300*11) + (400*12)] = -7100

2. 100*10.5 - [(300*11) + (400*12)] = -7050

3. 100*10 - [(300*11.5) + (400*12)] =  - 7250

4. 100*10.5 - [(300*11.5) + (400*12)]  = -7200

The best configuration will be the 2. because is a minimum local (next to zero)

[RobPratt]

The following does what I think you want:

data VData;
   input i v;
   datalines;
1 100
2 300
3 400
;

data PData;
   input i j p;
   datalines;
1 1 10
1 2 10.5
2 1 11
2 2 11.5
3 1 12
;

proc optmodel;
   set ISET;
   num v {ISET};
   read data VData into ISET=[i] v;
   set <num,num> IJ;
   num p {IJ};
   read data PData into IJ=[i j] p;

   var X {IJ} binary;

   min f = abs(sum {<i,j> in IJ} (if i = 1 then 1 else -1)*v[i]*p[i,j]*X[i,j]);

   con OneValue {i in ISET}:
      sum {<(i),j> in IJ} X[i,j] = 1;

   solve linearize;

   print X;
   print {i in ISET} (sum {<(i),j> in IJ} p[i,j]*X[i,j]);
quit;

The SAS System

The OPTMODEL Procedure

Problem Summary
Objective Sense	Minimization
Objective Function	f
Objective Type	Nonlinear
Number of Variables	5
Bounded Above	0
Bounded Below	0
Bounded Below and Above	5
Free	0
Fixed	0
Binary	5
Integer	0
Number of Constraints	3
Linear LE (<=)	0
Linear EQ (=)	3
Linear GE (>=)	0
Linear Range	0
Constraint Coefficients	5

---

The SAS System

The OPTMODEL Procedure

Solution Summary
Solver	MILP
Algorithm	Branch and Cut
Objective Function	f
Solution Status	Optimal
Objective Value	7050
Relative Gap	0
Absolute Gap	0
Primal Infeasibility	0
Bound Infeasibility	0
Integer Infeasibility	0
Best Bound	7050
Nodes	0
Solutions Found	2
Iterations	0
Presolve Time	0.02
Solution Time	0.02

X
	1	2
1	0	1
2	1	0
3	1

[1]
1	10.5
2	11.0
3	12.0

[harmonic]

Thank you,

Really nice solution, and if I have n variables do you think this solution will be acceptable?  

data have;
input cod $ vol pcs;
cards;
M1 100 10
M2 100 10.5

.......

Mn 150 11
S1 300 11
S1 300 11.5
S2 400 12
S2 400 .

......

Sn 200 10.5
;
run;

Where M got just one p for each observation and S different p. M and S are the variable vol.

[RobPratt]

The following will work even if you have more than one p for each M observation.  It is more efficient and less error-prone to avoid duplicating data, so I have retained two separate data sets:

data VData;
   input cod $ vol;
   datalines;
M1 100
M2 100
Mn 150
S1 300
S2 400
Sn 200
;

data PData;
   input cod $ j pcs;
   datalines;
M1 1 10
M2 1 10.5
Mn 1 11
S1 1 11
S1 2 11.5
S2 1 12
Sn 1 10.5
;

proc optmodel;
   set <str> ISET;
   num v {ISET};
   read data VData into ISET=[cod] v=vol;
   set <str,num> IJ;
   num p {IJ};
   read data PData into IJ=[cod j] p=pcs;

   var X {IJ} binary;

   min f = abs(sum {<i,j> in IJ} (if char(i,1) = 'M' then 1 else -1)*v[i]*p[i,j]*X[i,j]);

   con OneValue {i in ISET}:
      sum {<(i),j> in IJ} X[i,j] = 1;

   solve linearize;

   print X;
   print {i in ISET} (sum {<(i),j> in IJ} p[i,j]*X[i,j]);
quit;

[harmonic]

Thank you for you quick reply, I don't understand how to keep the value char, the final problem is more complex then this one, I created sample tables like that:

data VData;
input i v remi $ class $;
datalines;
1 100 205 MA
2 300 802 MB
3 600 605 PB
4 200 604 PB
5 300 725 SI
6 500 215 SI
;

data PData;
input i j p remi $ aop $ class $;
datalines;
1 1 10 205 15 MA
2 1 11.25 802 18 MB
3 1 10 605 26 PB
4 1 11.5 604 35 PB
5 1 10.25 725 44 SI
5 2 12 725 42 SI

5 3 110.75 725 46 SI
6 1 10.25 215 87 SI
6 2 11.25 215 89 SI

;

The objective is to output the configuration of "aop" and "pcs" codes to minimize abs(f)

tot_in=sum(sum(of M:), -1 * sum(of PB:));
tot_out=sum(sum(of P:), sum(of SI:), -1 * sum(of PB:));
f=sum(tot_in, -tot_out);

NB: just the points with class "SI" got different "pcs" for the same "remi" code.

[RobPratt · Accepted solution]

With your new data, you can omit i and j, and I think the following does what you want:

data VData;
   input v remi $ class $;
   datalines;
100 205 MA
300 802 MB
600 605 PB
200 604 PB
300 725 SI
500 215 SI
;

data PData;
   input p remi $ aop $ class $;
   datalines;
10 205 15 MA
11.25 802 18 MB
10 605 26 PB
11.5 604 35 PB
10.25 725 44 SI
12 725 42 SI
110.75 725 46 SI
10.25 215 87 SI
11.25 215 89 SI
;

proc optmodel;
   set <str,str> REMI_CLASS;
   num v {REMI_CLASS};
   read data VData into REMI_CLASS=[remi class] v;
   set <str,str,str> REMI_CLASS_AOP;
   num p {REMI_CLASS_AOP};
   read data PData into REMI_CLASS_AOP=[remi class aop] p;

   var X {REMI_CLASS_AOP} binary;

   min f = abs(sum {<r,c,a> in REMI_CLASS_AOP} (if char(c,1) = 'M' then 1 else -1)*v[r,c]*p[r,c,a]*X[r,c,a]);

   con OneValue {<r,c> in REMI_CLASS}:
      sum {<(r),(c),a> in REMI_CLASS_AOP} X[r,c,a] = 1;

   solve linearize;

   print X;
   print {<r,c> in REMI_CLASS} (sum {<(r),(c),a> in REMI_CLASS_AOP} p[r,c,a]*X[r,c,a]);
quit;

[harmonic]

Thank you,

and if I would like to add constraints for the p value like 9<p<12?

And try to change the solver with CLP or MILP?

[RobPratt]

To restrict the p values, you can either modify the input data or use a WHERE clause in the READ DATA statement:

   read data PData(where=(9 < p < 12)) into REMI_CLASS_AOP=[remi class aop] p;

You cannot use the CLP solver here because the variables needed to linearize the absolute value are continuous.

[harmonic]

1. My goal is to use the p values just with his correspond remi and aop number , how is possible that it will find everytime function next to zero with just a couple of p available for each couple of r and a.

2. And can I create a table output with r c a, for values of v p and x?

3. For some reason I get this error, I will upload the dataset input.

NOTE: Problem generation will use 4 threads.
NOTE: The problem has 22246 variables (0 free, 0 fixed).
NOTE: The problem has 22246 binary and 0 integer variables.
NOTE: The problem has 7215 linear constraints (0 LE, 7215 EQ, 0 GE, 0 range).
NOTE: The problem has 22246 linear constraint coefficients.
NOTE: The problem has 0 nonlinear constraints (0 LE, 0 EQ, 0 GE, 0 range).
NOTE: The OPTMODEL presolver removed 1567 variables, 1567 linear constraints, and 0 nonlinear constraints.
NOTE: The OPTMODEL presolver replaced 0 nonlinear constraints, 0 objectives, and 0 implicit variables.
NOTE: The OPTMODEL presolver added 0 variables and 0 linear constraints.
NOTE: The OPTMODEL presolved problem has 20679 variables, 5648 linear constraints, and 0 nonlinear constraints.
NOTE: The OPTMODEL presolver removed 1567 linear constraint coefficients, leaving 20679.
ERROR: The NLP solver does not allow integer variables. The black-box solver might be suitable for this problem.

4. Anyway I found this kind of bug for the X values, it should be 0 or 1.

5.Can I use bounds in constraint statement?

@RobPratt 

[RobPratt]

1. Yes, you can use an explicit constraint:

   con Pcon {<r,c,a> in REMI_CLASS_AOP: not(9 < p[r,c,a] < 12)}:
      X[r,c,a] = 0;

But it is more efficient to omit those <r,c,a> tuples.  See this doc example.

This is an answer to original question #1.  Sorry, I don't understand your edited question #1.

2. You can use the CREATE DATA statement.  Here are two alternatives:

   create data SolData1 from [r c a] v[r,c] p X;
   create data SolData2 from [r c a]={<r,c,a> in REMI_CLASS_AOP: X[r,c,a].sol > 0.5} v[r,c] p X;

3. What SOLVE statement did you use?

4. That is not a bug.  Those values are within the default integer tolerance of 1e-5.

5. Not sure what you mean.  Can you please clarify?