AnimaL

Tutorial

Documentation

---

Documentation

• Overview
• Modules
• Namespace List
• Class Hierarchy
• Alphabetical List
• Compound List
• File List
• Namespace Members
• Compound Members
• File Members
• Examples

---

linmin.h

Go to the documentation of this file.

#ifndef NR_LINMIN_H
#define NR_LINMIN_H

#include <iostream>
#include "brent.h"
#include "mnbrak.h"

namespace nr {


template < class TraitsT = Old_Traits<> >
struct linmin {

    typedef TraitsT Traits;
    typedef typename Traits::Real Real;
    typedef typename Traits::Vector Vector;

    
    template < class Func >
    linmin (
        Vector& p, 
        Vector& xi, 
        unsigned int n, 
        Real *fret, 
        Func (*func),
        const Real TOL = 2.0e-4
    )
    {
        unsigned int j;
        Real xx,xmin,fx,fb,fa,bx,ax;

        F1dim<Real,Vector,Func> f1dim;
        Traits::alloc_vector(f1dim.pcom,n);
        Traits::alloc_vector(f1dim.xicom,n);
        Traits::alloc_vector(f1dim.xt,n);
        f1dim.funct = func;
        f1dim.ncom = n;

        for (j=0;j<n;j++) {
            f1dim.pcom[j]=p[j];
            f1dim.xicom[j]=xi[j];
        }
        ax=0.0;
        xx=1.0;
//      std::cout<<"linmin, p init = "<<p<<std::endl;
//      std::cout<<"linmin, xi init = "<<xi<<std::endl;
        mnbrak(&ax,&xx,&bx,&fa,&fx,&fb,&f1dim);
        //std::cout<<"linmin: mnbrak returns "<<ax<<" "<<xx<<" "<<bx<<", values: "<<fa<<" "<<fx<<" "<<fb<<std::endl;
        *fret=brent(ax,xx,bx,&f1dim,TOL,&xmin);
//      std::cout<<"linmin, brent returns "<<*fret<<" at "<<xmin<<std::endl;
//      std::cout<<"linmin, p = "<<p<<std::endl;
//      std::cout<<"linmin, xi = "<<xi<<std::endl;
        for (j=0;j<n;j++) {
            xi[j] *= xmin;
            p[j] += xi[j];
        }
//      std::cout<<"linmin, new p = "<<p<<std::endl;
//      std::cout<<"linmin, new xi = "<<xi<<std::endl;
        Traits::free_vector( f1dim.pcom );
        Traits::free_vector( f1dim.xicom );
        Traits::free_vector( f1dim.xt );
//      std::cout<<"end linmin"<<std::endl<<std::endl<<std::endl<<std::endl<<std::endl;
    }

    


    template < class Real, class Vector, class Func >
    struct F1dim 
    {
        Real operator () (Real x)
        {
            //std::cout << "linmin::F1dim: " << x << std::endl;
            for ( unsigned int j=0;j<ncom;j++) {
                xt[j]=pcom[j]+x*xicom[j];
            }
            //std::cout<<" linmin: xt = "<<xt<<std::endl;
            return (*funct)( xt );
        }
        
        F1dim(){}

        Vector xt;    
        Vector pcom;    
        Vector xicom;    
        unsigned int ncom;  
        Func (*funct);  
    
    };

};

} // end namespace nr



#endif

---