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

## Help for Jacobian calculations Beginner
240 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
240 Views
>> the RCI_Request == -5. ??
This is an unexpected output task status - please show a complete example demonstrating the claimed bug. Beginner
240 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
240 Views
well, sounds better :), but somehow I don't see when thenext iteration values are different with both the cases. Beginner
240 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 