[ADOL-C] manual hessian parallelization

[Waldemar Hujo]

Dear all.

I want to compute the hessian matrix in parallel, where each step should give
just one component of the hessian matrix. That is, if my hessian matrix is a
6x6 matrix I want to get only one number from the entire matrix.
The background of this is, that I want to compute the hessian matrix on
different computers and assemble the matrix afterwards.

Up to know I can compute an entire row with the following code, which is just
taken from the original adol-c source code. In this code a loop is used to get
every row of the matrix. I could split up the loop and compute the rows
separately.

However, as I said, I want to compute just one element of the matrix. I stuck
with the "fos_forward" and the "hos_reverse" functions, where I suppose the
rows to be computed.

I cannot figure out how to split up this step. Has anybody a hint for me?

Thanks in advance.
Waldemar



///////////////////////////////////////////////////////
// Parallel Hessian evaluation
///////////////////////////////////////////////////////
  cout << "Parallel Hessian " << endl;

  double *p_argument = xp;
  double** p_hess   = new double*[n];
  for(i=0; i<n; i++)
      p_hess[i] = new double[n];
  for(int i=0; i<n; ++i){
   for(int j=0; j<i+1; ++j){
     p_hess[i][j] = 0.0 ;
   }
  }
  short tag = 1;

  double lagrange = 1.0;

    int rc= 3;
    double *tangent = myalloc1(n);
    double *result = myalloc1(n);
    for(i=0;i<n;i++) tangent[i] = 0;
    for(i=0;i<n;i++) {
        tangent[i] = 1;
/////
    int m=1;
    rc=-1;
    int degree = 1;
    int keep = degree+1;
    double **X, *y, *y_tangent;

    X = myalloc2(n,2);
    y = myalloc1(m);
    y_tangent = myalloc1(m);

    rc = fos_forward(tag, m, n, keep, p_argument, tangent, y, y_tangent);

    if(rc < 0){
      cout << "parallel hessian error 1" << endl;
      exit(1);
    }

    MINDEC(rc, hos_reverse(tag, m, n, degree, &lagrange, X));

    for(int k = 0; k < n; ++k){
        result[k] = X[k][1];
         // print the row
        cout << "r " << result[k] << endl;
    }
    cout << endl;

    myfree1(y_tangent);
    myfree1(y);
    myfree2(X);
/////
        if( rc < 0) {
            free((char *)tangent);
            free((char *) result);
            cout << "parallel hessian error 2" << endl;
            exit(1);
        }
        for(j=0;j<=i;j++){
            p_hess[i][j] = result[j];
        }
        tangent[i] = 0;
    }


    free((char *)tangent);
    free((char *) result);

/////////////////////////////////////////////////////////