Using Matlab coder with lsqcurvefit(), can not pass optimoptions as an argument to the target function

I want to convert my function to mex using matlab coder. The function fits data to an exponential model.

function ab_fitted = my_fit_ab_1(x0, xdata, ydata)
simple_exp_func = @(x, xdata) x(1)*exp(-x(2)*xdata);   % y=a*exp(-b)
Opt = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt','MaxIter', 100);
ab_fitted = lsqcurvefit(simple_exp_func, x0, xdata, ydata, [],[], Opt);
end

The script main_1.m that calls my_fit_ab_1() is as follows,

clear

a = rand();
b = rand();
xdata = rand(8,1);
ydata = a*exp(-b*xdata);

x0 = [rand, rand];
ab_fitted = my_fit_ab_1(x0, xdata, ydata);

With matlab coder App, I can easily convert my_fit_ab_1() to mex file my_fit_ab_1_mex.mexw64.

However, if I want to pass the variable Opt as an argument, i.e. my_fit_ab_2.m

function ab_fitted = my_fit_ab_2(x0, xdata, ydata, Opt)
simple_exp_func = @(x, xdata) x(1)*exp(-x(2)*xdata);   % y=a*exp(-b)
ab_fitted = lsqcurvefit(simple_exp_func, x0, xdata, ydata, [],[], Opt);
end

and main_2.m

clear

a = rand();
b = rand();
xdata = rand(8,1);
ydata = a*exp(-b*xdata);
Opt = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt','MaxIter', 100);

x0 = [rand, rand];
ab_fitted = my_fit_ab_2(x0, xdata, ydata, Opt);

the coder says " The 'optim.options.Lsqcurvefit' class does not support code generation. "

Why optimoptions can be set inside the function while can not be passed as an argument?

I've attached the four files mentioned above.

Answer 1

[Steve Grikschat]

Hi Xingwang,

Despite support for the optimoptions function, the objects it generates are not supported for code generation. In the generated code, optimoptions generates structs of options.

As a workaround, I can suggest 2 alternatives, both of which require you to call optimoptions in the function that is being codegened:

1. The most simple, which will work for this case but won't scale well is to pass the individual options you want set as inputs to your function. This won't be as easy to update as you add more options.

function ab_fitted = my_fit_ab_2(x0, xdata, ydata, maxiter)
Opt = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt','MaxIterations',maxiter);
simple_exp_func = @(x, xdata) x(1)*exp(-x(2)*xdata);   % y=a*exp(-b)
ab_fitted = lsqcurvefit(simple_exp_func, x0, xdata, ydata, [],[], Opt);
end

2. If you anticipate wanting to toggle different options in the future, then try passing a struct of options that can be copied over into the optimoptions:

function ab_fitted = my_fit_ab_2(x0, xdata, ydata, optIn)
f = fieldnames(optIn);
Opt = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');
for k = 1:numel(f)
    thisOption = f{k};
    Opt.(thisOption) = optIn.(thisOption);
end
simple_exp_func = @(x, xdata) x(1)*exp(-x(2)*xdata);   % y=a*exp(-b)
ab_fitted = lsqcurvefit(simple_exp_func, x0, xdata, ydata, [],[], Opt);
end

To generate code for this, try this script:

dblMatType = coder.typeof(double(1), [8, 1], [1, 1]);
dblVecType = coder.typeof(double(1), [2, 1], [1, 1]);
dblScalarType = coder.typeof(double(1), [1, 1], [1, 1]);
optsType = struct('Algorithm', 'levenberg-marquardt', 'MaxIterations', dblScalarType);
cfg = coder.config('mex');

%            x0          xdata      ydata      opts
cArgs = {dblVecType, dblMatType, dblMatType, optsType};
codegen -config cfg my_fit_ab_2 -args cArgs

Note that the 4th input is a struct of "typeof" objects that you can add to as needed.

[Xingwang Yong]

The 2nd solution seems useful, I'll try it. Thank you Steve!

[Xingwang Yong]

I tried, the 2nd solution does not work on Matlab2020b.

With the error "??? Error calling 'mtimes'. This call-site expects more outputs than this function can supply.

Error in ==> my_fit_ab_2 Line: 8 Column: 19"

[Steve Grikschat]

It turns out that MATLAB Coder doesn't like the anonymous function.

Make simple_exp_func into a sub-function and give lsqcurvefit a handle to it:

... % <rest of the code above goes here....>

ab_fitted = lsqcurvefit(@simple_exp_func, x0, xdata, ydata, [],[], Opt);

end % end of main function my_fit_ab_2

function F = simple_exp_func(x, xdata)

F = x(1)*exp(-x(2)*xdata); % y=a*exp(-b)

end

[Xingwang Yong]

Thanks Steve.

But this is not only related to anonymous function.

Sometimes we have to use anonymous functions, e.g. when we want to pass extra parameters to function in optimization. Here is a simple demo that can not use codegen.

function x_fit = my_fit_ab_3(x0, xdata, ydata, constant, Opt)
obj_fun = @(x,xdata) simple_exp_func_1(x, xdata, constant);

f = fieldnames(Opt);
options = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt'); % for codegen, must use levenberg-marquardt algorithm
% comment the for-loop, then this can be codegen
for k = 1:numel(f)
    thisOption = f{k};
    options.(thisOption) = Opt.(thisOption);
end

x_fit = lsqcurvefit(obj_fun,x0,xdata,ydata,[],[],options);
end


function y = simple_exp_func_1(x, xdata, constant)
y = x(1)*exp(-x(2)*xdata) + constant;
end

And the generation scripts like

%% Create configuration object of class 'coder.MexCodeConfig'.
cfg = coder.config('mex');
cfg.GenerateReport = true;
cfg.ReportPotentialDifferences = false;

%% Define argument types for entry-point 'my_fit_ab_3'.
ARGS = cell(1,1);
ARGS{1} = cell(5,1);
ARGS{1}{1} = coder.typeof(0,[1 2]);
ARGS{1}{2} = coder.typeof(0,[8 20],[1 1]);
ARGS{1}{3} = coder.typeof(0,[8 20],[1 1]);
ARGS{1}{4} = coder.typeof(0);
ARGS_1_5 = struct;
ARGS_1_5.Display = coder.typeof('X',[1 Inf],[0 1]);
ARGS_1_5.MaxIter = coder.typeof(0);
ARGS_1_5.FunctionTolerance = coder.typeof(0);
ARGS{1}{5} = coder.typeof(ARGS_1_5);

%% Invoke MATLAB Coder.
codegen -config cfg my_fit_ab_3 -args ARGS{1}

So the codegen would fail if we used a strcut to construct the options and used anonymous functions simultaneously.

Probably this peculiar behavior should be listed in the doc.

[Steve Grikschat]

Hi Xingwang,

Thanks for the follow up.

Have you been able to successfully get something working for the case above?

Another way to pass a constant or fixed parameter to a function is via the use of "nested" functions, which can share the workspace of the main/calling function. Perhaps that will work better?

I think some of the issue has to do with what Coder can infer from the code. My guess is that the anonymous function part of it. I'll take a look and see what I can come up with.

[Xingwang Yong]

I can codegen if I comment the for-loop, which converts the option struct to optimal options. In this way, although I can not change the options conveniently, it works at least.

Yes, nested function is another way to pass extra parameters, but it is less intuitive. For now, I prefer to use anonymous function.

Steve, thanks for you continuous help.

[Steve Grikschat]

I think there may be a bug in the codegen lsqcurvefit that prevents this from working (consistently) with anonymous functions.

I believe that the nested function approach might be the best workaround for this case until the bug is fixed. I've tested with nested functions and it works, even with the for-loop setting options.

There may be other workarounds, but I haven't been able to find them.

[Xingwang Yong]

Thanks Steve. I'll wait until it is fixed.

[Steve Grikschat]

Hi Xingwang,

It turns out that this may not be fixable in our code. The for-loop that sets options in a dynamic way breaks the codegen process for lsqcurvefit.

The for-loop that sets options in a dynamic way means that Coder cannot determine the value of SpecifyObjectiveGradient statically. This is a problem for generating the correct code since lsqcurvefit relies on this value to determine how to call the objective function (as we do in MATLAB).

Therefore, to deal with the uncertainty at codegen time, codegen wants to generate C code for a function with 2 output arguments (y and Jacobian) which is the maximum number of possible outputs.

Anyway, a workaround to add a Jacobian output with a dummy value. This works since Coder can safely make a 2 output syntax, and then the value of SpecifyObjectiveGradient = false will be checked at runtime to avoid ever using the dummy Jacobian value.

The modified simple_exp_func_1 looks like this

function [y,J] = simple_exp_func_1(x, xdata, constant)
y = x(1)*exp(-x(2)*xdata) + constant;
if nargout > 1
    J = 0;  % Dummy value, incorrect size, but doesn't matter since it never gets used at runtime
    % Real Jacobian looks like this BTW 
    %  J = [exp(-x(2)*xdata) -xdata*x(1).*exp(-x(2)*xdata)];
end 
end % simple_exp_func_1

[Xingwang Yong]

Steve, thanks for your expertise. This dummy Jacobian solution works perfectly.