[rge, obviously].
I am trying to learn and use proc optmodel (replacing an older program someone else wrote with proc nlp, which didn't really work as desired anyway.)
Am I doing this at all correctly?
Here is my code:
proc optmodel;
/*decision variables*/
var n {1..&n};
var denom;
var sum_n;
number pop {1..&n};
number p_rate {1..&n};
read data ced_rates into [_n_] pop p_rate;
denom = sum{k in 1..&n} n[k] * p_rate[k];
minimize sum_sqrs = sum {j in 1..&n.} ((n[j] * p_rate[j] / denom) - (pop[j] / &usa_population.))**2;
sum_n = sum{k in 1..&n} n[k];
/*problem constraints*/
con sum_n = &n_usa.,
p_rate[1] * n[1] >= &n_min_urban.,
p_rate[2] * n[2] >= &n_min_urban.,
p_rate[3] * n[3] >= &n_min_urban.,
p_rate[4] * n[4] >= &n_min_urban.,
p_rate[5] * n[5] >= &n_min_urban.,
p_rate[6] * n[6] >= &n_min_urban.,
p_rate[7] * n[7] >= &n_min_urban.,
p_rate[8] * n[8] >= &n_min_urban.,
p_rate[9] * n[9] >= &n_min_urban.,
p_rate[10]*n[10] >= &n_min_rural.,
p_rate[11]*n[11] >= &n_min_rural.,
p_rate[12]*n[12] >= &n_min_rural.,
p_rate[13]*n[13] >= &n_min_rural.,
p_rate[14]*n[14] >= &n_min_rural.,
p_rate[15]*n[15] >= &n_min_rural.,
p_rate[16]*n[16] >= &n_min_rural.,
p_rate[17]*n[17] >= &n_min_rural.,
p_rate[18]*n[18] >= &n_min_rural.,
p_rate[19]*n[19] >= &n_min_urban.,
p_rate[20]*n[20] >= &n_min_urban.,
p_rate[21]*n[21] >= &n_min_urban.,
p_rate[22]*n[22] >= &n_min_urban.,
p_rate[23]*n[23] >= &n_min_urban.,
p_rate[24]*n[24] >= &n_min_urban.,
p_rate[25]*n[25] >= &n_min_urban.,
p_rate[26]*n[26] >= &n_min_urban.,
p_rate[27]*n[27] >= &n_min_urban.,
p_rate[28]*n[28] >= &n_min_urban.,
p_rate[29]*n[29] >= &n_min_urban.,
p_rate[30]*n[30] >= &n_min_urban.,
p_rate[31]*n[31] >= &n_min_urban.,
p_rate[32]*n[32] >= &n_min_urban.,
p_rate[33]*n[33] >= &n_min_urban.,
p_rate[34]*n[34] >= &n_min_urban.,
p_rate[35]*n[35] >= &n_min_urban.,
p_rate[36]*n[36] >= &n_min_urban.,
p_rate[37]*n[37] >= &n_min_urban.,
p_rate[38]*n[38] >= &n_min_urban.,
p_rate[39]*n[39] >= &n_min_urban.,
p_rate[40]*n[40] >= &n_min_urban.,
p_rate[41]*n[41] >= &n_min_urban.;
solve;
print n;
print sum_sqrs;
quit;
(I wanted the number of constraints to be dynamic (&n. = 41), but I couldn't figure out a way to do that.)
I guess the major problem is that the log shows
WARNING: Objective function cannot be evaluated at starting point.
NOTE: Did not converge.
and then a bunch of "NOTE: Division by zero at line"... lines.
Another problem is that, when the program does print out the 41 n values, they do NOT add up to the desired total, &n_usa.
(They are all given values greater than their constraint minimums, so at least that's something.)
The sum_sqrs variable is output as "."
What am I doing wrong? (Am I doing anything right, even?)
Thanks in advance for any help or constructive input.
You are probably not solving the problem that you think you are. You might find the EXPAND statement useful in debugging your model. I suspect that you intended the following two lines to be constraints:
denom = sum{k in 1..&n} n[k] * p_rate[k];
sum_n = sum{k in 1..&n} n[k];
As written, they simply update the values of denom and sum_n, but then the solver is free to change those values. Instead, I think you want:
con denomCon:
denom = sum{k in 1..&n} n[k] * p_rate[k];
con sumCon:
sum_n = sum{k in 1..&n} n[k];
Also, you can write your family of 41 constraints more compactly as follows:
con pCon {j in 1..&n}:
p_rate[j] * n[j] >= if j in 10..18 then &n_min_rural. else &n_min_urban.;
Alternatively, if p_rate[j] is known to be positive, you can omit these explicit constraints and instead modify the lower bound of n[j]:
for {j in 1..&n}
n[j].lb = (if j in 10..18 then &n_min_rural. else &n_min_urban.) / p_rate[j];
OK, thank you for your help so far, but I still get the same [lack of] results:
Solution Summary
Solver NLP
Algorithm Interior Point
Objective Function sum_sqrs
Solution Status Failed
Objective Value .
Optimality Error .
Infeasibility 0
Iterations 0
Presolve Time 0.00
Solution Time 0.00
sum_sqrs
.
(%n n values omitted from results)
Excerpt from Log:
NOTE: Problem generation will use 4 threads.
NOTE: Previous errors might cause the problem to be resolved incorrectly.
NOTE: The problem has 43 variables (43 free, 0 fixed).
NOTE: The problem has 43 linear constraints (0 LE, 2 EQ, 41 GE, 0 range).
NOTE: The problem has 125 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, 41 linear constraints, and 0 nonlinear constraints.
NOTE: The OPTMODEL presolved problem has 43 variables, 2 linear constraints, and 0 nonlinear constraints.
NOTE: The OPTMODEL presolver removed 41 linear constraint coefficients, leaving 84.
NOTE: Using analytic derivatives for objective.
NOTE: Using 2 threads for nonlinear evaluation.
NOTE: The NLP solver is called.
NOTE: The Interior Point algorithm is used.
WARNING: Objective function cannot be evaluated at starting point.
NOTE: Did not converge.
15315 print n;
15316 print sum_sqrs;
NOTE: Division by zero at line 15249 column 63.
Why is it not converging?
I need more information to be able to diagnose this. But you have some error before the log that you posted:
NOTE: Previous errors might cause the problem to be resolved incorrectly.
Can you please share your code and data?
Oops, I didn't realize I had copied from an erroneous run when I tried to change a little bit of the code.
Here is my current PROC OPTMODEL and the results:
proc optmodel;
/*decision variables*/
var n {1..&n};
var denom;
var sum_n;
number pop {1..&n};
number p_rate {1..&n};
read data ced_rates into [_n_] pop p_rate;
/*objective function*/
minimize sum_sqrs = sum {j in 1..&n.} (((n[j] * p_rate[j] / denom) - (pop[j] / &usa_population.))**2);
/*problem constraints*/
con denomCon:
denom = sum{k in 1..&n} (n[k] * p_rate[k]);
con sumCon:
sum_n = sum{k in 1..&n} n[k],
sum_n = &n_usa.;
con pCon {j in 1..&n}:
p_rate[j] * n[j] >=
if j in 10..18 then &n_min_rural.
else &n_min_urban.;
solve;
print n;
print sum_sqrs;
quit;
The OPTMODEL Procedure
Solution Summary
Solver NLP
Algorithm Interior Point
Objective Function sum_sqrs
Solution Status Failed
Objective Value .
Optimality Error .
Infeasibility 12000
Iterations 0
Presolve Time 0.00
Solution Time 0.00
[1] n
1 119.271
2 140.021
(...)
40 135.444
41 140.301
sum_sqrs
.
18421 solve;
NOTE: Problem generation will use 4 threads.
NOTE: The problem has 43 variables (43 free, 0 fixed).
NOTE: The problem has 44 linear constraints (0 LE, 3 EQ, 41 GE, 0 range).
NOTE: The problem has 126 linear constraint coefficients.
NOTE: The problem has 0 nonlinear constraints (0 LE, 0 EQ, 0 GE, 0 range).
NOTE: The OPTMODEL presolver removed 1 variables, 42 linear constraints, and 0 nonlinear constraints.
NOTE: The OPTMODEL presolved problem has 42 variables, 2 linear constraints, and 0 nonlinear constraints.
NOTE: The OPTMODEL presolver removed 43 linear constraint coefficients, leaving 83.
NOTE: Using analytic derivatives for objective.
NOTE: Using 2 threads for nonlinear evaluation.
NOTE: The NLP solver is called.
NOTE: The Interior Point algorithm is used.
WARNING: Objective function cannot be evaluated at starting point.
NOTE: Did not converge.
18422 print n;
18423 print sum_sqrs;
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: Division by zero at line 18356 column 63.
NOTE: The maximum message limit was reached during execution of the statement block. 16 notes and warnings were not
displayed.
18424
18425 quit;
NOTE: PROCEDURE OPTMODEL used (Total process time):
real time 0.12 seconds
cpu time 0.06 seconds
The default initial value for each variable is 0. Because denom appears in a denominator of the objective function, you get a division by zero error. One thing you can try is initializing denom to something nonzero, say 1. You can do this during declaration:
var denom init 1;
Or with an assignment statement before the solve:
denom = 1;
If that doesn't help, please share your data set and macro variables so that I can diagnose further.
The program "works". That is, it compiles and presents a "feasible" solution (one which fits the constraints I gave it), but the solution is not at all realistic. Furthermore, it is highly subjective to initial constraints (such as changing the var sum_n line from the default of 0 to 1 or 10 or 308745538). It also provides a quite different solution if I specify PRESOLVE=BASIC rather than the default, but not one I like much better.
Shouldn't the sum of squares be a convex function, with only one local minimum?
Now I'm just not sure what to do, because my program seems to work syntactically, but I just don't like its output.
Unfortunately, I don't think I (should / am authorized to / legally can ?) publicly post my data set and macro variables.
Yes, the sum of squares of convex functions is convex, but the expression you are squaring is not convex. As a simpler example, consider (sqrt(f(x)))^2 + (sqrt(f(x)))^2, where f(x) is any nonconvex function.
To increase the likelihood of finding a globally optimal solution for a nonconvex minimization (or nonconcave maximization) problem, you can use the MULTISTART option:
solve with nlp / multistart;
Also, you can simplify the code a bit by omitting the sum_n variable and merging two of the constraints into one:
con sumCon:
sum{k in 1..&n} n[k] = &n_usa.;
From your log, it looks like the presolver is doing that anyway.
Registration is now open for SAS Innovate 2025 , our biggest and most exciting global event of the year! Join us in Orlando, FL, May 6-9.
Sign up by Dec. 31 to get the 2024 rate of just $495.
Register now!
Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.
Find more tutorials on the SAS Users YouTube channel.