double compute_ba_J()

in tools/ADOLC/main.cpp [406:506]


double compute_ba_J(int nruns, double time_limit, int n, int m, int p,
  double *cams, double *X, double *w, int *obs, double *feats,
  double *reproj_err, double *w_err, BASparseMat *J, double* t_sparsity)
{
  if (nruns == 0)
    return 0.;

  bool doRowCompression = false;
  high_resolution_clock::time_point start, end;

  int tapeTag = 1;

  start = high_resolution_clock::now();

  adouble *acams, *aX, *aw, *areproj_err,
    *af_prior_err, *aw_err;

  areproj_err = new adouble[2 * p];
  af_prior_err = new adouble[n - 2];
  aw_err = new adouble[p];

  // Record on a tape
  trace_on(tapeTag);

  acams = new adouble[BA_NCAMPARAMS*n];
  for (int i = 0; i < BA_NCAMPARAMS*n; i++)
    acams[i] <<= cams[i];
  aX = new adouble[3 * m];
  for (int i = 0; i < 3 * m; i++)
    aX[i] <<= X[i];
  aw = new adouble[p];
  for (int i = 0; i < p; i++)
    aw[i] <<= w[i];

  ba_objective(n, m, p, acams, aX, aw, obs, feats,
    areproj_err, aw_err);

  for (int i = 0; i < 2 * p; i++)
    areproj_err[i] >>= reproj_err[i];

  for (int i = 0; i < p; i++)
    aw_err[i] >>= w_err[i];

  trace_off();

  delete[] acams;
  delete[] aX;
  delete[] aw;
  delete[] areproj_err;
  delete[] af_prior_err;
  delete[] aw_err;

  end = high_resolution_clock::now();
  double t_tape = duration_cast<duration<double>>(end - start).count();

  //////// Compute J and compute sparsity always again and again
  vector<double> in(J->ncols);
  memcpy(&in[0], cams, BA_NCAMPARAMS*n*sizeof(double));
  int off = BA_NCAMPARAMS*n;
  memcpy(&in[off], X, 3 * m*sizeof(double));
  off += 3 * m;
  memcpy(&in[off], w, p*sizeof(double));

  int opt[4];
  opt[0] = 0; // default
  opt[1] = 0; // default
  opt[2] = 0; // 0=auto 1=F 2=R
  opt[3] = doRowCompression ? 1 : 0;
  int nnz;
  unsigned int *ridxs = nullptr, *cidxs = nullptr;
  double *nzvals = nullptr;

  int samePattern = 0;
  start = high_resolution_clock::now();
  sparse_jac(tapeTag, J->nrows, J->ncols, samePattern,
    in.data(), &nnz, &ridxs, &cidxs, &nzvals, opt);
  end = high_resolution_clock::now();
  *t_sparsity = duration_cast<duration<double>>(end - start).count();

  samePattern = 1;
  start = high_resolution_clock::now();
  for (int i = 0; i < nruns; i++)
  {
    sparse_jac(tapeTag, J->nrows, J->ncols, samePattern,
      in.data(), &nnz, &ridxs, &cidxs, &nzvals, opt);
  }
  end = high_resolution_clock::now();
  double t_J = duration_cast<duration<double>>(end - start).count() / nruns;

  convert_J(nnz, ridxs, cidxs, nzvals, J);

  delete[] ridxs;
  delete[] cidxs;
  delete[] nzvals;

  cout << "t_tape: " << t_tape << endl;
  cout << "t_sparsity: " << *t_sparsity << endl;
  cout << "t_J:" << t_J << endl;

  return t_J;
}