[ADOL-C] sparse_jac return unreasonable NaN's

何旭 mailhexu at gmail.com
Sun Aug 28 04:27:15 EDT 2011


Hello,

I made some simple test on sparse_jac, and I found that if my function is
with condassign() and if any branch of condassign() returns a nan, the
sparse_jac returns a nan,too.  I think it there might be a bug in condassign
or sparse_jac.

My function is very simple. y=x/(exp(x)-1) .   However if x is too large
(x>100 or so) ,exp(x) results in nan,  if fabs(x) is too small
(fabs(x)<1.0e-3),  exp(x)-1 ->0.  So I wrote some conditional statement to
avoid nan's.   But that seems to be useless with sparse_jac.  In my test,
x=2000, which means the function return 0 and the jacobian be 0.   But that
does not happen.  That's result of test1.

I tried to substitute y=x/(exp(x)-1) with y=x, which never return nan.  And
the result is right.    (test 3)

So I guess there might be some logical error in sparse_jac .

I put the test results here.  test2 seems to be  of no use. Just ignore it.

Result of test1:
x= 20000:
0.000000e+00
dense:0.000000e+00
number of sparse matrix elements: 1
sparse: -nan

Result of test2:
x= 20000:
0.000000e+00
dense:0.000000e+00
number of sparse matrix elements: 0
sparse: -3.848143e-42

Result of test3:
x= 20000:
0.000000e+00
dense:0.000000e+00
number of sparse matrix elements: 1
sparse: 0.000000e+00



And here's my testing code.
---------------------------------------------------------------------------------------
#include <iostream>
#include <adolc/adolc.h>
#include <adolc/adolc_sparse.h>
using namespace std;

adouble myfunc(adouble x)
{
 adouble y;
 condassign(y,1.0e-03-fabs(x),1.0);
 condassign(y,x-100.0,0.0);
 condassign(y,(1.0e-3-fabs(x))*(x-100.0),x/(exp(x)-1.0));     //  TEST 1
//              NOTHING HERE
//TEST 2
//condassign(y,(1.0e-3-fabs(x))*(x-100.0),x);                         //TEST
3
 return y;
}


int main()
{
  adouble ax,ay;
  double x,y;
  x=200.0;
  trace_on(1);
  ax<<=x;
  ay=myfunc(ax);
  ay>>=y;
  trace_off();


  double **jac;
  jac=(double **)malloc(sizeof(double*));
  jac[0]=(double *)malloc(sizeof(double));

  x=20000.0;
  cout<<"x= "<<x<<":"<<endl;
  function(1,1,1,&x,&y);
  cout<<scientific<<y<<endl;


//   jacobian
  jacobian(1,1,1,&x,jac);
  cout<<"dense:"<<jac[0][0]<<endl;

  free(jac[0]);
  free(jac);


//sparse_jac
  int nnz;
  unsigned int *rindex=NULL;
  unsigned int *cindex=NULL;
  double *sparse_Jacobi_value=NULL;

  int options[4]={0,0,0,1};
  sparse_jac(1,1,1,0,&x,&nnz,&rindex,&cindex,&sparse_Jacobi_value,options);

   cout<<"number of sparse matrix elements: "<<nnz<<'\n'<<"sparse:
"<<sparse_Jacobi_value[0]<<endl;

  free(rindex);
  free(cindex);
  free(sparse_Jacobi_value);
  return 0;

}
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://list.coin-or.org/pipermail/adol-c/attachments/20110828/5fdbec71/attachment.html 


More information about the ADOL-C mailing list