File indexing completed on 2025-02-21 10:03:51
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #ifndef __GSL_MULTIFIT_NLIN_H__
0021 #define __GSL_MULTIFIT_NLIN_H__
0022
0023 #include <stdlib.h>
0024 #include <gsl/gsl_types.h>
0025 #include <gsl/gsl_math.h>
0026 #include <gsl/gsl_vector.h>
0027 #include <gsl/gsl_matrix.h>
0028 #include <gsl/gsl_permutation.h>
0029
0030 #undef __BEGIN_DECLS
0031 #undef __END_DECLS
0032 #ifdef __cplusplus
0033 # define __BEGIN_DECLS extern "C" {
0034 # define __END_DECLS }
0035 #else
0036 # define __BEGIN_DECLS
0037 # define __END_DECLS
0038 #endif
0039
0040 __BEGIN_DECLS
0041
0042 int gsl_multifit_gradient (const gsl_matrix * J, const gsl_vector * f,
0043 gsl_vector * g);
0044
0045 int gsl_multifit_covar (const gsl_matrix * J, const double epsrel, gsl_matrix * covar);
0046 int gsl_multifit_covar_QRPT (gsl_matrix * r, gsl_permutation * perm,
0047 const double epsrel, gsl_matrix * covar);
0048
0049
0050
0051
0052 struct gsl_multifit_function_struct
0053 {
0054 int (* f) (const gsl_vector * x, void * params, gsl_vector * f);
0055 size_t n;
0056 size_t p;
0057 void * params;
0058 };
0059
0060 typedef struct gsl_multifit_function_struct gsl_multifit_function ;
0061
0062 #define GSL_MULTIFIT_FN_EVAL(F,x,y) (*((F)->f))(x,(F)->params,(y))
0063
0064 typedef struct
0065 {
0066 const char *name;
0067 size_t size;
0068 int (*alloc) (void *state, size_t n, size_t p);
0069 int (*set) (void *state, gsl_multifit_function * function, gsl_vector * x, gsl_vector * f, gsl_vector * dx);
0070 int (*iterate) (void *state, gsl_multifit_function * function, gsl_vector * x, gsl_vector * f, gsl_vector * dx);
0071 void (*free) (void *state);
0072 }
0073 gsl_multifit_fsolver_type;
0074
0075 typedef struct
0076 {
0077 const gsl_multifit_fsolver_type * type;
0078 gsl_multifit_function * function ;
0079 gsl_vector * x ;
0080 gsl_vector * f ;
0081 gsl_vector * dx ;
0082 void *state;
0083 }
0084 gsl_multifit_fsolver;
0085
0086 gsl_multifit_fsolver *
0087 gsl_multifit_fsolver_alloc (const gsl_multifit_fsolver_type * T,
0088 size_t n, size_t p);
0089
0090 void gsl_multifit_fsolver_free (gsl_multifit_fsolver * s);
0091
0092 int gsl_multifit_fsolver_set (gsl_multifit_fsolver * s,
0093 gsl_multifit_function * f,
0094 const gsl_vector * x);
0095
0096 int gsl_multifit_fsolver_iterate (gsl_multifit_fsolver * s);
0097
0098 int gsl_multifit_fsolver_driver (gsl_multifit_fsolver * s,
0099 const size_t maxiter,
0100 const double epsabs, const double epsrel);
0101
0102 const char * gsl_multifit_fsolver_name (const gsl_multifit_fsolver * s);
0103 gsl_vector * gsl_multifit_fsolver_position (const gsl_multifit_fsolver * s);
0104
0105
0106
0107
0108 struct gsl_multifit_function_fdf_struct
0109 {
0110 int (* f) (const gsl_vector * x, void * params, gsl_vector * f);
0111 int (* df) (const gsl_vector * x, void * params, gsl_matrix * df);
0112 int (* fdf) (const gsl_vector * x, void * params, gsl_vector * f, gsl_matrix *df);
0113 size_t n;
0114 size_t p;
0115 void * params;
0116 size_t nevalf;
0117 size_t nevaldf;
0118 };
0119
0120 typedef struct gsl_multifit_function_fdf_struct gsl_multifit_function_fdf ;
0121
0122 typedef struct
0123 {
0124 const char *name;
0125 size_t size;
0126 int (*alloc) (void *state, size_t n, size_t p);
0127 int (*set) (void *state, const gsl_vector * wts,
0128 gsl_multifit_function_fdf * fdf, gsl_vector * x,
0129 gsl_vector * f, gsl_vector * dx);
0130 int (*iterate) (void *state, const gsl_vector * wts,
0131 gsl_multifit_function_fdf * fdf, gsl_vector * x,
0132 gsl_vector * f, gsl_vector * dx);
0133 int (*gradient) (void *state, gsl_vector * g);
0134 int (*jac) (void *state, gsl_matrix * J);
0135 void (*free) (void *state);
0136 }
0137 gsl_multifit_fdfsolver_type;
0138
0139 typedef struct
0140 {
0141 const gsl_multifit_fdfsolver_type * type;
0142 gsl_multifit_function_fdf * fdf ;
0143 gsl_vector * x;
0144 gsl_vector * f;
0145 gsl_vector * dx;
0146 gsl_vector * g;
0147 gsl_vector * sqrt_wts;
0148 size_t niter;
0149 void *state;
0150 }
0151 gsl_multifit_fdfsolver;
0152
0153
0154 gsl_multifit_fdfsolver *
0155 gsl_multifit_fdfsolver_alloc (const gsl_multifit_fdfsolver_type * T,
0156 size_t n, size_t p);
0157
0158 int
0159 gsl_multifit_fdfsolver_set (gsl_multifit_fdfsolver * s,
0160 gsl_multifit_function_fdf * fdf,
0161 const gsl_vector * x);
0162 int gsl_multifit_fdfsolver_wset (gsl_multifit_fdfsolver * s,
0163 gsl_multifit_function_fdf * f,
0164 const gsl_vector * x,
0165 const gsl_vector * wts);
0166
0167 int
0168 gsl_multifit_fdfsolver_iterate (gsl_multifit_fdfsolver * s);
0169
0170 int gsl_multifit_fdfsolver_driver (gsl_multifit_fdfsolver * s,
0171 const size_t maxiter,
0172 const double xtol,
0173 const double gtol,
0174 const double ftol,
0175 int *info);
0176
0177 int gsl_multifit_fdfsolver_jac (gsl_multifit_fdfsolver * s,
0178 gsl_matrix * J);
0179
0180 void
0181 gsl_multifit_fdfsolver_free (gsl_multifit_fdfsolver * s);
0182
0183 const char * gsl_multifit_fdfsolver_name (const gsl_multifit_fdfsolver * s);
0184 gsl_vector * gsl_multifit_fdfsolver_position (const gsl_multifit_fdfsolver * s);
0185 gsl_vector * gsl_multifit_fdfsolver_residual (const gsl_multifit_fdfsolver * s);
0186 size_t gsl_multifit_fdfsolver_niter (const gsl_multifit_fdfsolver * s);
0187 int gsl_multifit_eval_wf(gsl_multifit_function_fdf *fdf,
0188 const gsl_vector *x, const gsl_vector *wts,
0189 gsl_vector *y);
0190 int gsl_multifit_eval_wdf(gsl_multifit_function_fdf *fdf,
0191 const gsl_vector *x, const gsl_vector *wts,
0192 gsl_matrix *dy);
0193
0194 int gsl_multifit_fdfsolver_test (const gsl_multifit_fdfsolver * s,
0195 const double xtol,
0196 const double gtol,
0197 const double ftol, int *info);
0198 int gsl_multifit_test_delta (const gsl_vector * dx, const gsl_vector * x,
0199 double epsabs, double epsrel);
0200
0201 int gsl_multifit_test_gradient (const gsl_vector * g, double epsabs);
0202
0203 int gsl_multifit_fdfsolver_dif_df(const gsl_vector *x,
0204 const gsl_vector *wts,
0205 gsl_multifit_function_fdf *fdf,
0206 const gsl_vector *f, gsl_matrix *J);
0207 int gsl_multifit_fdfsolver_dif_fdf(const gsl_vector *x, gsl_multifit_function_fdf *fdf,
0208 gsl_vector *f, gsl_matrix *J);
0209
0210 typedef struct
0211 {
0212 size_t n;
0213 size_t p;
0214 double lambda;
0215 const gsl_vector *L_diag;
0216 const gsl_matrix *L;
0217 gsl_vector *f;
0218 gsl_vector *wts;
0219 gsl_multifit_fdfsolver * s;
0220 gsl_multifit_function_fdf *fdf;
0221 gsl_multifit_function_fdf fdftik;
0222 } gsl_multifit_fdfridge;
0223
0224 gsl_multifit_fdfridge *
0225 gsl_multifit_fdfridge_alloc (const gsl_multifit_fdfsolver_type * T,
0226 const size_t n, const size_t p);
0227 void gsl_multifit_fdfridge_free(gsl_multifit_fdfridge *work);
0228 const char *gsl_multifit_fdfridge_name(const gsl_multifit_fdfridge * w);
0229 gsl_vector *gsl_multifit_fdfridge_position (const gsl_multifit_fdfridge * w);
0230 gsl_vector *gsl_multifit_fdfridge_residual (const gsl_multifit_fdfridge * w);
0231 size_t gsl_multifit_fdfridge_niter (const gsl_multifit_fdfridge * w);
0232 int gsl_multifit_fdfridge_set (gsl_multifit_fdfridge * w,
0233 gsl_multifit_function_fdf * f,
0234 const gsl_vector * x,
0235 const double lambda);
0236 int gsl_multifit_fdfridge_wset (gsl_multifit_fdfridge * w,
0237 gsl_multifit_function_fdf * f,
0238 const gsl_vector * x,
0239 const double lambda,
0240 const gsl_vector * wts);
0241 int gsl_multifit_fdfridge_set2 (gsl_multifit_fdfridge * w,
0242 gsl_multifit_function_fdf * f,
0243 const gsl_vector * x,
0244 const gsl_vector * lambda);
0245 int gsl_multifit_fdfridge_wset2 (gsl_multifit_fdfridge * w,
0246 gsl_multifit_function_fdf * f,
0247 const gsl_vector * x,
0248 const gsl_vector * lambda,
0249 const gsl_vector * wts);
0250 int gsl_multifit_fdfridge_set3 (gsl_multifit_fdfridge * w,
0251 gsl_multifit_function_fdf * f,
0252 const gsl_vector * x,
0253 const gsl_matrix * L);
0254 int gsl_multifit_fdfridge_wset3 (gsl_multifit_fdfridge * w,
0255 gsl_multifit_function_fdf * f,
0256 const gsl_vector * x,
0257 const gsl_matrix * L,
0258 const gsl_vector * wts);
0259 int gsl_multifit_fdfridge_iterate (gsl_multifit_fdfridge * w);
0260 int gsl_multifit_fdfridge_driver (gsl_multifit_fdfridge * w,
0261 const size_t maxiter,
0262 const double xtol,
0263 const double gtol,
0264 const double ftol,
0265 int *info);
0266
0267
0268
0269 GSL_VAR const gsl_multifit_fdfsolver_type * gsl_multifit_fdfsolver_lmsder;
0270 GSL_VAR const gsl_multifit_fdfsolver_type * gsl_multifit_fdfsolver_lmder;
0271 GSL_VAR const gsl_multifit_fdfsolver_type * gsl_multifit_fdfsolver_lmniel;
0272
0273 __END_DECLS
0274
0275 #endif