\$title Test of correctness for levels & marginals of QCP (QCP01,SEQ=74) \$ontext Test of correctness of the levels and marginals returned. We have QP terms in the objective only - all QCP solvers accept this. All cases are considered, e.g. 1) =L=, =G=, =E= constraints (should we add =N=?) 2) variables at lower and upper bound, and not at bound 3) min [convex obj] or max [concave obj] 4) special attention paid to the form of the obj. constraint, i.e. cz * z = xQx + cx + b where cz and b take different values \$offtext \$if not set MTYPE \$set MTYPE qcp \$if not set TESTTOL \$set TESTTOL 1e-6 scalar mchecks / 0 /; \$if not %QPMCHECKS% == 0 mchecks = 1; scalar failed / 0 /; \$escape = \$echo if{%=1, display '%=1 failed', '%=2'; failed=1}; > gtest.gms \$escape % set J / j1 * j7 /; set GT / gt1 * gt2 /; set LT / lt1 * lt2 /; set EQ / eq1 /; scalars cz, cb; parameters c(J) / j1 0 j2 0 j3 1 j4 2 j5 -1 j6 -1 j7 3 / bgt(GT) / gt1 2 gt2 4 /, blt(LT) / lt1 -17 lt2 10 /, beq(EQ) / eq1 9 /; table Agt(GT,J) j1 j2 gt1 1 1 gt2 4 1 ; table Alt(LT,J) j3 j4 lt1 -6 1 lt2 1 2 ; parameter Aeq(EQ,J); Aeq('eq1',j)\$(ord(j) le 7) = 1; variable x(J) z; * these bounds should never be active - but they help the global solvers x.lo(J) = -2; x.up(J) = 5; * these are active and tested to be so x.lo('j5') = 1; x.up('j6') = 3; equations gte(GT), lte(LT), eqe(EQ), obj; obj.. cz*z =E= cb + sum{J,c(J)*x(J)} + sqr(x('j1')) + sqr(x('j2')) - x('j1')*x('j2') + sqr(x('j3')-x('j4')) + sqr(x('j5')+1) + sqr(x('j6')-4) + sqr(x('j7')-1); gte(GT).. sum{J, Agt(GT,J)*x(J)} =G= bgt(GT); lte(LT).. sum{J, Alt(LT,J)*x(J)} =L= blt(LT); eqe(EQ).. sum{J, Aeq(EQ,J)*x(J)} =E= beq(EQ); model m / all /; set czvals 'obj multipliers' / 'cz=1', 'cz=0.5' 'cz=3' 'cz=-1' 'cz=-0.5' 'cz=-3' /; parameter czv(czvals) / 'cz=1' 1 'cz=0.5' 0.5 'cz=3' 3 'cz=-1' -1 'cz=-0.5' -0.5 'cz=-3' -3 /; set cbvals 'obj constants' / 'cb=0' 'cb=2' 'cb=-2' /; parameter cbv(cbvals) / 'cb=0' 0 'cb=2' 2 'cb=-2' -2 /; scalars isMin / 0 /, tol / %TESTTOL% /, obj_l / 0.0 / obj_m / 1.0 / objval / 12.0 /; parameters x_l(J) / j1 1.0 j2 1.0 j3 3.0 j4 1.0 j5 1.0 j6 3.0 j7 -1.0 / x_m(J) / j1 0.0 j2 0.0 j3 0.0 j4 0.0 j5 4.0 j6 -2.0 j7 0.0 / gte_l(GT) / gt1 2.0 gt2 5.0 / gte_m(GT) / gt1 2.0 gt2 0.0 / lte_l(LT) / lt1 -17.0 lt2 5.0 / lte_m(LT) / lt1 -1.0 lt2 0.0 / eqe_l(EQ) / eq1 9.0 / eqe_m(EQ) / eq1 -1.0 /; loop {czvals\$[ord(czvals) <= INF], cz = czv(czvals); loop {cbvals\$[ord(cbvals) <= INF], cb = cbv(cbvals); if {(cz > 0), isMin = 1; solve m using %MTYPE% min z; else isMin = -1; solve m using %MTYPE% max z; } \$ batinclude gtest "( m.solvestat <> %solvestat.NormalCompletion% or (m.modelstat > %modelstat.LocallyOptimal% and m.modelstat <> %modelstat.FeasibleSolution%))" "wrong status codes" \$ batinclude gtest "(abs(cz*z.l-cb-objval) > tol)" "bad z.l" \$ batinclude gtest "(abs(obj.l-cb) > tol)" "bad obj.l" \$ batinclude gtest "(smax(J, abs(x.l(j)-x_l(j))) > tol)" "bad x.l" \$ batinclude gtest "(smax(GT,abs(gte.l(GT)-gte_l(GT))) > tol)" "bad gte.l" \$ batinclude gtest "(smax(LT,abs(lte.l(LT)-lte_l(LT))) > tol)" "bad lte.l" \$ batinclude gtest "(smax(EQ,abs(eqe.l(EQ)-eqe_l(EQ))) > tol)" "bad eqe.l" if {mchecks, \$ batinclude gtest "(abs(z.m) > tol)" "bad z.m" \$ batinclude gtest "(abs(cz*obj.m-obj_m) > tol)" "bad obj.m" \$ batinclude gtest "(smax(J, abs(cz*x.m(j)-x_m(j))) > tol)" "bad x.m" \$ batinclude gtest "(smax(GT,abs(cz*gte.m(GT)-gte_m(GT))) > tol)" "bad gte.m" \$ batinclude gtest "(smax(LT,abs(cz*lte.m(LT)-lte_m(LT))) > tol)" "bad lte.m" \$ batinclude gtest "(smax(EQ,abs(cz*eqe.m(EQ)-eqe_m(EQ))) > tol)" "bad eqe.m" }; \$if set doscalar if (failed, option %MTYPE%=convert; if {(cz > 0), solve m using %MTYPE% min z; else solve m using %MTYPE% max z;}) abort\$failed 'test failed', cz, cb; }; };