g = @(x) ([x(1).^2+2*x(2).^2-5*x(1)+7*x(2)-40, 3*x(1).^2-x(2).^2+4*x(1)+2*x(2)-28]) fsolve(g, [2,3])