Computing Partial derivates of anonymous functions comprised of sum of standard functions

Hi there,

I am implementing a conjugate gradient based optimisation algorithm, and to achieve this i need to calculate the grad of the input function. I currently pass in an anonymous function, which is comprised as the sum of other anonymous functions, of which those are comprised of standard matlab functions. Here is how the objective function is defined

% one constraint for example - all 6 are similar form
c2 = @(x) mat.minSF_yield - mat.yield ...
        / axial_stress(mat.F, mat.E, x(1), x(2));


transf_obj = @(x) obj(x) - r * (1 / (c1(x) + c2(x) + c3(x) + c4(x) + c5(x) + c6(x)));

the partial derivative function calculates the partial derivative anon functions and adds them into a cell. This allows the partials to be calculated before the optimisation loop saving computation time

function [grad_cell] = partials(func, x)
    % Get dimensionality of inputs
    n = numel(x);
    
    % Convert the anonymous function to a symbolic expression
    syms a [1 n]
    funcSym = func(a);
    
    % Initialize the array to hold partial derivative functions
    grad_vect = zeros(1, n);
    
    grad_cell = cell(1,n);
    % Compute partial derivatives
    for i = 1:n
        % Compute symbolic partial derivative with respect to the i-th variable
        partialDerivSym = diff(funcSym, a(i));
        
        % Convert symbolic partial derivative back to an anon function
        grad = matlabFunction(partialDerivSym, 'Vars', {a});
        
        % Add anon function to cell
        grad_cell{i} = grad;
    end
end

To evalute the gradient at a specific point the following function is used

function [nabla] = compute_grad(cell, x)
    nabla = zeros(1,numel(x))
    for i =1:numel(x)
        nabla(i) = cell{i}(x);
    end
end

When evaluting any partial, for some values such as x = [14,3] something within the feasable domain, an empty array is returned. And for values such as x=[7.9,7], a value close to the results i was getting from using fmincon an error of the following is produced.

Error using diff
Difference order N must be a positive integer scalar.
Error in sym/matlabFunction>@(in1)(in1(:,1).*pi.*(3.0./2.5e+...

My current theory is that the conversion to and from symbolic variables to compute the derivative is not being achieved successfuly given the combination of anon functions and locally defined functions. I have tested the partials function with simple functions directly defined and it appears to work. This program is split across 3 files, the first which defines variables and calls the conj-grad m-file. This file then calls the partials.m file which is just the partials function.

%% Variable Definition ...
...
%% Call 'Conj_grad.m'
    %% Call 'partials.m'
        %% Returns cell
    %% compute at specific point, indexing from cell using compute_grad function
        %% Error

Any help is much appreciated.

Answer 1

[Torsten]

Stay numerical in your computation and approximate the gradient of your objective function "obj" by finite-difference quotients:

compute_gradient(@objective,[3 5 -1], 1e-8)

ans = 3x1

2.0000 0.2837 -60.2566

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

[2 ;cos(5); -3*exp(-3*(-1))]

ans = 3x1

2.0000 0.2837 -60.2566

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

function value = objective(x)
  value = 2*x(1) + sin(x(2)) + exp(-3*x(3));
end

function gradobj = compute_gradient(obj,x,h)
  n = numel(x);
  gradobj = zeros(n,1);
  f = obj(x);
  for i = 1:n
    xh = x(i);
    x(i) = x(i) + h;
    fh = obj(x);
    gradobj(i) = (fh-f)/h;
    x(i) = xh;
  end  
end

[Paul]

Not disagreeing with @Torsten's recommendation.

Just showing that the OPs approach may be feasible. Can't say for sure w/o having the full code in order to recreate the problem. Could just as easily have used gradient, but it seems like the jacobian is probably what's really needed.

J = matlabFunction(jacobian(objective(sym('x',[3 1]))),'Var',sym('x',[3 1]))

J = function_handle with value:
    @(x1,x2,x3)[2.0,cos(x2),exp(x3.*-3.0).*-3.0]

J(3,5,-1)

ans = 1x3

2.0000 0.2837 -60.2566

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

function value = objective(x)
  value = 2*x(1) + sin(x(2)) + exp(-3*x(3));
end

Answer 2

[Torsten]

Your code works for me. Maybe you have a MATLAB version where the second argument to "diff" is only interpreted as the order of the derivative instead of the differentiation variable ?

func = @(x) 2*x(1) + sin(x(2)) + exp(-3*x(3));
x = [3 5 -1];
grad_cell = partials(func,x)

grad_cell = 1x3 cell array
    {@(in1)2.0}    {@(in1)cos(in1(:,2))}    {@(in1)exp(in1(:,3).*-3.0).*-3.0}

compute_grad(grad_cell,x)

ans = 1x3

2.0000 0.2837 -60.2566

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

function [grad_cell] = partials(func, x)
    % Get dimensionality of inputs
    n = numel(x);
    
    % Convert the anonymous function to a symbolic expression
    syms a [1 n]
    funcSym = func(a);
    
    % Initialize the array to hold partial derivative functions
    grad_vect = zeros(1, n);
    
    grad_cell = cell(1,n);
    % Compute partial derivatives
    for i = 1:n
        % Compute symbolic partial derivative with respect to the i-th variable
        partialDerivSym = diff(funcSym, a(i));
        
        % Convert symbolic partial derivative back to an anon function
        grad = matlabFunction(partialDerivSym, 'Vars', {a});
        
        % Add anon function to cell
        grad_cell{i} = grad;
    end
end

function [nabla] = compute_grad(cell, x)
    nabla = zeros(1,numel(x));
    for i =1:numel(x)
        nabla(i) = cell{i}(x);
    end
end