%% Optimize a Nonsmooth Function
% This example shows the importance of choosing an appropriate solver. It also shows that a single point of non-smoothness can cause problems for Optimization Toolbox™ solvers.
%% A Function with a Single Nonsmooth Point
% The function f(x) = ||x||^{1/2} is nonsmooth only at the point 0, which is the minimizing point. Here is a 2-D plot using the Euclidean norm.

figure
fsurf(@(x,y)sqrt(x.^2+y.^2).^(1/2))
xlabel('x')
ylabel('y')
title('Norm([x,y])^{1/2}')

% Instead of Euclidean norm, this example uses matrix norm for a 2-by-6 matrix x. The matrix norm relates to the singular value decomposition, which is not as smooth as the Euclidean norm. See 2-Norm of Matrix.
%
%% Minimize Using fminsearch
% Start fminsearch from a nonzero 2-by-6 matrix x0, and attempt to locate the minimum of f. The documentation states that fminsearch sometimes can handle discontinuities, so this is a reasonable first try.

fun = @(x)norm([x(1:6);x(7:12)])^(1/2);
x0 = [1:6;7:12];
rng default
x0 = x0 + rand(size(x0))

[xfms,fvalfms,eflagfms,outputfms] = fminsearch(fun,x0)

% fminsearch runs out of function evaluations before it converges to a solution. The reported solution is poor.
%
%% Use patternsearch
% patternsearch is the recommended first solver to try for nonsmooth problems. See Table for Choosing a Solver.

[xps,fvalps,eflagps,outputps] = patternsearch(fun,x0)

% patternsearch reaches a good solution, but takes over 10,000 function evaluations to do so.
%
%% Use particleswarm
% particleswarm is recommended as the next solver to try. See Choosing Between Solvers for Nonsmooth Problems.

[xpsw,fvalpsw,eflagpsw,outputpsw] = particleswarm(fun,12)

% particleswarm finds an even more accurate solution than patternsearch, but takes over 35,000 function evaluations.
%
%% Use ga
% ga is a popular solver, but is not recommended as the first solver to try. See how well it works on this problem.

[xga,fvalga,eflagga,outputga] = ga(fun,12)

% ga does not find as good a solution as patternsearch or particleswarm, and takes about twice as many function evaluations as particleswarm.
%
%% Use fminunc from Optimization Toolbox
% fminunc is not recommended for nonsmooth functions. See how well it does on this one.

[xfmu,fvalfmu,eflagfmu,outputfmu] = fminunc(fun,x0)

% The fminunc solution is not as good as the ga solution. However, fminunc reaches the rather poor solution in relatively few function evaluations.
%
%% Use fmincon from Optimization Toolbox
% fmincon can sometimes minimize nonsmooth functions. See how well it works on this one.

[xfmc,fvalfmc,eflagfmc,outputfmc] = fmincon(fun,x0)

% fmincon with default options produces an accurate solution after fewer than 1000 function evaluations.
%
%% Summary of Results
% Choosing the appropriate solver leads to better, faster results. This summary shows how disparate the results can be.

Solver = {'fminsearch';'patternsearch';'particleswarm';'ga';'fminunc';'fmincon'};
SolutionQuality = {'Poor';'Good';'Good';'Mediocre';'Poor';'Good'};
FVal = [fvalfms,fvalps,fvalpsw,fvalga,fvalfmu,fvalfmc]';
NumEval = [outputfms.funcCount,outputps.funccount,outputpsw.funccount,...
    outputga.funccount,outputfmu.funcCount,outputfmc.funcCount]';
results = table(Solver,SolutionQuality,FVal,NumEval)

% Copyright 2012 The MathWorks, Inc.