samer@0: /*
samer@0:  *	Copyright (c) 2000, Samer Abdallah, King's College London.
samer@0:  *	All rights reserved.
samer@0:  *
samer@0:  *	This software is provided AS iS and WITHOUT ANY WARRANTY; 
samer@0:  *	without even the implied warranty of MERCHANTABILITY or 
samer@0:  *	FITNESS FOR A PARTICULAR PURPOSE.
samer@0:  */
samer@0: 
samer@0: package samer.maths.opt;
samer@0: import  samer.maths.*;
samer@0: import  samer.core.*;
samer@0: 
samer@0: public class State
samer@0: {
samer@0: 	public Functionx	func;			// function
samer@0: 	public int			n;				// number of dimensions
samer@0: 	public double[]	h;				// search direction
samer@0: 	public Datum		P1, P2;		// two bracket points and a new point
samer@0: 	public double		alpha;		// step length
samer@0: 	public double		normh;		// some norm of h
samer@0: 
samer@0: 	/** x must have an accessible array */
samer@0: 	public State( Vec x, Functionx func)
samer@0: 	{
samer@0: 		this.n = x.size();
samer@0: 		this.func = func;
samer@0: 		P1 = new Datum(n,x.array());
samer@0: 		P2 = new Datum(n);
samer@0: 		h = new double[n];
samer@0: 	}
samer@0: 
samer@0: 	public void dispose() 
samer@0: 	{ 
samer@0: 		func.dispose(); 
samer@0: 		P1.dispose(); 
samer@0: 		P2.dispose(); 
samer@0: 		h=null; 
samer@0: 	}
samer@0: 
samer@0: 	/** estimated step to minimum using cubic interpolation of current state */
samer@0: 	public final double cubic() { return Util.cubici(P1.f,P1.s,P2.f,P2.s,alpha); }
samer@0: 
samer@0: 	/** evaluate at P1 */
samer@0: 	public final void evaluate() { func.evaluate(P1);	}
samer@0: 
samer@0: 	/** move P1 to where P2 is now */
samer@0: 	public final void move() { P1.copy(P2); } // assert that P1.s<0 ?
samer@0: 
samer@0: 	/** take step of length lambda from 
samer@0: 	    P1 in direction of h, evaluate and store in P2 */
samer@0: 	public void step(double lambda)
samer@0: 	{
samer@0: 		for (int i=0; i<n; i++) {
samer@0: 			P2.x[i] = P1.x[i] + lambda*h[i];
samer@0: 		}
samer@0: 
samer@0: 		func.evaluate(P2);
samer@0: 		P2.s = Mathx.dot(P2.g,h);
samer@0: 		alpha=lambda;
samer@0: 	}
samer@0: 
samer@0: 	/** evaluate at P2 */
samer@0: 	public final void eval2()
samer@0: 	{
samer@0: 		func.evaluate(P2);
samer@0: 		P2.s = Mathx.dot(P2.g,h);
samer@0: 	}
samer@0: 
samer@0: 	/** evaluate slope at P1 using gradient at P1 */	
samer@0: 	public final void setSlope() { 
samer@0: 		double s1 = Mathx.dot(P1.g,h); 
samer@0: 		if (Double.isNaN(s1)) throw new Error( "slope 1 is NaN");
samer@0: 		if (s1>0) {
samer@0: 			Shell.trace(toString());
samer@0: 			Shell.trace("bad new slope="+s1);
samer@0: 		}
samer@0: 		P1.s = s1;
samer@0: 	}
samer@0: 
samer@0: 	/** set slope at P2, assert new slope is negative */	
samer@0: 	public final void setSlope2(double s2)
samer@0: 	{
samer@0: 		if (Double.isNaN(s2)) throw new Error( "slope 2 is NaN");
samer@0: 
samer@0: 		if (s2>0) {
samer@0: 			Shell.trace(" ***** BAD DIRECTION UPDATE **** ");
samer@0: 			Shell.trace(toString());
samer@0: 			Shell.trace("-- slope = "+s2);
samer@0: 		}
samer@0: 	
samer@0: 		P2.s=s2;
samer@0: 	}
samer@0: 
samer@0: 	public String toString() {
samer@0: 		return append(new StringBuffer()).toString();
samer@0: 	}
samer@0: 	
samer@0: 	public double initialStep() {
samer@0: 		double step=2*Math.abs(P1.f/P1.s);
samer@0: 		if (Double.isNaN(step)) return 1;
samer@0: 		return Util.clip(0.01,step,1);
samer@0: 	}
samer@0: 
samer@0: 	/** estimate initial step for line search after turn */
samer@0: 	public double nextStep()
samer@0: 	{
samer@0: 		double beta;
samer@0: 
samer@0: 		if (P2.s<=P1.s) return 1;		
samer@0: 		if (P2.s>0) {
samer@0: 			// signal: 'overshoot, updhess'
samer@0: 			// expect beta between 0 and 1
samer@0: 			if (alpha<0.9) beta=alpha;
samer@0: 			else beta=cubic();
samer@0: 			if (Double.isNaN(beta) || Double.isInfinite(beta)) beta=1;
samer@0: 		} else {
samer@0: 			// signal: 'undershoot, updhess'
samer@0: 			double gamma=1.5*alpha; if (gamma<2) gamma=2;
samer@0: 			double delta=alpha*(P2.s-P1.s)-P2.s+(alpha>1?alpha:1);
samer@0: 			// double delta=uv-P2.s+(alpha>1?alpha:1);
samer@0: 			beta=cubic();
samer@0: 			if (Double.isNaN(beta) || Double.isInfinite(beta) || beta<alpha) beta=2*alpha+1e-5;
samer@0: 			beta*=1.2;
samer@0: 			if (gamma<beta) beta=gamma;
samer@0: 			if (delta<beta) beta=delta;
samer@0: 		}
samer@0: 
samer@0: 		return beta;
samer@0: 	}
samer@0: 
samer@0: 	public StringBuffer append(StringBuffer info) 
samer@0: 	{
samer@0: 		info.append(DoubleFormat.format(P1.f,4))
samer@0: 			.append(", ")
samer@0: 			.append(DoubleFormat.format(P1.s,4))
samer@0: 			.append("\t")
samer@0: 			.append(DoubleFormat.format(alpha,4))
samer@0: 			.append("\t")
samer@0: 			.append(DoubleFormat.format(P2.f,4))
samer@0: 			.append(", ")
samer@0: 			.append(DoubleFormat.format(P2.s,4))
samer@0: 			.append(" : "); // .append(this.info);
samer@0: 		return info;		
samer@0: 	}
samer@0: }