Exchanging information between function, gradient, and Hessian

Hello there

I have a constrained nonlinear optimzation problem. Denote the objective, constraints, and their respective gradient vectors and Hessian matricies as

.

The Hessian of the Lagrangian function is .

The gradient vectors of the objective and constraints functions are given in the (conditionalized) output as

function [f,df] = myObj(x)
f; %compute f
    if nargout > 1
    df; %compute df
    end
end

function [g,dg] = myConst(x)
g(k); % compute g_k
   if nargout > 1
   dg(k); % compute dg_k    
   end
end

NOTE: The gradient vectors are only calculated when fmincon needs to evaluate them

Fmicon requests that the Hessian is supplied via a seperate function,

function [H] = myHess(x,lam)
H; %compute H
end

However, in the case of my problem, the computation of H requires the gradient vectors Obviously, I don't want to recompute them inside myHess. Instead, I want to pass the information to myHess as

function [H] = myHess(x,lam,df,dg)
H; %compute H using df and dg
end

My Question: How can I efficiently use the results of df and dg inside myHess?

I have seen the page on using nested functions: https://ch.mathworks.com/help/optim/ug/objective-and-nonlinear-constraints-in-the-same-function.html. However, here they did not show what to do when I have conditionalized outputs.

Thanks a lot!

Answer 1

[Marty]

See my reply to Sai

Answer 2

[Sai Kiran]

Hi,

As per my understanding you dont want to recompute the df,dg. You can directly pass the functions myObj & myConst in the myHess function to get the df, dg .

Please follow the below example for better understanding.

myHess function uses the func1 and func2 outputs.

function [H] = myHess(a,b)
    H=func1(a)+func2(b);
end

function [a]=func1(x)
a=x*x;
end

function [b]=func2(y)
b=y*y;
end

Name the matlab file as myHess.m and use the myHess function.

I hope it helps!

Thanks.

[Marty]

Hi Sai

Thanks for your reply!

Mmm, this is not exactly what I was looking for---unless I am not understanding your code correctly. The computations of a and b may already be done, so it would be redudant to do them inside of myHess. The Hessian is not called at every iteration.

However, I did solve this problem using nested functions. I expanded the example shown here: https://ch.mathworks.com/help/optim/ug/objective-and-nonlinear-constraints-in-the-same-function.html. Below, you can find an example of what I did. NOTE: I had to remove and re-write some parts, since it contained confidential material, so code may be incomplete.

function [x,f_out,flag,output]  = solve_problem(x0)

%% nested variables (shared between functions)  
x_last_obj = []; % last iteration used to compute objective
x_last_nonlincons = []; % last iteration used to compute nonlinear constraints

f_obj = []; % value of objective
df_obj = []; % gradient of objective

c_nonlin_ineq = []; % values of inequality constraints
dc_nonlin_ineq = []; % gradient ofinequality constraints
    
%% handles for fmincon 
obj = @(x) my_obj(x)
nonlin = @(x) my_constr(x)
HessFcn = @(x,lam) my_hessian(x,lam)

%% options for fmincon 
options = optimoptions('fmincon','Algorithm','interior-point',...
        "SpecifyConstraintGradient",true,"SpecifyObjectiveGradient",true,...
        'HessianFcn',@HessFcn);

%% solve problem
[x,f_out,flag,output] = fmincon(obj,x0,[],[],[],[],[],[],nonlin,options);
      
%% Nested functions 
    function [f,df] = my_obj(x)      
        if nargout > 1 % need gradient
            [f,df] =  objective_function(x);
            f_obj = f;
            df_obj = df;
        else % need only value
            f =  objective_function(x);
            f_obj = f;
            df_obj = [];
        end
         x_last_obj = x;
    end


    function [c,ceq,dc,deq] = my_constr(x)
        if nargout > 2 % need gradient
            [c,~,dc,~]  = nonlinear_constraints(x);
            c_nonlin_ineq = c;
            dc_nonlin_ineq = dc;
        else % need only value
            c  =  nonlinear_constraints(x);
            c_nonlin_ineq = c;
            dc_nonlin_ineq = [];     
        end
        ceq = [];
        deq = [];
        x_last_nonlincons = x;
    end

    function H = my_hessian(x,lam)
        % check if gradients have been calculated and if the values are up
        % to date; if false, calculate them
        if isempty(dc_nonlin_ineq) ||  isempty(df_obj) ||  ~isequal(x_last_obj,x) || ~isequal(x_last_nonlincons,x)
           [c_nonlin_ineq,~,dc_nonlin_ineq,~]  = nonlinear_constraints(x);
           [f_obj,df_obj] =  objective_function(x);
        end
        
        H = hessian(x,lam,f_obj,df_obj,c_nonlin_ineq,dc_nonlin_ineq);
       
    end
  
end