Gromacs
5.1.5
|
#include "gmxpre.h"
#include "expfit.h"
#include <math.h>
#include <string.h>
#include <algorithm>
#include "external/lmfit/lmcurve.h"
#include "gromacs/correlationfunctions/integrate.h"
#include "gromacs/fileio/xvgr.h"
#include "gromacs/legacyheaders/macros.h"
#include "gromacs/math/vec.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/futil.h"
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/real.h"
#include "gromacs/utility/smalloc.h"
Implements routine for fitting a data set to a curve.
Typedefs | |
typedef double(* | t_lmcurve )(double x, const double *a) |
function type for passing to fitting routine | |
Functions | |
int | effnNparams (int effn) |
Returns number of function parameters associated with a fitting function. More... | |
const char * | effnDescription (int effn) |
Returns description corresponding to the enum above, or NULL if out of range. More... | |
int | sffn2effn (const char **sffn) |
Returns corresponding to the selected enum option in sffn. More... | |
static double | myexp (double x, double A, double tau) |
Compute exponential function A exp(-x/tau) | |
static double | lmc_erffit (double x, const double *a) |
Compute y=(a0+a1)/2-(a0-a1)/2*erf((x-a2)/a3^2) | |
static double | safe_exp (double x) |
Exponent function that prevents overflow. | |
static double | safe_expm1 (double x) |
Exponent minus 1 function that prevents overflow. | |
static double | lmc_exp_one_parm (double x, const double *a) |
Compute y = exp(-x/|a0|) | |
static double | lmc_exp_two_parm (double x, const double *a) |
Compute y = a1 exp(-x/|a0|) | |
static double | lmc_exp_exp (double x, const double *a) |
Compute y = a1 exp(-x/|a0|) + (1-a1) exp(-x/|a2|) | |
static double | lmc_exp_5_parm (double x, const double *a) |
Compute y = a0 exp(-x/|a1|) + a2 exp(-x/(|a1|+|a3|)) + a4. | |
static double | lmc_exp_7_parm (double x, const double *a) |
Compute 7 parameter exponential function value. More... | |
static double | lmc_exp_9_parm (double x, const double *a) |
Compute 9 parameter exponential function value. More... | |
static double | lmc_pres_6_parm (double x, const double *a) |
Compute y = (1-a0)*exp(-(x/|a2|)^|a3|)*cos(x*|a1|) + a0*exp(-(x/|a4|)^|a5|) | |
static double | lmc_vac_2_parm (double x, const double *a) |
Compute vac function. | |
static double | lmc_errest_3_parm (double x, const double *a) |
Compute error estimate. | |
double | fit_function (const int eFitFn, const double parm[], const double x) |
Returns the value of fit function eFitFn at x. More... | |
static gmx_bool | lmfit_exp (int nfit, const double x[], const double y[], const double dy[], double parm[], gmx_bool bVerbose, int eFitFn, int nfix) |
lmfit_exp supports fitting of different functions More... | |
static void | initiate_fit_params (int eFitFn, double params[]) |
Ensure the fitting parameters are well-behaved. More... | |
static void | extract_fit_params (int eFitFn, double params[]) |
Process the fitting parameters to get output parameters. More... | |
static void | print_chi2_params (FILE *fp, const int eFitFn, const double fitparms[], const char *label, const int nfitpnts, const double x[], const double y[]) |
Print chi-squared value and the parameters. | |
real | do_lmfit (int ndata, real c1[], real sig[], real dt, real x0[], real begintimefit, real endtimefit, const output_env_t oenv, gmx_bool bVerbose, int eFitFn, double fitparms[], int fix, const char *fn_fitted) |
See description in header file. | |
real | fit_acf (int ncorr, int fitfn, const output_env_t oenv, gmx_bool bVerbose, real tbeginfit, real tendfit, real dt, real c1[], real *fit) |
Fit an autocorrelation function to a pre-defined functional form. More... | |
Variables | |
static const int | nfp_ffn [effnNR] = { 0, 1, 2, 3, 5, 7, 9, 2, 4, 3, 6 } |
Number of parameters for each fitting function. | |
const char * | s_ffn [effnNR+2] |
Short name of each function type. This is exported for now in order to use when calling parse_common_args. More... | |
static const char * | longs_ffn [effnNR] |
Long description for each fitting function type. More... | |
t_lmcurve | lmcurves [effnNR+1] |
array of fitting functions corresponding to the pre-defined types More... | |
const char* effnDescription | ( | int | effn | ) |
Returns description corresponding to the enum above, or NULL if out of range.
[in] | effn | Index |
int effnNparams | ( | int | effn | ) |
Returns number of function parameters associated with a fitting function.
[in] | effn | Index |
|
static |
Process the fitting parameters to get output parameters.
See comment at the previous function.
real fit_acf | ( | int | ncorr, |
int | fitfn, | ||
const output_env_t | oenv, | ||
gmx_bool | bVerbose, | ||
real | tbeginfit, | ||
real | tendfit, | ||
real | dt, | ||
real | c1[], | ||
real * | fit | ||
) |
Fit an autocorrelation function to a pre-defined functional form.
See description in header file.
double fit_function | ( | const int | eFitFn, |
const double | parm[], | ||
const double | x | ||
) |
Returns the value of fit function eFitFn at x.
[in] | eFitFn | the index to the fitting function (0 .. effnNR) |
[in] | parm | Array of parameters, the length of which depends on eFitFn |
[in] | x | The value of x |
|
static |
Ensure the fitting parameters are well-behaved.
In order to make sure that fitting behaves according to the constraint that time constants are positive and increasing we enforce this here by setting all time constants to their absolute value and by adding e.g. |a_0| to |a_2|. This is done in conjunction with the fitting functions themselves. When there are multiple time constants we make sure that the successive time constants are at least double the previous ones and during fitting we enforce the they remain larger. It may very well help the convergence of the fitting routine.
|
static |
Compute 7 parameter exponential function value.
Compute y = a0 exp(-x/|a1|) + a2 exp(-x/(|a1|+|a3|)) + a4 exp(-x/(|a1|+|a3|+|a5|)) + a6
|
static |
Compute 9 parameter exponential function value.
Compute y = a0 exp(-x/|a1|) + a2 exp(-x/(|a1|+|a3|)) + a4 exp(-x/(|a1|+|a3|+|a5|)) + a6 exp(-x/(|a1|+|a3|+|a5|+|a7|)) + a8
|
static |
lmfit_exp supports fitting of different functions
This routine calls the Levenberg-Marquardt non-linear fitting routine for fitting a data set with errors to a target function. Fitting routines included in gromacs in src/external/lmfit.
int sffn2effn | ( | const char ** | sffn | ) |
Returns corresponding to the selected enum option in sffn.
[in] | sffn | Two dimensional string array coming from parse_common_args |
t_lmcurve lmcurves[effnNR+1] |
array of fitting functions corresponding to the pre-defined types
|
static |
Long description for each fitting function type.
const char* s_ffn[effnNR+2] |
Short name of each function type. This is exported for now in order to use when calling parse_common_args.