00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040 #include "LAStreams.h"
00041 #include "math_num.h"
00042 #include "std.h"
00043
00044 namespace linalg
00045 {
00046 using namespace linalg;
00047
00048
00049
00050
00051
00052
00053 class FPoint
00054 {
00055 Vector& x;
00056 double fval;
00057 MultivariateFunctor& fproc;
00058 const bool free_x_on_destructing;
00059
00060
00061
00062 public:
00063 FPoint(Vector& b, MultivariateFunctor& f);
00064 FPoint(const FPoint& fp);
00065 ~FPoint(void);
00066
00067 FPoint& operator = (const FPoint& fp);
00068
00069 double f(void) const { return fval; }
00070
00071 double fiddle_around(LAStreamIn h);
00072
00073
00074
00075
00076
00077 friend void update_in_direction(FPoint& from, FPoint& to);
00078
00079
00080
00081 bool is_step_relatively_small(LAStreamIn hs, const double tau);
00082 };
00083
00084
00085 inline FPoint::FPoint(Vector& b, MultivariateFunctor& f)
00086 : x(b), fval(f(b)), fproc(f), free_x_on_destructing(false)
00087 {
00088 }
00089
00090
00091 FPoint::FPoint(const FPoint& fp)
00092 : x(*(new Vector(fp.x))),
00093 fval(fp.fval), fproc(fp.fproc),
00094 free_x_on_destructing(true)
00095 {
00096 }
00097
00098 FPoint::~FPoint(void)
00099 {
00100 if( free_x_on_destructing )
00101 delete &x;
00102 }
00103
00104
00105
00106 inline FPoint& FPoint::operator = (const FPoint& fp)
00107 {
00108 x = fp.x;
00109 fval = fp.fval;
00110 return *this;
00111 }
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122 double FPoint::fiddle_around(LAStreamIn hs)
00123 {
00124
00125 for(LAStreamOut xs(x); !xs.eof(); xs.get() )
00126 {
00127 const double hi = hs.get();
00128 register const double xi_old = xs.peek();
00129 register double fnew;
00130
00131 if( xs.peek() = xi_old + hi, (fnew = fproc(x)) < fval )
00132 fval = fnew;
00133 else if( xs.peek() = xi_old - hi, (fnew = fproc(x)) < fval )
00134 fval = fnew;
00135 else
00136 xs.peek() = xi_old;
00137 }
00138 return fval;
00139 }
00140
00141
00142
00143
00144
00145 void update_in_direction(FPoint& from, FPoint& to)
00146 {
00147 for(LAStreamOut tos(to.x), froms(from.x); !tos.eof(); )
00148 {
00149 register const double t = tos.peek();
00150 tos.get() += (t - froms.peek());
00151 froms.get() = t;
00152 }
00153 from.fval = to.fval;
00154 to.fval = (to.fproc)(to.x);
00155 }
00156
00157
00158
00159 bool FPoint::is_step_relatively_small(LAStreamIn hs, const double tau)
00160 {
00161 for(LAStreamIn xs(x); !hs.eof(); )
00162 if( !(hs.get() /(1 + fabs(xs.get())) < tau) )
00163 return false;
00164 return true;
00165 }
00166
00167
00168
00169
00170
00171 #if 0
00172 double hjmin(Vector& b, Vector& h, MultivariateFunctor& ff)
00173 {
00174
00175 const double tau = 10*DBL_EPSILON;
00176 const double threshold = 1e-8;
00177
00178
00179 const double step_reduce_factor = 10;
00180
00181 are_compatible(b,h);
00182
00183 FPoint pmin(b,ff);
00184 FPoint pbase(pmin);
00185
00186 for(;;)
00187 {
00188 if( pbase.fiddle_around(h) < pmin.f() - threshold )
00189 {
00190 do
00191 update_in_direction(pmin,pbase);
00192 while( pbase.fiddle_around(h) < pmin.f() - threshold );
00193 pbase = pmin;
00194 }
00195 else
00196 if(
00197 h *= 1/step_reduce_factor,
00198 pbase.is_step_relatively_small(h,tau) )
00199 return pmin.f();
00200 }
00201 }
00202
00203
00204
00205
00206
00207 double hjmin(Vector& b, const double h0, MultivariateFunctor& f)
00208 {
00209 Vector h(b.q_lwb(),b.q_upb()); h = h0;
00210 return hjmin(b,h,f);
00211 }
00212 #endif // of 0
00213 }