I've used NR3 interpolation classes to construct the following error function:
// function calculates error function as a function of sigma and dielectric
// result is square of differences
Doub error_sq_func(VecDoub x)
{
int i, j;
double sigma, dielectric, result, sum = 0.0;

sigma = x[0]; dielectric = x[1];
for (i=0; i<func_parameters.n_splines; i++) {
if (eval_spline((unsigned long *) func_parameters.spline+i, sigma, dielectric, &result, func_parameters.error) < 0)
return NaN;
sum += sqr(func_parameters.vals[i] - result);
}
Of course, most of the code isn't shown here, but it has been tested and works well.

I just can't figure out how to call the Powell minimizer, all my attempts at building a constructor come up empty -- I've been trying to construct the second instance mentioned on page 513 of the book (with simple C++ function instead of functor). :confused:

Does anyone have an example of Powell minimization? Is there another place for NR3 examples?

...Dan

Doub error_sq_func(VecDoub x)
{
.
.
.
}

...attempts at building a constructor come up empty
With that function, the constructor could be like this:

Powell <Doub(VecDoub)> powell(error_sq_func); // Make sure there is no '&' in the template function argument!


Although, I think many (most?) C++ programmers would prefer the form shown on page 513 (Pass a const reference to an object rather than passing a copy of the object):


Doub error_sq_func(VecDoub_I & x)
{
.
.
.
}
In which case the constructor could be like this:

Powell <Doub(VecDoub_I &)> powell(error_sq_func);


The function implementation code would be exactly the same, whichever style you use.



Does anyone have an example of Powell minimization?

A program based on the one in the Numerical Recipes 2 C++ Examples:


#include "../code/nr3.h"
#include "../code/bessel.h"
#include "../code/mins.h"
#include "../code/mins_ndim.h"

using namespace std;

// Driver for Powell minimization
//
// Find minimum value of 0.5 - J0[(x-1)^2 +(y-2)^2 +(z-3)^2]
//

Doub func(VecDoub_I & x)
{
Bessjy bessjy;
return 0.5 - bessjy.j0(SQR(x[0]-1.0) + SQR(x[1]-2.0) + SQR(x[2]-3.0));
}

int main()
{

int NDIM = 3;
Doub initx[] = { // Initial guess: has NDIM elements
1.5, 1.5, 2.5
};
VecDoub p(NDIM, initx); // Initialize the vector with the guess
VecDoub pp;

// Feed the function to the constructor
Powell <Doub(VecDoub_I &)> powell(func);

// Perform the minimization
pp = powell.minimize(p);

cout << "Number of iterations = " << powell.iter << endl << endl;

cout << scientific;

cout << "Minimum function value = ";
cout << powell.fret
<< ", found at: " << endl;
for (int i = 0; i < NDIM; i++) {
cout << setw(15) << pp[i];
}
cout << endl << endl;


cout << "Actual minimum value = "
<< -0.5
<< " at:" << endl;
cout << setw(15) << 1.0
<< setw(15) << 2.0
<< setw(15) << 3.0
<< endl;
return 0;
}


Output

Number of iterations = 1

Minimum function value = -5.000000e-01, found at:
1.000000e+00 2.000000e+00 2.999954e+00

Actual minimum value = -5.000000e-01 at:
1.000000e+00 2.000000e+00 3.000000e+00


Regards,

Dave

davekw7x;

I modeled my code based on your recommendations and it compiles and appears to work (I haven't yet rigorously tested it and still have to obtain experimental data for verification). Anyway, a million thanks! :)

C++ is difficult enough, the added layers of templates, typedefs, and polymorphism of the library make this old FORTRAN programmer (30 years, 20 years for C) strain to make sense of it. :confused:

Again, thanks.

...Dan