topic Re: Pls help me check this code in Mathematical Optimization, Discrete-Event Simulation, and OR

Pls help me check this code

violet2386 — Thu, 20 Jun 2013 03:49:10 GMT

I am using SAS to solve a nonlinear problem.

I once tried "solve with QP", but the calculated result is wrong. (I use matlab fmincon to solve the same problem, and the two results are different. The matlab result is much more convincing, and I believe the result given by matlab is correct). However, when I tried "solve;" in my code, the SAS just can not solve it. In the note, it indicated that The SAS System stopped processing this step because of errors. Can someone help me check my code? Thanks so much.

data data_read;

input x1 x2 y1 dmu_id;

datalines;

20.22 1 85 1

95.51 1 25 2

1662.98 5 22 3

;

/*---------solve QP--------*/

proc optmodel printlevel=0;

set x_num=1..2; * set dimension of X ;

set y_num=1..1;

set <num>DMU_ID;

num X{DMU_ID,x_num};

num Y{DMU_ID,y_num};

read data data_read into DMU_ID=[dmu_id]

{r in x_num}<X[dmu_id,r]=col("x"||r)>

{s in y_num}<Y[dmu_id,s]=col("y"||s)>;

var epi{DMU_ID}>=0;

var alpha{DMU_ID};

var beta{DMU_ID,x_num}>=0;

num epi_sol{DMU_ID};

num alpha_sol{DMU_ID};

num beta_sol{DMU_ID,x_num};

min obj=sum{q in DMU_ID}epi**2;

con regression{q in DMU_ID,s in y_num}:Y[q,s]=alpha-epi+sum{r in x_num}X[q,r]*beta[q,r];

con concave{q in DMU_ID,h in DMU_ID}:alpha+sum{r in x_num}X[q,r]*beta[q,r]<=alpha+sum{r in x_num}X[q,r]*beta[h,r];

solve;

for{q in DMU_ID}epi_sol=epi.sol;

for{q in DMU_ID}alpha_sol=alpha.sol;

for{q in DMU_ID,r in x_num}beta_sol[q,r]=beta[q,r].sol;

create data epi_sol from [dmu_id]epi_sol;

create data alpha_sol from [dmu_id]alpha_sol;

create data beta_sol from [dmu_id x_num]beta_sol;

quit;

run;

Re: Pls help me check this code

LeoLopes — Thu, 20 Jun 2013 13:28:02 GMT

I can't reproduce the error on SAS 9.3. There are a few minor messages that can be addressed, but no errors. Could you please attach your log file?

Here is what I get:

------

172 data data_read;

173 input x1 x2 y1 dmu_id;

174 datalines;

NOTE: The data set WORK.DATA_READ has 3 observations and 4 variables.

NOTE: DATA statement used (Total process time):

real time 0.00 seconds

cpu time 0.01 seconds

178 ;

179 /*---------solve QP--------*/

180 proc optmodel printlevel=0;

181 set x_num=1..2;

181! * set dimension of X ;

182 set y_num=1..1;

183 set <num>DMU_ID;

184 num X{DMU_ID,x_num};

185 num Y{DMU_ID,y_num};

186 read data data_read into DMU_ID=[dmu_id]

777

NOTE 777-782: The name 'dmu_id' hides an outer declaration.

187 {r in x_num} <X[dmu_id,r]=col("x"||r)>

188 {s in y_num} <Y[dmu_id,s]=col("y"||s)>;

NOTE: There were 3 observations read from the data set WORK.DATA_READ.

189

190 var epi{DMU_ID}>=0;

191 var alpha{DMU_ID};

192 var beta{DMU_ID,x_num}>=0;

193

194 min obj= sum{q in DMU_ID} epi**2;

195 con regression{q in DMU_ID,s in y_num}:

196 Y[q,s] = alpha

197 - epi

198 + sum{r in x_num} X[q,r]*beta[q,r];

199 con concave{q in DMU_ID,h in DMU_ID}:

200 alpha + sum{r in x_num} X[q,r]*beta[q,r]

201 <= alpha + sum{r in x_num} X[q,r]*beta[h,r];

202 solve;

NOTE: Problem generation will use 4 threads.

NOTE: The constraint 'concave[1,1]' is empty.

NOTE: The constraint 'concave[2,2]' is empty.

NOTE: The constraint 'concave[3,3]' is empty.

NOTE: The problem has 12 variables (3 free, 0 fixed).

NOTE: The problem has 12 linear constraints (9 LE, 3 EQ, 0 GE, 0 range).

NOTE: The problem has 48 linear constraint coefficients.

NOTE: The problem has 0 nonlinear constraints (0 LE, 0 EQ, 0 GE, 0 range).

NOTE: The OPTMODEL presolver removed 0 variables, 3 linear constraints, and 0 nonlinear

constraints.

NOTE: The OPTMODEL presolved problem has 12 variables, 9 linear constraints, and 0 nonlinear

constraints.

NOTE: The QP presolver value AUTOMATIC is applied.

NOTE: The QP presolver removed 0 variables and 0 constraints.

NOTE: The QP presolver removed 0 constraint coefficients.

NOTE: The presolved problem has 12 variables, 9 constraints, and 48 constraint coefficients.

NOTE: The QP solver is called.

NOTE: The Interior Point algorithm is used.

NOTE: The deterministic parallel mode is enabled.

NOTE: The Interior Point algorithm is using up to 4 threads.

Primal Bound Dual

Iter Complement Duality Gap Infeas Infeas Infeas

0 2.91929e+03 5.84738e-01 5.84520e-01 6.60197e+01 3.94884e+04

1 2.16097e+03 1.42333e+00 3.60634e-01 4.07325e+01 2.43634e+04

2 1.37176e+03 5.09543e+00 7.69786e-02 8.69449e+00 5.20044e+03

3 3.89368e+02 2.27760e+00 2.22052e-02 2.50801e+00 1.50012e+03

4 3.92930e+02 2.51196e+00 2.18544e-02 2.46839e+00 1.47642e+03

5 4.47907e+02 1.21485e+01 2.08696e-02 2.35716e+00 1.40989e+03

6 2.97149e+03 1.80112e-01 1.05151e-02 1.18765e+00 7.10371e+02

7 2.27681e+03 1.25190e-01 7.63588e-03 8.62449e-01 5.15858e+02

8 1.16576e+03 4.90437e-02 3.61949e-03 4.08810e-01 2.44522e+02

9 4.99641e+02 2.16038e-02 1.44582e-03 1.63301e-01 9.76752e+01

10 4.49887e+02 3.53944e-02 1.20773e-03 1.36409e-01 8.15905e+01

11 4.09736e+02 1.85557e-01 1.00743e-03 1.13786e-01 6.80590e+01

12 6.36467e+03 1.65790e-01 5.39117e-04 6.08916e-02 3.64211e+01

13 2.43454e+02 2.21331e-02 8.60910e-06 9.72371e-04 5.81605e-01

14 2.58773e+00 2.46941e-04 8.67933e-08 9.80304e-06 5.86350e-03

15 1.77261e-02 2.46950e-06 7.33576e-07 8.64549e-08 2.56299e-05

16 1.77261e-04 2.46950e-08 7.33576e-09 8.64559e-10 2.56299e-07

17 1.77269e-06 2.46962e-10 7.33592e-11 8.64698e-12 2.56296e-09

18 1.77269e-08 2.47041e-12 7.37477e-13 1.37854e-13 2.59279e-11

NOTE: Optimal.

NOTE: Objective = 7569.

NOTE: The Interior Point solve time is 0.02 seconds.

203

204 *num epi_sol{q in DMU_ID} = epi.sol;

205 num alpha_sol{DMU_ID};

206 num beta_sol{DMU_ID,x_num};

207

208 *for{q in DMU_ID}epi_sol= epi.sol;

209 for{q in DMU_ID}alpha_sol=alpha.sol;

210 for{q in DMU_ID,r in x_num}beta_sol[q,r]=beta[q,r].sol;

211 create data epi_sol from [DMU_ID] epi;

NOTE: The data set WORK.EPI_SOL has 3 observations and 2 variables.

212 create data alpha_sol from [DMU_ID]alpha_sol;

NOTE: The data set WORK.ALPHA_SOL has 3 observations and 2 variables.

213 create data beta_sol from [DMU_ID x_num]beta_sol;

NOTE: The data set WORK.BETA_SOL has 6 observations and 3 variables.

214 quit;

NOTE: PROCEDURE OPTMODEL used (Total process time):

real time 0.06 seconds

cpu time 0.06 seconds

------

Thanks,

Leo.

Re: Pls help me check this code

violet2386 — Thu, 20 Jun 2013 14:15:28 GMT

Thanks so much. I use SAS 9.2 and below is the log file. It seems that my SAS version does have some problem. However, the optimal obj calculated from matlab is zero. I hope to find what is something wrong in my code. Thanks again.

---------------------------------------------------------------------------------------------

385 data data_read;

386 input x1 x2 y1 dmu_id;

387 datalines;

NOTE: The data set WORK.DATA_READ has 3 observations and 4 variables.

NOTE: DATA statement used (Total process time):

real time 0.03 seconds

cpu time 0.01 seconds

391 ;

392 proc optmodel printlevel=0;

393 set x_num=1..2;

393! * set dimension of X ;

394 set y_num=1..1;

395 set <num>DMU_ID;

396 num X{DMU_ID,x_num};

397 num Y{DMU_ID,y_num};

398 read data data_read into DMU_ID=[dmu_id]

399 {r in x_num}<X[dmu_id,r]=col("x"||r)>

400 {s in y_num}<Y[dmu_id,s]=col("y"||s)>;

NOTE: There were 3 observations read from the data set WORK.DATA_READ.

401 var epi{DMU_ID}>=0;

402 var alpha{DMU_ID};

403 var beta{DMU_ID,x_num}>=0;

404 num epi_sol{DMU_ID};

405 num alpha_sol{DMU_ID};

406 num beta_sol{DMU_ID,x_num};

407 min obj=sum{q in DMU_ID}epi**2;

408 con regression{q in DMU_ID,s in y_num}:Y[q,s]=alpha-epi+sum{r in x_num}X[q,r]*beta[q,r];

409 con concave{q in DMU_ID,h in DMU_ID}:alpha+sum{r in x_num}X[q,r]*beta[q,r]<=alpha+sum{r in x_num}X[q,r]*beta[h,r];

410 solve;

NOTE: The constraint 'concave[1,1]' is empty.

NOTE: The constraint 'concave[2,2]' is empty.

NOTE: The constraint 'concave[3,3]' is empty.

NOTE: The problem has 12 variables (3 free, 0 fixed).

NOTE: The problem has 12 linear constraints (9 LE, 3 EQ, 0 GE, 0 range).

NOTE: The problem has 48 linear constraint coefficients.

NOTE: The problem has 0 nonlinear constraints (0 LE, 0 EQ, 0 GE, 0 range).

NOTE: The OPTMODEL presolver removed 0 variables, 3 linear constraints, and 0 nonlinear constraints.

NOTE: The OPTMODEL presolved problem has 12 variables, 9 linear constraints, and 0 nonlinear constraints.

NOTE: Using analytic derivatives for objective.

NOTE: Initial point was changed to be feasible to bound and linear constraints.

NOTE: The trust region method is used.

NOTE: The SAS System stopped processing this step because of errors.

NOTE: PROCEDURE OPTMODEL used (Total process time):

real time 0.03 seconds

cpu time 0.03 seconds

411 for{q in DMU_ID}epi_sol=epi.sol;

412 for{q in DMU_ID}alpha_sol=alpha.sol;

413 for{q in DMU_ID,r in x_num}beta_sol[q,r]=beta[q,r].sol;

414 create data epi_sol from [dmu_id]epi_sol;

415 create data alpha_sol from [dmu_id]alpha_sol;

416 create data beta_sol from [dmu_id x_num]beta_sol;

417 quit;

418 run;

Re: Pls help me check this code

LeoLopes — Thu, 20 Jun 2013 18:09:48 GMT

I am glad to help.

I have forwarded your model to someone who might be able to reproduce your problem.

There are a couple of things that may help further, if you are able to provide them:

1. Your Matlab code

2. The very top of the log file, where the detailed version information is.

Cheers,

Leo.

Re: Pls help me check this code

LeoLopes — Thu, 20 Jun 2013 21:07:27 GMT

Hi Violet,

Our folks suggested that if you can, please upgrade to at least 9.3. This problem does not occur there.

Technical Support may also be able to help: http://support.sas.com/techsup/contact/

Good luck and thanks for sharing!

Leo.

Re: Pls help me check this code

violet2386 — Fri, 21 Jun 2013 03:55:09 GMT

Thanks. Yes, I should upgrade my SAS version.

Re: Pls help me check this code

LeoLopes — Sat, 22 Jun 2013 02:12:56 GMT

You're welcome!