git-svn-id: file:///home/svn/opt_utilities@118 ed2142bd-67ad-457f-ba7c-d818d4011675

5 files changed, 1 insertions, 364 deletions

diff --git a/methods/powell/bas_util.hpp b/methods/powell/bas_util.hpp
deleted file mode 100644
index ccdce9e..0000000
--- a/methods/powell/bas_util.hpp
+++ /dev/null
@@ -1,65 +0,0 @@
-#ifndef BAS_UTIL
-#define BAS_UTIL
-#define OPT_HEADER
-#include <core/opt_traits.hpp>
-#include <algorithm>
-namespace opt_utilities
-{
-  template <typename T>
-  T tabs(T x)
-  {
-    return T(x)<T(0)?T(-x):T(x);
-  }
-
-  template <typename T>
-  T sqr(T x)
-  {
-    return x*x;
-  }
-  
-  
-  template <typename T>
-  void shft3(T&a,T& b,T& c,T d)
-  {
-    opt_eq(a,b);
-    opt_eq(b,c);
-    opt_eq(c,d);
-  }
-
-  template <typename T>
-  void shft(T& a,T& b,T& c,T d)
-  {
-    opt_eq(a,b);
-    opt_eq(b,c);
-    opt_eq(c,d);
-  }
-  //  template <typename T>  
-  //  void swap(T& ax,T& bx)
-  //{
-    //  swap(ax,bx);
-    //    T temp;
-    //opt_eq(temp,ax);
-    //opt_eq(ax,bx);
-    //opt_eq(bx,temp);
-    //}
-  
-  template <typename T>
-  T sign(const T& a,const T& b)
-  {
-    return b>=0?T(a>=0?T(a):T(-a)):T(a>=0?T(-a):T(a));
-  }
-  
-  template <typename T>
-  T max(T a,T b)
-  {
-    return b>a?T(b):T(a);
-  }
-
-  template <typename T>
-  void mov3(T& a,T& b,T& c, T& d,T& e,T& f)
-  {
-    opt_eq(a,d);opt_eq(b,e);opt_eq(c,f);
-  }
-}
-
-#endif
diff --git a/methods/powell/brent.hpp b/methods/powell/brent.hpp
deleted file mode 100644
index 1a0ee39..0000000
--- a/methods/powell/brent.hpp
+++ /dev/null
@@ -1,111 +0,0 @@
-#ifndef BRENT_HPP
-#define BRENT_HPP
-#define OPT_HEADER
-#include <iostream>
-#include "bas_util.hpp"
-//#include "optimizer.hpp"
-namespace opt_utilities
-{
-  template<typename T>
-  T brent(T ax,T bx,T cx,func_obj<T,T>& f,T tol,T& xmin)
-  {
-    const int ITMAX=100;
-    const T CGOLD=0.3819660;
-    const T ZEPS=std::numeric_limits<T>::epsilon()*1.e-3;
-    
-    int iter;
-    T a=0,b=0,d(0),etemp=0,fu=0,fv=0,fw=0,fx=0,p=0,q=0
-      ,r=0,tol1=0,tol2=0,u=0,v=0,w=0,x=0,xm=0;
-    T e=0.;
-    a=(ax<cx?ax:cx);
-    b=(ax>cx?ax:cx);
-    x=w=v=bx;
-    fw=fv=fx=f.eval(x);
-    for(iter=0;iter<ITMAX;++iter)
-      {
-	xm=.5*(a+b);
-	tol2=2.*(tol1=tol*tabs(x)+ZEPS);
-	if(tabs(T(x-xm))<=(tol2-.5*(b-a)))
-	  {
-	    xmin=x;
-	    return fx;
-	  }
-	if(tabs(e)>tol1)
-	  {
-	    r=(x-w)*(fx-fv);
-	    q=(x-v)*(fx-fw);
-	    p=(x-v)*q-(x-w)*r;
-	    q=2.*(q-r);
-	    if(q>0.)
-	      {
-		p=-p;
-	      }
-	    q=tabs(q);
-	    etemp=e;
-	    e=d;
-	    if(tabs(p)>=tabs(T(T(.5)*p*etemp))||p<=q*(a-x)||p>=q*(b-x))
-	      {
-		d=CGOLD*(e=(x>=xm?a-x:b-x));
-	      }
-	    else
-	      {
-		d=p/q;
-		u=x+d;
-		if(u-a<tol2||b-u<tol2)
-		  {
-		    d=sign(tol1,T(xm-x));
-		  }
-	      }
-	    
-	  }
-	else
-	  {
-	    d=CGOLD*(e=(x>=xm?a-x:b-x));
-	  }
-	u=(tabs(d)>=tol1?x+d:x+sign(tol1,d));
-	fu=f.eval(u);
-	if(fu<=fx)
-	  {
-	    if(u>=x)
-	      {
-		a=x;
-	      }
-	    else
-	      {
-		b=x;
-	      }
-	    shft3(v,w,x,u);
-	    shft3(fv,fw,fx,fu);
-	  }
-	else
-	  {
-	    if(u<x)
-	      {
-		a=u;
-	      }
-	    else
-	      {
-		b=u;
-	      }
-	    if(fu<=fw||w==x)
-	      {
-		v=w;
-		w=u;
-		fv=fw;
-		fw=fu;
-	      }
-	    else if(fu<=fv||v==x||v==w)
-	      {
-		v=u;
-		fv=fu;
-	      }
-	  }
-      }
-    std::cerr<<"Too many iterations in brent"<<std::endl;
-    xmin=x;
-    return fx;
-    
-  }
-}
-
-#endif
diff --git a/methods/powell/linmin.hpp b/methods/powell/linmin.hpp
deleted file mode 100644
index 962b052..0000000
--- a/methods/powell/linmin.hpp
+++ /dev/null
@@ -1,105 +0,0 @@
-#ifndef LINMIN_HPP
-#define LINMIN_HPP
-#define OPT_HEADER
-#include "mnbrak.hpp"
-#include "brent.hpp"
-#include <core/opt_traits.hpp>
-
-namespace opt_utilities
-{
-  template <typename rT,typename pT>
-  class func_adaptor
-    :public func_obj<rT,rT>
-  {
-  private:
-    const pT p1,xi1;
-    const func_obj<rT,pT>* pfoo;
-    func_adaptor(){}
-    func_adaptor(const func_adaptor&)
-      :func_obj<rT,rT>(),p1(),xi1(),pfoo(0)
-    {}
-    
-  public:
-    /*
-    void set_origin(pT& p2)
-    {
-      p1=p2;
-    }
-    
-    void set_direction(pT& xi2)
-    {
-      xi1=xi2;
-    }
-    
-    void set_func_obj(func_obj<rT,pT>& foo)
-    {
-      pfoo=&foo;
-      }*/
-  public:
-    func_adaptor(const pT& _p,const pT& _xi,const func_obj<rT,pT>& pf)
-      :p1(_p),xi1(_xi),pfoo(&pf)
-    {}
-    
-  private:
-    func_obj<rT,rT>* do_clone()const
-    {
-      return new func_adaptor(*this);
-    }
-      
-    rT do_eval(const rT& x)
-    {
-      //assert(p1.size()==xi1.size());
-      
-      pT xt;
-      opt_eq(xt,p1);
-      for(size_t i=0;i<get_size(xt);++i)
-	{
-	  //get_element(xt,i)+=x*get_element((pT)xi1,i);
-	  set_element(xt,i,
-		      get_element(xt,i)+x*get_element((pT)xi1,i));
-	  //get_element((pT)xi1,i);
-	}
-      return const_cast<func_obj<rT,pT>&>(*pfoo).eval(xt);
-      //return x;
-    }
-  };
-  
-
-  template<typename rT,typename pT>
-  void linmin(pT& p,pT& xi,rT& fret,func_obj<rT,pT>& func)
-  {
-
-    //  assert(p.size()==10);
-    //assert(xi.size()==10);
-    func_adaptor<rT,pT> fadpt(p,xi,func);
-
-    int j=0;
-    const rT TOL=sqrt(std::numeric_limits<rT>::epsilon());
-    rT xx=0,xmin=0,fx=0,fb=0,fa=0,bx=0,ax=0;
-    int n=(int)get_size(p);
-
-
-    ax=0.;
-    xx=1.;
-    
-
-    mnbrak(ax,xx,bx,fa,fx,fb,fadpt);
-    //cout<<xx<<endl;
-    fret=brent(ax,xx,bx,fadpt,TOL,xmin);
-    //cout<<xmin<<endl;
-    for(j=0;j<n;++j)
-      {
-	//get_element(xi,j)*=xmin;
-	set_element(xi,j,
-		    get_element(xi,j)*xmin);
-	//get_element(p,j)+=get_element(xi,j);
-	set_element(p,j,
-		    get_element(p,j)+get_element(xi,j));
-      }
-    //  delete xicom_p;
-    //delete pcom_p;
-  }
-}
-
-
-#endif
diff --git a/methods/powell/mnbrak.hpp b/methods/powell/mnbrak.hpp
deleted file mode 100644
index 9871473..0000000
--- a/methods/powell/mnbrak.hpp
+++ /dev/null
@@ -1,82 +0,0 @@
-#ifndef MNBRAK_HPP
-#define MNBRAK_HPP 
-#define OPT_HEADER
-//#include "optimizer.hpp"
-#include "bas_util.hpp"
-namespace opt_utilities
-{
-
-  
-  template <typename T>
-  void mnbrak(T& ax,T& bx,T& cx,T& fa,T& fb,T& fc,func_obj<T,T>& func)
-  {
-    const T GOLD=1.618034;
-    const T GLIMIT=100;
-    const T TINY=std::numeric_limits<T>::epsilon();
-    T ulim,u,r,q,fu;
-    fa=func.eval(ax);
-    fb=func.eval(bx);
-    
-    if(fb>fa)
-      {
-	//shft(dum,ax,bx,dum);
-	//shft(dum,fb,fa,dum);
-	std::swap(ax,bx);
-	std::swap(fa,fb);
-      }
-
-    cx=bx+GOLD*(bx-ax);
-    fc=func.eval(cx);
-    while(fb>fc)
-      {
-	r=(bx-ax)*(fb-fc);
-	q=(bx-cx)*(fb-fa);
-	u=bx-T((bx-cx)*q-(bx-ax)*r)/
-	  T(T(2.)*sign(T(max(T(tabs(T(q-r))),T(TINY))),T(q-r)));
-	ulim=bx+GLIMIT*(cx-bx);
-	if((bx-u)*(u-cx)>0.)
-	  {
-	    fu=func.eval(u);
-	    if(fu<fc)
-	      {
-		ax=bx;
-		bx=u;
-		fa=fb;
-		fb=fu;
-		return;
-	      }
-	    else if(fu>fb)
-	      {
-		cx=u;
-		fc=fu;
-		return;
-	      }
-	    u=cx+GOLD*(cx-bx);
-	    fu=func.eval(u);
-	  }
-	else if((cx-u)*(u-ulim)>0.)
-	  {
-	    fu=func.eval(u);
-	    if(fu<fc)
-	      {
-		shft3(bx,cx,u,T(cx+GOLD*(cx-bx)));
-		shft3(fb,fc,fu,func.eval(u));
-	      }
-	  }
-	else if((u-ulim)*(ulim-cx)>=0)
-	  {
-	    u=ulim;
-	    fu=func.eval(u);
-	  }
-	else
-	  {
-	    u=cx+GOLD*(cx-bx);
-	    fu=func.eval(u);
-	  }
-	shft3(ax,bx,cx,u);
-	shft3(fa,fb,fc,fu);
-      }
-  }
-}
-
-#endif
diff --git a/methods/powell/powell_method.hpp b/methods/powell/powell_method.hpp
index 1e135b0..270c2b2 100644
--- a/methods/powell/powell_method.hpp
+++ b/methods/powell/powell_method.hpp
@@ -10,7 +10,7 @@
 #include <limits>
 #include <cassert>
 #include <cmath>
-#include "linmin.hpp"
+#include "../linmin/linmin.hpp"
 #include <algorithm>
 #include <iostream>