Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.

## Help for Jacobian calculations Beginner
146 Views
Hi,

I was trying to use trust region solver for one of my code. I found that the jacobian given by the solver is calculating only half of thevalues (i.e. n*m/2)whilefor rest of thecalculation it gives zeros and brakes out with the RCI_Request == -5. ??
I have tested solver with some small eqation where i have calculated the jacobian(i.e. (f(x) - f(x'))/(x - x')) and compared the same with jacobian given by the solver which are same but not for my original code.
One more thing is i found that, though i am giving initial guess as upper bounds ,while perturbing the variable in costfun calculation goes beyond upper bounds .??

I know without code it will be difficul to find the errors, hope someone might be able to hint the possible reasons.

followingare thecalls which i use. it may help

if (dtrnlspbc_init (&handle, &n, &m, x, LW, UP, eps, &iter1, &iter2,&rs) != TR_SUCCESS) {}
if (dtrnlspbc_solve (&handle, fvec, fjac, &RCI_Request) != TR_SUCCESS) {}
if (RCI_Request == 1)
{
CostFun(&m,&n,x,fvec);
}
if (RCI_Request == 2)
{
if (djacobi (CostFun, &n, &m, fjac, x, eps) != TR_SUCCESS)
{
return 0;
}
fp_temp= fopen("solverResFjac.dat","w");
for(i=0;i fprintf(fp_temp,"%f\\n",fjac);

Jacobian(fjac1,n,m,x,fvec); // user defined jacobian
fp_temp= fopen("solverResFjac1.dat","w");
for(i=0;i fprintf(fp_temp,"%f\\n",fjac1);
}

Thanks ,
Suresh

4 Replies Moderator
146 Views
>> the RCI_Request == -5. ??
This is an unexpected output task status - please show a complete example demonstrating the claimed bug. Beginner
146 Views

Thanks for the reply. I found there are someproblem with my code..:) .

But still if my initial guess is one of the boundary condition, in jacobian it crosses the boundary condition which i dont want to happen ,So i am calculating jacobian.

I found the values are same in solver calculated and user defined but still the next iteartion values are different with both the cases.

Please let me know if i am missing something here..

```[bash]#include
#include
#include
#include

void myJacobi(double *c,double *y,double *fjec,int n,int m,double *fvec);
void myCostFun(int *m, int *n, double *x, double *f);
double y;

int main()
{

//double y;
double xin;
double c = {3 ,2 ,3};

/* n - number of function variables
m - dimension of function value */
int n = 3, m = 7;
/* precisions for stop-criteria (see manual for more details) */
double eps;
/* solution vector. contains values x for f(x) */
double *x;
/* iter1 - maximum number of iterations
iter2 - maximum number of iterations of calculation of trial-step */
int iter1 = 1000, iter2 = 100;
/* initial step bound */
double rs = 100.0;
/* reverse communication interface parameter */
int RCI_Request;

/* controls of rci cycle */
int successful;
/* function (f(x)) value vector */
double *fvec;
/* jacobi matrix */
double *fjac,*fjac1;
/* lower and upper bounds */
double *LW, *UP;
/* number of iterations */
int iter;
/* number of stop-criterion */
int st_cr;
/* initial and final residuals */
double r1, r2;
/* TR solver handle */
_TRNSPBC_HANDLE_t handle;
/* cycles counter */
int i;

/* memory allocation */
x = (double*) malloc (sizeof (double)*n);
fvec = (double*) malloc (sizeof (double)*m);
LW = (double*) malloc (sizeof (double)*n);
UP = (double*) malloc (sizeof (double)*n);
fjac = (double*) malloc (sizeof (double)*m*n);
fjac1 = (double*) calloc (m*n,sizeof (double));
/* set precisions for stop-criteria */
for (i = 0; i < 6; i++)
{
eps  = 0.00001;
}

/* set the initial guess */
for (i = 0; i < n; i++)
{
x  = 4.2;
}
/* set the initial values */
for (i = 0; i < m; i++)
fvec  = 0.0;
for (i = 0; i < m*n; i++)
fjac  = 0.0;
/* set bounds */
for (i = 0; i < n; i++)
{
LW  = 0;
UP  = 4;
}

for(i=0;i<7;i++)
{
xin = i+1;
y = c*xin +c*xin*xin +c*xin*xin;
}

if (dtrnlspbc_init (&handle, &n, &m, x, LW, UP, eps, &iter1, &iter2,
&rs) != TR_SUCCESS)
{
printf ("| error in dtrnlspbc_initn");
return 0;
}

/* set initial rci cycle variables */
RCI_Request = 0;
successful = 0;
/* rci cycle */
while (successful == 0)
{
if (dtrnlspbc_solve (&handle, fvec, fjac, &RCI_Request) != TR_SUCCESS)
{
printf ("error in dtrnlspbc_solven");
return 0;
}
/* according to rci_request value we do next step */
if (RCI_Request == -1 ||
RCI_Request == -2 ||
RCI_Request == -3 ||
RCI_Request == -4 ||
RCI_Request == -5 ||
RCI_Request == -6)
/* exit rci cycle */
successful = 1;

if (RCI_Request == 1)
{
myCostFun(&m, &n, x, fvec);
printf("Cost function n");
for (i=0;i);
}
}
if (RCI_Request == 2)
{
/* compute jacobi matrix*/
if (djacobi(myCostFun, &n, &m, fjac, x, eps) != TR_SUCCESS)
{
printf ("error in djacobin");
return 0;
}

myJacobi(x,y,fjac1,n,m,fvec);
printf("Jacobiann MKLtt User definedn");
for (i=0;i,fjac1);
}
}
}

if (dtrnlspbc_get (&handle, &iter, &st_cr, &r1, &r2) != TR_SUCCESS)
{
return 0;
}

if (dtrnlspbc_delete (&handle) != TR_SUCCESS)
{
return 0;
}

printf ("The result isn");
for(i=0;i);
//printf ("n The Cost is %fn");

/* free allocated memory */
free (x);
free (fvec);
free (fjac);
free (LW);
free (UP);
/* if final residual is less than required precision then print pass */
if (r2 < 0.1)
printf ("| dtrnlspbc ..........PASSn");
/* else print failed */
else
printf ("| dtrnlspbc .........FAILEDn");

getch();

return 0;

}

void myCostFun(int *m, int *n, double *x, double *f) //  (double* c,double *y,double *fvec)
{
int  i= 0;
double fx;
double xin;

for(i=0;i<7;i++)
{
xin = i+1;
fx = x*xin +x*xin*xin +x*xin*xin;
f = (y-fx);
}
//return resCost;
}

void myJacobi(double *c1,double *y,double *fjec,int n,int m,double *fvec)
{
int i =0,j=0,k=0;
float x;
double fx;
double fv;
float var = 0.01;
float del;
double c;

for(j = 0;j = c1;
}
del = (c+var*c)- c;
c= c+var*c;
for(i=0;i = i+1;

fx = c*x +c*x*x +c*x*x;
fv = (y-fx);

fjec[j*m+i] = (fv - fvec)/del;

}
}

}[/bash]```

Thanks,
Suresh Moderator
146 Views
well, sounds better :), but somehow I don't see when thenext iteration values are different with both the cases.
please show these values which you received. Beginner
146 Views
Hi,

The output is as follows..

Cost function in 1st iteration -

-4.600000

-16.000000

-34.200000

-59.200000

-91.000000

-129.600000

-175.000000

Jacobian -

MKL User Defined

-1.000000 -1.000000

-2.000000 -2.000000

-3.000000 -3.000000

-4.000000 -4.000000

-5.000000 -5.000000

-6.000000 -6.000000

-7.000000 -7.000000

-1.000000 -1.000000

-4.000000 -4.000000

-9.000000 -9.000000

-16.000000 -16.000000

-25.000000 -25.000000

-36.000000 -36.000000

-49.000000 -49.000000

-1.000000 -1.000000

-4.000000 -4.000000

-9.000000 -9.000000

-16.000000 -16.000000

-25.000000 -25.000000

-36.000000 -36.000000

-49.000000 -49.000000

These are the cost function and jacobian for 1st iteration which are same..

But the cost function( new values of variable ) with djacobi is ...
Cost function in 2nd iteration

-0.022444

-0.036071

-0.040883

-0.036880

-0.024060

-0.002425

0.028026

whereas afterusing myjacobi the output is,

Cost function in 2nd iteartion

0.780000

3.120000

7.020000

12.480000

19.500000

28.080000

38.220000

if i feed the same data as cost and jacobian why there should be difference in next iteration?

Thanks,
Suresh 