1 |
|
|
2 |
from Numeric import * |
from Numeric import * |
3 |
|
import RandomArray |
4 |
from Tkinter import * |
from Tkinter import * |
5 |
import math |
import math |
6 |
import os |
import os |
41 |
geom += "%s %s %s\n" % (x[i], y[i], z[i]) |
geom += "%s %s %s\n" % (x[i], y[i], z[i]) |
42 |
return geom |
return geom |
43 |
|
|
44 |
|
def flatsum(x): |
45 |
|
return add.reduce(ravel(x)) |
46 |
|
|
47 |
|
class LineMin: |
48 |
|
def __init__(self): |
49 |
|
self.l = 0.01 |
50 |
|
def linemin(self, func, x, dx): |
51 |
|
(u, g) = func(x) |
52 |
|
origu = u |
53 |
|
origg = g |
54 |
|
# print "Start linemin at ",u |
55 |
|
dot = -flatsum(g * dx) |
56 |
|
if dot < 0: |
57 |
|
print "Inverting: ",dot |
58 |
|
dx = dx * -1 |
59 |
|
dot = -dot |
60 |
|
l = 0 |
61 |
|
while dot > 0: |
62 |
|
prevl = l |
63 |
|
l += self.l |
64 |
|
l *= 2 |
65 |
|
(uc, gc) = func(x + l * dx) |
66 |
|
dot = -flatsum(gc * dx) |
67 |
|
# print "Lextend: ",l,dot |
68 |
|
curmax = l # Point where gradient dot was known to be positive |
69 |
|
curmin = prevl # Point where gradient dot was known to be negative |
70 |
|
while (curmax-curmin)/(curmax+curmin) > 0.05: |
71 |
|
testl = 0.5*(curmax+curmin) |
72 |
|
(uc, gc) = func(x + testl * dx) |
73 |
|
# print "Lloop: ",curmin,curmax,uc |
74 |
|
dot = -flatsum(gc * dx) |
75 |
|
if dot < 0: |
76 |
|
curmax = testl |
77 |
|
else: |
78 |
|
curmin = testl |
79 |
|
final = (0.5 * (curmax +curmin)) |
80 |
|
(tu, tg) = func(x + final * dx) |
81 |
|
self.l = final |
82 |
|
if tu > (u * 1.03): # Allow small jiggle to help out of local minima |
83 |
|
print "ARGH! LINE SEARCH FAILURE: ",u, tu, self.l, final, g, tg |
84 |
|
return x, origu, origg |
85 |
|
return x + final * dx, tu, tg |
86 |
|
|
87 |
|
class ConjGrad: |
88 |
|
def __init__(self, func, x): |
89 |
|
self.func = func |
90 |
|
self.x = x |
91 |
|
self.lin = LineMin() |
92 |
|
(u,g) = func(x) |
93 |
|
self.grad = g |
94 |
|
self.dir = -g |
95 |
|
def round(self, step, brake): |
96 |
|
(self.x, up, gp) = self.lin.linemin(self.func, self.x, self.dir) |
97 |
|
# Polak-Ribiere |
98 |
|
beta = flatsum(gp * (gp-self.grad)) / flatsum(self.grad**2); |
99 |
|
self.dir = -gp + beta * self.dir |
100 |
|
self.grad = gp |
101 |
|
return up |
102 |
|
def kick(self): |
103 |
|
(u,g) = self.func(self.x) |
104 |
|
self.x += 0.01 * (RandomArray.random(g.shape)-0.5) * g |
105 |
|
(u,g) = self.func(self.x) |
106 |
|
self.grad = g |
107 |
|
self.dir = -g |
108 |
|
|
109 |
|
|
110 |
|
|
111 |
|
|
112 |
|
|
113 |
class MomGrad: |
class MomGrad: |
114 |
def __init__(self, func, x): |
def __init__(self, func, x): |
115 |
self.func = func |
self.func = func |
116 |
self.mom = x * 0 |
self.mom = x * 0 |
117 |
self.x = x |
self.x = x |
118 |
|
self.lin = LineMin() |
119 |
def round(self, step, brake): |
def round(self, step, brake): |
120 |
|
|
121 |
(u, g) = self.func(self.x) |
(u, g) = self.func(self.x) |
122 |
self.u = u |
self.u = u |
123 |
self.mom += g |
self.mom -= g |
124 |
l = add.reduce(ravel(self.mom*self.mom)) |
# l = flatsum(self.mom*self.mom) |
125 |
stp = -step / sqrt(l) * self.mom |
# stp = -step / sqrt(l) * self.mom |
126 |
self.x += stp |
# self.x += stp |
|
# self.x = maximum(self.x, 0) |
|
127 |
self.mom *= brake |
self.mom *= brake |
128 |
|
|
129 |
|
self.x, foo = self.lin.linemin(self.func, self.x, self.mom) |
130 |
# print "Mom, l:", self.mom, l |
# print "Mom, l:", self.mom, l |
131 |
return self.u |
return self.u |
132 |
|
def kick(self): |
133 |
|
(u,g) = self.func(self.x) |
134 |
|
self.x += 0.01 * (RandomArray.random(g.shape)-0.5) * g |
135 |
|
|
136 |
|
|
137 |
|
|
138 |
def slider(obj, name, var, default, from_, to, resolution=0.001): |
def slider(obj, name, var, default, from_, to, resolution=0.0001, parent=None): |
139 |
|
if parent == None: parent = obj |
140 |
setattr(obj, name, var) |
setattr(obj, name, var) |
141 |
var.set(default) |
var.set(default) |
142 |
Scale(obj, label=name, orient="horizontal", |
Scale(parent, label=name, orient="horizontal", |
143 |
resolution=resolution, |
resolution=resolution, |
144 |
from_=from_, to=to, variable=var).pack( |
from_=from_, to=to, variable=var).pack( |
145 |
fill="both", expand=1) |
fill="both", expand=1) |
151 |
Button(self, text="RESET", command=self.reset).pack() |
Button(self, text="RESET", command=self.reset).pack() |
152 |
self.pack(expand=1, fill="both") |
self.pack(expand=1, fill="both") |
153 |
Button(self, text="Start", command = self.idle).pack() |
Button(self, text="Start", command = self.idle).pack() |
154 |
|
Button(self, text="kick", command = self.kick).pack() |
155 |
# Button(self, text="Kick", command = self.stretch.kick).pack() |
# Button(self, text="Kick", command = self.stretch.kick).pack() |
156 |
slider(self, "dice", IntVar(), 10, 1, 100, resolution=1); |
slider(self, "dice", IntVar(), 20, 1, 100, resolution=1); |
157 |
slider(self, "niter", IntVar(), 5, 1, 50, resolution=1); |
slider(self, "niter", IntVar(), 1, 1, 50, resolution=1); |
158 |
|
|
159 |
slider(self, "mom_step", DoubleVar(), 0.03, 0, 0.2) |
slider(self, "mom_step", DoubleVar(), 0.03, 0, 0.2) |
160 |
slider(self, "mom_brake", DoubleVar(), 0.97, 0.8, 1.0) |
slider(self, "mom_brake", DoubleVar(), 0.97, 0.8, 1.0) |
161 |
|
|
162 |
|
self.problem.addWidgets(self) |
163 |
|
|
164 |
self.pot = DoubleVar(); |
self.pot = DoubleVar(); |
165 |
Label(self, textvariable=self.pot).pack() |
Label(self, textvariable=self.pot).pack() |
166 |
|
|
167 |
|
|
168 |
self.reset() |
self.reset() |
169 |
self.i = 0 |
self.i = 0 |
170 |
|
|
171 |
def reset(self): |
def reset(self): |
172 |
self.problem.reset(self.dice.get()) |
self.problem.reset(self.dice.get()) |
173 |
self.opt = MomGrad(self.problem, self.problem.defaultX()) |
self.opt = ConjGrad(self.problem, self.problem.defaultX()) |
174 |
|
def kick(self): |
175 |
|
self.opt.kick() |
176 |
def idle(self): |
def idle(self): |
177 |
# print "IDLE" |
# print "IDLE" |
178 |
u = self.opt.round(step = self.mom_step.get(), brake = self.mom_brake.get()) |
u = self.opt.round(step = self.mom_step.get(), brake = self.mom_brake.get()) |