samer@0
|
1 /*
|
samer@0
|
2 * Copyright (c) 2000, Samer Abdallah, King's College London.
|
samer@0
|
3 * All rights reserved.
|
samer@0
|
4 *
|
samer@0
|
5 * This software is provided AS iS and WITHOUT ANY WARRANTY;
|
samer@0
|
6 * without even the implied warranty of MERCHANTABILITY or
|
samer@0
|
7 * FITNESS FOR A PARTICULAR PURPOSE.
|
samer@0
|
8 */
|
samer@0
|
9
|
samer@0
|
10 package samer.maths.opt;
|
samer@0
|
11 import samer.maths.*;
|
samer@0
|
12 import samer.core.*;
|
samer@0
|
13
|
samer@0
|
14 public class State
|
samer@0
|
15 {
|
samer@0
|
16 public Functionx func; // function
|
samer@0
|
17 public int n; // number of dimensions
|
samer@0
|
18 public double[] h; // search direction
|
samer@0
|
19 public Datum P1, P2; // two bracket points and a new point
|
samer@0
|
20 public double alpha; // step length
|
samer@0
|
21 public double normh; // some norm of h
|
samer@0
|
22
|
samer@0
|
23 /** x must have an accessible array */
|
samer@0
|
24 public State( Vec x, Functionx func)
|
samer@0
|
25 {
|
samer@0
|
26 this.n = x.size();
|
samer@0
|
27 this.func = func;
|
samer@0
|
28 P1 = new Datum(n,x.array());
|
samer@0
|
29 P2 = new Datum(n);
|
samer@0
|
30 h = new double[n];
|
samer@0
|
31 }
|
samer@0
|
32
|
samer@0
|
33 public void dispose()
|
samer@0
|
34 {
|
samer@0
|
35 func.dispose();
|
samer@0
|
36 P1.dispose();
|
samer@0
|
37 P2.dispose();
|
samer@0
|
38 h=null;
|
samer@0
|
39 }
|
samer@0
|
40
|
samer@0
|
41 /** estimated step to minimum using cubic interpolation of current state */
|
samer@0
|
42 public final double cubic() { return Util.cubici(P1.f,P1.s,P2.f,P2.s,alpha); }
|
samer@0
|
43
|
samer@0
|
44 /** evaluate at P1 */
|
samer@0
|
45 public final void evaluate() { func.evaluate(P1); }
|
samer@0
|
46
|
samer@0
|
47 /** move P1 to where P2 is now */
|
samer@0
|
48 public final void move() { P1.copy(P2); } // assert that P1.s<0 ?
|
samer@0
|
49
|
samer@0
|
50 /** take step of length lambda from
|
samer@0
|
51 P1 in direction of h, evaluate and store in P2 */
|
samer@0
|
52 public void step(double lambda)
|
samer@0
|
53 {
|
samer@0
|
54 for (int i=0; i<n; i++) {
|
samer@0
|
55 P2.x[i] = P1.x[i] + lambda*h[i];
|
samer@0
|
56 }
|
samer@0
|
57
|
samer@0
|
58 func.evaluate(P2);
|
samer@0
|
59 P2.s = Mathx.dot(P2.g,h);
|
samer@0
|
60 alpha=lambda;
|
samer@0
|
61 }
|
samer@0
|
62
|
samer@0
|
63 /** evaluate at P2 */
|
samer@0
|
64 public final void eval2()
|
samer@0
|
65 {
|
samer@0
|
66 func.evaluate(P2);
|
samer@0
|
67 P2.s = Mathx.dot(P2.g,h);
|
samer@0
|
68 }
|
samer@0
|
69
|
samer@0
|
70 /** evaluate slope at P1 using gradient at P1 */
|
samer@0
|
71 public final void setSlope() {
|
samer@0
|
72 double s1 = Mathx.dot(P1.g,h);
|
samer@0
|
73 if (Double.isNaN(s1)) throw new Error( "slope 1 is NaN");
|
samer@0
|
74 if (s1>0) {
|
samer@0
|
75 Shell.trace(toString());
|
samer@0
|
76 Shell.trace("bad new slope="+s1);
|
samer@0
|
77 }
|
samer@0
|
78 P1.s = s1;
|
samer@0
|
79 }
|
samer@0
|
80
|
samer@0
|
81 /** set slope at P2, assert new slope is negative */
|
samer@0
|
82 public final void setSlope2(double s2)
|
samer@0
|
83 {
|
samer@0
|
84 if (Double.isNaN(s2)) throw new Error( "slope 2 is NaN");
|
samer@0
|
85
|
samer@0
|
86 if (s2>0) {
|
samer@0
|
87 Shell.trace(" ***** BAD DIRECTION UPDATE **** ");
|
samer@0
|
88 Shell.trace(toString());
|
samer@0
|
89 Shell.trace("-- slope = "+s2);
|
samer@0
|
90 }
|
samer@0
|
91
|
samer@0
|
92 P2.s=s2;
|
samer@0
|
93 }
|
samer@0
|
94
|
samer@0
|
95 public String toString() {
|
samer@0
|
96 return append(new StringBuffer()).toString();
|
samer@0
|
97 }
|
samer@0
|
98
|
samer@0
|
99 public double initialStep() {
|
samer@0
|
100 double step=2*Math.abs(P1.f/P1.s);
|
samer@0
|
101 if (Double.isNaN(step)) return 1;
|
samer@0
|
102 return Util.clip(0.01,step,1);
|
samer@0
|
103 }
|
samer@0
|
104
|
samer@0
|
105 /** estimate initial step for line search after turn */
|
samer@0
|
106 public double nextStep()
|
samer@0
|
107 {
|
samer@0
|
108 double beta;
|
samer@0
|
109
|
samer@0
|
110 if (P2.s<=P1.s) return 1;
|
samer@0
|
111 if (P2.s>0) {
|
samer@0
|
112 // signal: 'overshoot, updhess'
|
samer@0
|
113 // expect beta between 0 and 1
|
samer@0
|
114 if (alpha<0.9) beta=alpha;
|
samer@0
|
115 else beta=cubic();
|
samer@0
|
116 if (Double.isNaN(beta) || Double.isInfinite(beta)) beta=1;
|
samer@0
|
117 } else {
|
samer@0
|
118 // signal: 'undershoot, updhess'
|
samer@0
|
119 double gamma=1.5*alpha; if (gamma<2) gamma=2;
|
samer@0
|
120 double delta=alpha*(P2.s-P1.s)-P2.s+(alpha>1?alpha:1);
|
samer@0
|
121 // double delta=uv-P2.s+(alpha>1?alpha:1);
|
samer@0
|
122 beta=cubic();
|
samer@0
|
123 if (Double.isNaN(beta) || Double.isInfinite(beta) || beta<alpha) beta=2*alpha+1e-5;
|
samer@0
|
124 beta*=1.2;
|
samer@0
|
125 if (gamma<beta) beta=gamma;
|
samer@0
|
126 if (delta<beta) beta=delta;
|
samer@0
|
127 }
|
samer@0
|
128
|
samer@0
|
129 return beta;
|
samer@0
|
130 }
|
samer@0
|
131
|
samer@0
|
132 public StringBuffer append(StringBuffer info)
|
samer@0
|
133 {
|
samer@0
|
134 info.append(DoubleFormat.format(P1.f,4))
|
samer@0
|
135 .append(", ")
|
samer@0
|
136 .append(DoubleFormat.format(P1.s,4))
|
samer@0
|
137 .append("\t")
|
samer@0
|
138 .append(DoubleFormat.format(alpha,4))
|
samer@0
|
139 .append("\t")
|
samer@0
|
140 .append(DoubleFormat.format(P2.f,4))
|
samer@0
|
141 .append(", ")
|
samer@0
|
142 .append(DoubleFormat.format(P2.s,4))
|
samer@0
|
143 .append(" : "); // .append(this.info);
|
samer@0
|
144 return info;
|
samer@0
|
145 }
|
samer@0
|
146 }
|