clc
clear
syms x y z r1 r2 w
f=x^2+y^2+z^2+w^2;
g1=2*x-y+z-w-1;
g2=x+y-z+w-1;
h=f-r1*g1 -r2*g2;
hx=diff(h,x);
hy=diff(h,y);
hz=diff(h,z);
hw=diff(h,w);
hr1=diff(h,r1);
hr2=diff(h,r2);
r=solve([hx==0,hy==0,hz==0,hw==0,hr1==0,hr2==0],[x,y,z,w,r1,r2]);
arr_x=double(r.x)
arr_y=double(r.y)
arr_z=double(r.z)
arr_w=double(r.w)
arr_fv=[];
for i=1 : length(arr_x)
xv=arr_x(i);
yv=arr_y(i);
zv=arr_z(i);
wv=arr_w(i);
fv=subs(f,[x,y,z,w],[xv,yv,zv,wv]);
arr_fv(i)= double(fv);
g1v=subs(g1,[x,y,z,w],[xv,yv,zv,wv])
g2v=subs(g2,[x,y,z,w],[xv,yv,zv,wv])
end
arr_fv
disp('after sort:')
sort(arr_fv)
