Bi Level Optimization (solve unable to obtain explicit solution)

I am having sove command issues with a bi level optimization problem. Vpasolve did not address the issue. I can solve the problem of Maple so I know that solution EXISTS.

Bi Level Optimization with High Degree Polynomials

Code:

syms x y z t s; %x y z upper level optmization variables and t s are lower levels
% Upper Level Optimization
f=-0.008*x^2*y - 0.0315*x^2*z - 0.002*x*y + 0.0665*x*z + s*z - t*y - y^2 - 1.4*y*z - 0.7*z^2 + 0.91*y + 0.6650*z;
f_x=diff(f,x);
f_y=diff(f,y);
f_z=diff(f,z);
eqns = [f_x == 0, f_y == 0,f_z==0];
S = solve(eqns,[x y z]);
%Lower Level Optimization
f2= -0.7*y*z - 0.5000000000*y^2 - 0.3500000000*z^2 + 0.4776*y + 0.2176*z + 0.1822*x*z - 0.008*x^2*y - 0.0315*x^2*z - 0.002*x*y;
f3=subs(f2,S);
f_t=diff(f3,t);
f_s=diff(f3,s);
eqnsl = [f_t ==0, f_s == 0];
SL = solve([f_t==0,f_s==0],[t s]);

Warning: Unable to find explicit solution. For options, see help.

[Shumail Mazahir]

Thanks Walter. You are a life saver.

Answer 1

[Walter Roberson]

syms x y z t s; %x y z upper level optmization variables and t s are lower levels
% Upper Level Optimization
f=-0.008*x^2*y - 0.0315*x^2*z - 0.002*x*y + 0.0665*x*z + s*z - t*y - y^2 - 1.4*y*z - 0.7*z^2 + 0.91*y + 0.6650*z;
f_x=diff(f,x);
f_y=diff(f,y);
f_z=diff(f,z);
eqns = [f_x == 0, f_y == 0,f_z==0];

%it is a cubic in eqns, let us find explicit solutions
S = solve(eqns, [x y z], 'MaxDegree', 3)

S = struct with fields:
    x: [3×1 sym]
    y: [3×1 sym]
    z: [3×1 sym]

%Lower Level Optimization
f2= -0.7*y*z - 0.5000000000*y^2 - 0.3500000000*z^2 + 0.4776*y + 0.2176*z + 0.1822*x*z - 0.008*x^2*y - 0.0315*x^2*z - 0.002*x*y;
f3=subs(f2,S);
f_t=diff(f3,t);
f_s=diff(f3,s);
eqnsl = [f_t ==0; f_s == 0];

size(eqnsl)

ans = 1×2

6 1

symvar(eqnsl)

ans = 

%look at that: you have 6 equations in two variables. You will never be
%able to solve that



%but perhaps you can solve them in pairs

%but that takes too long for the Answers system to respond too. For your
%real work, remove the conditional test. And be prepared to wait a fair
%while
if false
    alternative_SL = arrayfun(@(E1,E2) vpasolve([E1,E2], [t s]), f_t == 0, f_s == 0, 'uniform', 0);
    celldisp(alternative_SL)
end

[Walter Roberson]

The results come out as

t is

[0.18750535865719207556855836628588, 0.18750535865719207556855836628588, -0.0050814631775379818251855413854812]

s is

[1.4822538940960940679193398562996, 1.4822538940960940679193398562996, 0.022317295109322082845757335468151]

The first two results are identical to within 10^-37 .

Answer 2

[Torsten]

I couldn't test it:

syms x y z t s %x y z upper level optmization variables and t s are lower levels
% Upper Level Optimization
f=-0.008*x^2*y - 0.0315*x^2*z - 0.002*x*y + 0.0665*x*z + s*z - t*y - y^2 - 1.4*y*z - 0.7*z^2 + 0.91*y + 0.6650*z;
f_x=diff(f,x);
f_y=diff(f,y);
f_z=diff(f,z);
eqns = [f_x == 0, f_y == 0,f_z==0];
S = solve(eqns,[x y z],'MaxDegree',3)
S.x
S.y
S.z
%Lower Level Optimization
f2 = -0.7*y*z - 0.5000000000*y^2 - 0.3500000000*z^2 + 0.4776*y + 0.2176*z + 0.1822*x*z - 0.008*x^2*y - 0.0315*x^2*z - 0.002*x*y;
f31 = subs(f2,[x y z],[S.x(1) S.y(1) S.z(1)])
f32 = subs(f2,[x y z],[S.x(2) S.y(2) S.z(2)])
f33 = subs(f2,[x y z],[S.x(3) S.y(3) S.z(3)])
f1_t=diff(f31,t)
f1_s=diff(f31,s)
f2_t=diff(f32,t)
f2_s=diff(f32,s)
f3_t=diff(f33,t)
f3_s=diff(f33,s)
%eqnsl = [f_t ==0, f_s == 0]
SL1 = vpasolve([f1_t ==0, f1_s == 0],[t s])
SL2 = vpasolve([f2_t ==0, f2_s == 0],[t s])
SL3 = vpasolve([f3_t ==0, f3_s == 0],[t s])