%nested1
%outer done
%inner done

clear;clc
sdpvar x k m y n
sdpvar pc0 pc1 pc2 pc3 pc4 pc5 pc11 pc12 pc13 pc14 pc15 pc22 pc23 pc24 pc25 pc33 pc34 pc35 pc44 pc45 pc55
sdpvar pd0 pd1 pd2 pd3 pd4 pd5 pd11 pd12 pd13 pd14 pd15 pd22 pd23 pd24 pd25 pd33 pd34 pd35 pd44 pd45 pd55
sdpvar pe0 pe1 pe2 pe3 pe4 pe5 pe11 pe12 pe13 pe14 pe15 pe22 pe23 pe24 pe25 pe33 pe34 pe35 pe44 pe45 pe55
sdpvar v1 v2 v3 v4 v5 v6 v7 v8 v9
sdpvar u1 u2 u3 u4 u5 u6 u7 u8 u9


%Inv=pc0+pc1*x+pc2*k+pc3*m+pc4*y+pc5*n+pc11*x^2+pc12*x*k+pc13*x*m+pc14*x*y+pc15*x*n+pc22*k^2+pc23*m*k+pc24*k*y+pc25*k*n
%+pc33*m^2+pc34*m*y+pc35*m*n+pc44*y*y+pc45*y*n+pc55*n*n
%Inv2=pd0+pd1*x+pd2*k+pd3*m+pd4*y+pd5*n+pd11*x^2+pd12*x*k+pd13*x*m+pd14*x*y+pd15*x*n+pd22*k^2+pd23*m*k+pd24*k*y+pd25*k*n
%+pd33*m^2+pd34*m*y+pd35*m*n+pd44*y*y+pd45*y*n+pd55*n*n
%preEinn = pe0+pe1*x+pe2*k+pe3*m+pe4*y+pe5*n+pe11*x^2+pe12*x*k+pe13*x*m+pe14*x*y+pe15*x*n+pe22*k^2+pe23*m*k+pe24*k*y+pe25*k*n
%+pe33*m^2+pe34*m*y+pe35*m*n+pe44*y*y+pe45*y*n+pe55*n*n
FeasibilityConstraints=[sos(pc0+pc1*x+pc2*k+pc3*m+pc4*y+pc5*n+pc11*x^2+pc12*x*k+pc13*x*m+pc14*x*y+pc15*x*n+pc22*k^2+pc23*m*k+pc24*k*y+pc25*k*n+pc33*m^2+pc34*m*y+pc35*m*n+pc44*y*y+pc45*y*n+pc55*n*n-20*(m-x)-v1*k-v2*x*k-v3*k*m-u5*m), ...
    sos(k-(pc0+pc1*x+pc2*k+pc3*m+pc4*y+pc5*n+pc11*x^2+pc12*x*k+pc13*x*m+pc14*x*y+pc15*x*n+pc22*k^2+pc23*m*k+pc24*k*y+pc25*k*n+pc33*m^2+pc34*m*y+pc35*m*n+pc44*y*y+pc45*y*n+pc55*n*n)-v4*(m-x)-u4*m), ...
    sos((pe0+pe1*x+pe2*k+pe3*m+pe4*y+pe5*n+pe11*x^2+pe12*x*k+pe13*x*m+pe14*x*y+pe15*x*n+pe22*k^2+pe23*m*k+pe24*k*y+pe25*k*n+pe33*m^2+pe34*m*y+pe35*m*n+pe44*y*y+pe45*y*n+pe55*n*n)-...
    (pc0+pc1*x+pc2*k+pc3*m+pc4*y+pc5*n+pc11*x^2+pc12*x*k+pc13*x*m+pc14*x*y+pc15*x*n+pc22*k^2+pc23*m*k+pc24*k*y+pc25*k*n+pc33*m^2+pc34*m*y+pc35*m*n+pc44*y*y+pc45*y*n+pc55*n*n)-u2*m-u3*(m-x)), ...
    sos((pe0+pe1*x+pe2*k+pe3*m+pe4*y+pe5*n+pe11*x^2+pe12*x*k+pe13*x*m+pe14*x*y+pe15*x*n+pe22*k^2+pe23*m*k+pe24*k*y+pe25*k*n+pe33*m^2+pe34*m*y+pe35*m*n+pe44*y*y+pe45*y*n+pe55*n*n)-...
    (pc0+pc1*x+pc2*k+pc3*m+pc4*y+pc5*n+pc11*x^2+pc12*x*k+pc13*x*m+pc14*x*y+pc15*x*n+pc22*k^2+pc23*m*k+pc24*k*y+pc25*k*n+pc33*m^2+pc34*m*y+pc35*m*n+pc44*y*y+pc45*y*n+pc55*n*n)),...
    sos((pd0+pd1*x+pd2*k+pd3*m+pd4*y+pd5*n+pd11*x^2+pd12*x*k+pd13*x*m+pd14*x*y+pd15*x*n+pd22*k^2+pd23*m*k+pd24*k*y+pd25*k*n+pd33*m^2+pd34*m*y+pd35*m*n+pd44*y*y+pd45*y*n+pd55*n*n)-...
    (pc0+pc1*(x+0.05)+pc2*(k+1)+pc3*m+pc4*n+pc5*n+pc11*(x+0.05)^2+pc12*(x+0.05)*(k+1)+pc13*(x+0.05)*m+pc14*(x+0.05)*n+pc15*(x+0.05)*n+pc22*(k+1)^2+pc23*m*k+pc24*(k+1)*n+pc25*(k+1)*n+pc33*m^2+pc34*m*n+pc35*m*n+pc44*n*n+pc45*n*n+pc55*n*n)-v5*y-v6*y*n-u6*n), ...
    sos((pc0+pc1*x+pc2*k+pc3*m+pc4*0+pc5*n+pc11*x^2+pc12*x*k+pc13*x*m+pc14*x*0+pc15*x*n+pc22*k^2+pc23*m*k+pc24*k*0+pc25*k*n+pc33*m^2+pc34*m*0+pc35*m*n+pc44*0*0+pc45*0*n+pc55*n*n)-...
    (pd0+pd1*x+pd2*k+pd3*m+pd4*y+pd5*n+pd11*x^2+pd12*x*k+pd13*x*m+pd14*x*y+pd15*x*n+pd22*k^2+pd23*m*k+pd24*k*y+pd25*k*n+pd33*m^2+pd34*m*y+pd35*m*n+pd44*y*y+pd45*y*n+pd55*n*n)-v7*(n-y)-v8*(n-y)*n-v9*(n-y)*y-u7*n),...
    sos((pd0+pd1*x+pd2*k+pd3*m+pd4*(y+0.05)+pd5*n+pd11*x^2+pd12*x*k+pd13*x*m+pd14*x*(y+0.05)+pd15*x*n+pd22*k^2+pd23*m*k+pd24*k*(y+0.05)+pd25*k*n+pd33*m^2+pd34*m*(y+0.05)+pd35*m*n+pd44*(y+0.05)*(y+0.05)+pd45*(y+0.05)*n+pd55*n*n)-...
    (pd0+pd1*x+pd2*k+pd3*m+pd4*y+pd5*n+pd11*x^2+pd12*x*k+pd13*x*m+pd14*x*y+pd15*x*n+pd22*k^2+pd23*m*k+pd24*k*y+pd25*k*n+pd33*m^2+pd34*m*y+pd35*m*n+pd44*y*y+pd45*y*n+pd55*n*n)-u8*(n-y)-u9*n),...
    u1>=0, u2>=0, u3>=0, u4>=0, u5>=0, u6>=0, u7>=0, u8>=0, u9>=0];

options=sdpsettings('verbose',1,'solver','sedumi');

solvesos(FeasibilityConstraints,[ ],options,[pc0; pc1; pc2; pc3; pc4; pc5; pc11; pc12; pc13; pc14; pc15; pc22; pc23; pc24; pc25; pc33; pc34; pc35; pc44; pc45; pc55;...
    pd0; pd1; pd2; pd3; pd4; pd5; pd11; pd12; pd13; pd14; pd15; pd22; pd23; pd24; pd25; pd33; pd34; pd35; pd44; pd45; pd55;...
    pe0; pe1; pe2; pe3; pe4; pe5; pe11; pe12; pe13; pe14; pe15; pe22; pe23; pe24; pe25; pe33; pe34; pe35; pe44; pe45; pe55;
    v1; v2; v3; v4; v5; v6; v7; v8; v9; u1; u2; u3; u4; u5; u6; u7; u8; u9])

vpc0=value(pc0)
vpc1=value(pc1)
vpc2=value(pc2)
vpc3=value(pc3)
vpc11=value(pc11)
vpc12=value(pc12)
vpc13=value(pc13)
vpc22=value(pc22)
vpc23=value(pc23)
vpc33=value(pc33)
L=vpc0+vpc1*x+vpc2*k+vpc3*m+vpc11*x^2+vpc12*x*k+vpc13*x*m+vpc22*k^2+vpc23*m*k+vpc33*m^2;


