htswanalysis/MACS/lib/gsl/gsl-1.11/multifit/lmset.c

   1 static int
   2 set (void *vstate, gsl_multifit_function_fdf * fdf, gsl_vector * x, gsl_vector * f, gsl_matrix * J, gsl_vector * dx, int scale)
   3 {
   4   lmder_state_t *state = (lmder_state_t *) vstate;
   5
   6   gsl_matrix *r = state->r;
   7   gsl_vector *tau = state->tau;
   8   gsl_vector *diag = state->diag;
   9   gsl_vector *work1 = state->work1;
  10   gsl_permutation *perm = state->perm;
  11
  12   int signum;
  13
  14   GSL_MULTIFIT_FN_EVAL_F_DF (fdf, x, f, J);
  15
  16   state->par = 0;
  17   state->iter = 1;
  18   state->fnorm = enorm (f);
  19
  20   gsl_vector_set_all (dx, 0.0);
  21
  22   /* store column norms in diag */
  23
  24   if (scale)
  25     {
  26       compute_diag (J, diag);
  27     }
  28   else
  29     {
  30       gsl_vector_set_all (diag, 1.0);
  31     }
  32
  33   /* set delta to 100 |D x| or to 100 if |D x| is zero */
  34
  35   state->xnorm = scaled_enorm (diag, x);
  36   state->delta = compute_delta (diag, x);
  37
  38   /* Factorize J into QR decomposition */
  39
  40   gsl_matrix_memcpy (r, J);
  41   gsl_linalg_QRPT_decomp (r, tau, perm, &signum, work1);
  42
  43   gsl_vector_set_zero (state->rptdx);
  44   gsl_vector_set_zero (state->w);
  45
  46   /* Zero the trial vector, as in the alloc function */
  47
  48   gsl_vector_set_zero (state->f_trial);
  49
  50 #ifdef DEBUG
  51   printf("r = "); gsl_matrix_fprintf(stdout, r, "%g");
  52   printf("perm = "); gsl_permutation_fprintf(stdout, perm, "%d");
  53   printf("tau = "); gsl_vector_fprintf(stdout, tau, "%g");
  54 #endif
  55
  56   return GSL_SUCCESS;
  57 }