In this example the Hessian H is calculated and given to fmincon:

x0 = [-2; 1];

function [f, gf, H] = objf(x)
   r = x(1).^2 + x(2).^2;
   s = exp(-r);

   f = x(1).*s + r/20;
   gf = [(1-2.*x(1).^2).*s+x(1)./10; -2.*x(1).*x(2).*s+x(2)/10];
   H = [2*x(1)*(2*x(1)^2-3)*s+1/10, 2*x(2)*(2*x(1)^2-1)*s;
        2*x(2)*(2*x(1)^2-1)*s,      2*x(1)*(2*x(2)^2-1)*s+1/10];
end

function [c,ceq,gc,gceq] = constraints(x)
   c = x(1).*x(2)/2 + (x(1)+2).^2 + (x(2)-2).^2/2 - 2;
   ceq = [];

   % Note that the gradients have to be row-vectors!!!
   gc = [x(2)/2+2*(x(1)+2), x(1)/2+x(2)-2];
   gceq = [];
end

opt = optimset('GradObj', 'on', 'GradConstr', 'on', 'HessianFcn', 'objective');

[x, fval, exitflag, output] = fmincon(@objf, x0, [], [], [], [], [], [], @constraints, opt)

I get the following error message:

error: structure has no member 'inverse_hessian'
error: called from
    __lm_feasible__ at line 49 column 9
    fmincon at line 457 column 24
    Test_optimization7 at line 29 column 29

If I add "hook.inverse_hessian = false;" after line 434 in fmincon.m

  ## passed function for hessian of objective function
  hook.hessian = f.hessian;
  hook.inverse_hessian = false;

it fixes the problem, but I'm not sure if this is really the correct way to solve it.

Another strange thing is that in contrast to Matlab, the constraint's gradient has to be a row vector. For compatibility it should also accept a column vector. Surprisingly both types are already accepted for the function gradient.

I'm using optim package version 1.6.2.