Problems Optimization

  1. Find, using one of the minimization routines from GSL, the minimum of the Rosenbrock function,
    f(x,y) = (1-x)² + 100(y-x²)²
    
  2. Least-squares fit.

    Suppose an experiment measuring the activity of a radioactive substance as function of time gives the following result,

    # t_i   y_i   e_i
     0.47  5.49  0.26
     1.41  4.08  0.12
     2.36  3.54  0.27
     3.30  2.61  0.10
     4.24  2.09  0.15
     5.18  1.91  0.11
     6.13  1.55  0.13
     7.07  1.47  0.07
     8.01  1.45  0.15
     8.95  1.25  0.09
    
    or
    double t[]= {0.47,1.41,2.36,3.30,4.24,5.18,6.13,7.07,8.01,8.95};
    double y[]= {5.49,4.08,3.54,2.61,2.09,1.91,1.55,1.47,1.45,1.25};
    double e[]= {0.26,0.12,0.27,0.10,0.15,0.11,0.13,0.07,0.15,0.09};
    int n = sizeof(t)/sizeof(t[0]);
    
    where the time ti, activity yi, and the experimental error σi are given in some rescaled units.

    Do the least-squares fit of the function

    f(t)=A*exp(-t/T)+B
    
    (where the initial amount A, the lifetime T, and the background noise B are the fitting parameters) to the data and estimate the lifetime, T, of the substance.

    Plot the experimental data and the fitting function.

    Hints:

    1. It is probably easier to use the algorithm without derivatives gsl_multimin_fminimizer_nmsimplex2.

    2. The (master) function to minimize is F(A,T,B)=∑i(f(ti)-yi)²/σi². That is, you have a three-dimensional minimization problem in the space of the parameters A, T, and B.

      It can be programmed like this,

      struct experimental_data {int n; double *t,*y,*e;};
      
      double function_to_minimize (const gsl_vector *x, void *params) {
      	double  A = gsl_vector_get(x,0);
      	double  T = gsl_vector_get(x,1);
      	double  B = gsl_vector_get(x,2);
      	struct experimental_data *p = (struct experimental_data*) params;
      	int     n = p->n;
      	double *t = p->t;
      	double *y = p->y;
      	double *e = p->e;
      	double sum=0;
      	#define f(t) A*exp(-(t)/T) + B
      	for(int i=0;i<n;i++) sum+ = pow( f(t[i]) - y[i] )/e[i] ,2);
      	return sum;
      }
      
  3. (Optional) Do the same least-squares fit using one of the algorithms with derivatives. It might be more effective, that is, it might arrive at the optimum with fewer evaluations of the function-to-minimize.