Having problems solving for critical values (Multi variable function)

Greetings, I'm trying to manually solve for critical points on a function I have, but Matlab seems to time out everytime it tries to solve. Here's my code:

if true
  % syms f(x,y)
f(x,y)= (100/((3*x+13)^2+(3*y+2)^2+(13+2))) - ((2*100)/((3*x+5)^2+(3*y-5)^2+10)-((3*25)/((3*x+13)^2+(3*y-5)^2+18))+((4*100)/(3*x-5)^2+(3*y-2)^2+7));

%ezsurf(f,[-10,10]) colormap(jet) hold off

dx= diff(f,x); display(dx);

dy= diff(f,y); display(dy);

dx = dx==0; dy= dy==0;

sol = solve([dx, dy], [x, y]); xSol = sol.x; ySol = sol.y;

end

Answer 1

[Stephan]

Hi,

my machine also doesnt stop being busy, when trying to solve this problem symbollically, so here is an numeric approach:

syms f(x,y)
fun = f == (100/((3*x+13)^2+(3*y+2)^2+(13+2))) - ((2*100)/((3*x+5)^2+(3*y-5)^2+10)-((3*25)/((3*x+13)^2+(3*y-5)^2+18))+((4*100)/(3*x-5)^2+(3*y-2)^2+7)); 

fsurf(rhs(fun),[-10,10]) 
colormap(jet) 
hold on

% define the function f as function handle (same with the derivates at the
% function below. x0 is a guess from looking at the function plot.
func = @(x)-(1.0./(x(1).*3.0-5.0).^2.*-4.0e2-(x(2).*3.0-2.0).^2-2.0e2./((x(1).*3.0+5.0).^2+(x(2).*3.0-5.0).^2+1.0e1)+1.0e2./((x(1).*3.0+1.3e1).^2+(x(2).*3.0+2.0).^2+1.5e1)+7.5e1./((x(1).*3.0+1.3e1).^2+(x(2).*3.0-5.0).^2+1.8e1)-7.0);
sol = fmincon(func,[-5 0.5],[],[],[],[],[-10 -10], [10 10], @derivates)

% show the solution on the plot
max_value = -func(sol);
scatter3(sol(1), sol(2), max_value, 'og', 'LineWidth',3)
hold off

% Make sure that the derivates of f with respect to x and y
% are equal to zero
function [c,ceq] = derivates(x)
c = [];
ceq = [(x(1).*1.8e1+3.0e1).*1.0./((x(1).*3.0+5.0).^2+(x(2).*3.0-5.0).^2+1.0e1).^2.*2.0e2-(x(1).*1.8e1+7.8e1).*1.0./((x(1).*3.0+1.3e1).^2+(x(2).*3.0+2.0).^2+1.5e1).^2.*1.0e2-(x(1).*1.8e1+7.8e1).*1.0./((x(1).*3.0+1.3e1).^2+(x(2).*3.0-5.0).^2+1.8e1).^2.*7.5e1+1.0./(x(1).*3.0-5.0).^3.*2.4e3, x(2).*-1.8e1-(x(2).*1.8e1+1.2e1).*1.0./((x(1).*3.0+1.3e1).^2+(x(2).*3.0+2.0).^2+1.5e1).^2.*1.0e2+(x(2).*1.8e1-3.0e1).*1.0./((x(1).*3.0+5.0).^2+(x(2).*3.0-5.0).^2+1.0e1).^2.*2.0e2-(x(2).*1.8e1-3.0e1).*1.0./((x(1).*3.0+1.3e1).^2+(x(2).*3.0-5.0).^2+1.8e1).^2.*7.5e1+1.2e1]; 
end

By using this i got:

sol =

    -4.727519221259621   0.605820153518793

which gives this as result when calculating the function value:

>> max_value

max_value =

     -4.266365719653375

or checked against the result fro your symbolic function:

>> check = double(subs(rhs(fun),[x y],sol))

check =

    -4.266365719653376

As far as in the plot can be seen, this should be (one) correct solution to your function. Note that it could be a local extremum.

Best regards

Stephan