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

   1 /* multifit/fdfsolver.c
   2  *
   3  * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Brian Gough
   4  *
   5  * This program is free software; you can redistribute it and/or modify
   6  * it under the terms of the GNU General Public License as published by
   7  * the Free Software Foundation; either version 3 of the License, or (at
   8  * your option) any later version.
   9  *
  10  * This program is distributed in the hope that it will be useful, but
  11  * WITHOUT ANY WARRANTY; without even the implied warranty of
  12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13  * General Public License for more details.
  14  *
  15  * You should have received a copy of the GNU General Public License
  16  * along with this program; if not, write to the Free Software
  17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  18  */
  19
  20 #include <config.h>
  21 #include <stdlib.h>
  22 #include <string.h>
  23 #include <gsl/gsl_errno.h>
  24 #include <gsl/gsl_multifit_nlin.h>
  25
  26 gsl_multifit_fdfsolver *
  27 gsl_multifit_fdfsolver_alloc (const gsl_multifit_fdfsolver_type * T,
  28                               size_t n, size_t p)
  29 {
  30   int status;
  31
  32   gsl_multifit_fdfsolver * s;
  33
  34   if (n < p)
  35     {
  36       GSL_ERROR_VAL ("insufficient data points, n < p", GSL_EINVAL, 0);
  37     }
  38
  39   s = (gsl_multifit_fdfsolver *) malloc (sizeof (gsl_multifit_fdfsolver));
  40
  41   if (s == 0)
  42     {
  43       GSL_ERROR_VAL ("failed to allocate space for multifit solver struct",
  44                      GSL_ENOMEM, 0);
  45     }
  46
  47   s->x = gsl_vector_calloc (p);
  48
  49   if (s->x == 0)
  50     {
  51       free (s);
  52       GSL_ERROR_VAL ("failed to allocate space for x", GSL_ENOMEM, 0);
  53     }
  54
  55   s->f = gsl_vector_calloc (n);
  56
  57   if (s->f == 0)
  58     {
  59       gsl_vector_free (s->x);
  60       free (s);
  61       GSL_ERROR_VAL ("failed to allocate space for f", GSL_ENOMEM, 0);
  62     }
  63
  64   s->J = gsl_matrix_calloc (n,p);
  65
  66   if (s->J == 0)
  67     {
  68       gsl_vector_free (s->x);
  69       gsl_vector_free (s->f);
  70       free (s);
  71       GSL_ERROR_VAL ("failed to allocate space for g", GSL_ENOMEM, 0);
  72     }
  73
  74   s->dx = gsl_vector_calloc (p);
  75
  76   if (s->dx == 0)
  77     {
  78       gsl_matrix_free (s->J);
  79       gsl_vector_free (s->x);
  80       gsl_vector_free (s->f);
  81       free (s);
  82       GSL_ERROR_VAL ("failed to allocate space for dx", GSL_ENOMEM, 0);
  83     }
  84
  85   s->state = malloc (T->size);
  86
  87   if (s->state == 0)
  88     {
  89       gsl_vector_free (s->dx);
  90       gsl_vector_free (s->x);
  91       gsl_vector_free (s->f);
  92       gsl_matrix_free (s->J);
  93       free (s);         /* exception in constructor, avoid memory leak */
  94
  95       GSL_ERROR_VAL ("failed to allocate space for multifit solver state",
  96                         GSL_ENOMEM, 0);
  97     }
  98
  99   s->type = T ;
 100
 101   status = (s->type->alloc)(s->state, n, p);
 102
 103   if (status != GSL_SUCCESS)
 104     {
 105       free (s->state);
 106       gsl_vector_free (s->dx);
 107       gsl_vector_free (s->x);
 108       gsl_vector_free (s->f);
 109       gsl_matrix_free (s->J);
 110       free (s);         /* exception in constructor, avoid memory leak */
 111
 112       GSL_ERROR_VAL ("failed to set solver", status, 0);
 113     }
 114
 115   s->fdf = NULL;
 116
 117   return s;
 118 }
 119
 120 int
 121 gsl_multifit_fdfsolver_set (gsl_multifit_fdfsolver * s,
 122                              gsl_multifit_function_fdf * f,
 123                              const gsl_vector * x)
 124 {
 125   if (s->f->size != f->n)
 126     {
 127       GSL_ERROR ("function size does not match solver", GSL_EBADLEN);
 128     }
 129
 130   if (s->x->size != x->size)
 131     {
 132       GSL_ERROR ("vector length does not match solver", GSL_EBADLEN);
 133     }
 134
 135   s->fdf = f;
 136   gsl_vector_memcpy(s->x,x);
 137
 138   return (s->type->set) (s->state, s->fdf, s->x, s->f, s->J, s->dx);
 139 }
 140
 141 int
 142 gsl_multifit_fdfsolver_iterate (gsl_multifit_fdfsolver * s)
 143 {
 144   return (s->type->iterate) (s->state, s->fdf, s->x, s->f, s->J, s->dx);
 145 }
 146
 147 void
 148 gsl_multifit_fdfsolver_free (gsl_multifit_fdfsolver * s)
 149 {
 150   (s->type->free) (s->state);
 151   free (s->state);
 152   gsl_vector_free (s->dx);
 153   gsl_vector_free (s->x);
 154   gsl_vector_free (s->f);
 155   gsl_matrix_free (s->J);
 156   free (s);
 157 }
 158
 159 const char *
 160 gsl_multifit_fdfsolver_name (const gsl_multifit_fdfsolver * s)
 161 {
 162   return s->type->name;
 163 }
 164
 165 gsl_vector *
 166 gsl_multifit_fdfsolver_position (const gsl_multifit_fdfsolver * s)
 167 {
 168   return s->x;
 169 }
 170