topic Nonlinear Least Squares Problem with Linear (Bound) Constraints in Intel® oneAPI Math Kernel Library

Nonlinear Least Squares Problem with Linear (Bound) Constraints

Shokriafra__Mohammad — Thu, 18 Jun 2020 08:25:48 GMT

Hi,

I noticed that the nonlinear solver uses central difference to calculate the Jacobi matrix, I was wondering is it possible to pass the analytical Jacobi elements to the code or not? I mean, in my problem, I want to use analytical Jacobi instead of numerical Jacobi. How can I do that?

Thanks
Javad

The nonlinear solver uses the

mecej4 — Fri, 19 Jun 2020 13:48:00 GMT

The nonlinear solver uses the reverse call interface. When RCI_REQUEST is equal to 2, you have to update the Jacobian matrix FJAC and call the solver again. The method used for calculating the Jacobian is entirely in your hands. You can enter formulae for the elements of the Jacobian in line, or in a separate subroutine that you write and call; you can use finite difference approximations for some or all elements, or obtain the derivatives from a finite-field extension using complex variables.

You can, for example, take the furnished testcase ex_nlsqp_f.f and modify it so that FJAC is calculated analytically instead of calling DJACOBI.

Re: The nonlinear solver uses the

Shokriafra__Mohammad — Tue, 30 Jun 2020 12:42:03 GMT

Thanks for your reply,

I am using this solver in C++, I used the "ex_nlsqp_bc_c.c" example and tried to implement what you mentioned as the following.

I defined the function as

extern void jacobian(MKL_INT *, MKL_INT *, double *, double *);

----

I called the Jacobian function as

if (RCI_Request == 2)
{
jacobian(&m, &n, x, fjac);
}

----

and defined the Jacobian like the following:

void jacobian(MKL_INT * m, MKL_INT * n, double *x, double *fj)
{
fj[0] = 2 * x[0];
fj[1] = 2 * x[1];
fj[2] = 1;
fj[3] = -1;

return;
}

----

The function was defined as:

void extended_powell(MKL_INT * m, MKL_INT * n, double *x, double *f)
{
f[0] = pow(x[0],2) + pow(x[1],2) - 9;
f[1] = x[0] - x[1] + 3;

return;
}

----

When I call the jacobian, unfortunately it does not work, but if I call djacobi it works.

I tracked the iteration of jacobian, it does around 20 iterations and becomes close to the solution, but I do not know why it fails after several iteration while it is getting closer to the correct solution.

The correct solution for this problem is "zero" and "3", and in the following I printed the two last iterations of the solution when I am using jacobian.

x1= 0.10554969128064530 x2= 2.9935339376683032

x1= 0.10554965354656229 x2= 2.9935339407828785

Can anyone help me?

Re: The nonlinear solver uses the

mecej4 — Wed, 01 Jul 2020 02:29:42 GMT

I wrote a detailed reply, but the forum software marked it as spam and removed it. Here is a short replacement:

Your Jacobian vector arranges the elements of the 2-dimensional matrix in row-major order. The MKL routine may expect it to be in column-major (as in Fortran) order. Try printing out the values of fjac from the finite-difference calculation and compare those with the analytically calculated values for accuracy and matching order.