1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38 package ffx.algorithms.optimize.anneal;
39
40 import ffx.algorithms.AlgorithmListener;
41 import ffx.algorithms.Terminatable;
42 import ffx.algorithms.dynamics.MolecularDynamics;
43 import ffx.algorithms.dynamics.integrators.IntegratorEnum;
44 import ffx.algorithms.dynamics.thermostats.ThermostatEnum;
45 import ffx.numerics.Potential;
46 import ffx.potential.MolecularAssembly;
47
48 import java.io.File;
49 import java.util.Set;
50 import java.util.logging.Level;
51 import java.util.logging.Logger;
52
53
54
55
56
57
58
59
60 public class SimulatedAnnealing implements Runnable, Terminatable {
61
62 private static final Logger logger = Logger.getLogger(SimulatedAnnealing.class.getName());
63
64 private final MolecularDynamics molecularDynamics;
65
66 private final AnnealingSchedule schedule;
67
68 private final long mdSteps;
69
70 private final double timeStep;
71
72 private final boolean reInitVelocity;
73
74 private final int trajSteps = 1;
75
76 private double printInterval = 0.01;
77
78 private boolean done = true;
79
80 private boolean terminate;
81
82 private double saveFrequency = 0.1;
83
84 private File dynFile;
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100 public SimulatedAnnealing(MolecularAssembly molecularAssembly, Potential potentialEnergy,
101 AlgorithmListener algorithmListener, ThermostatEnum requestedThermostat,
102 IntegratorEnum requestedIntegrator, AnnealingSchedule annealingSchedule, long mdSteps,
103 double timeStep, boolean reInitVelocity, File dynFile) {
104
105 molecularDynamics = MolecularDynamics.dynamicsFactory(molecularAssembly, potentialEnergy,
106 algorithmListener, requestedThermostat, requestedIntegrator);
107 this.schedule = annealingSchedule;
108 this.mdSteps = mdSteps;
109 this.timeStep = timeStep;
110 this.reInitVelocity = reInitVelocity;
111 this.dynFile = dynFile;
112 }
113
114
115 public void anneal() {
116
117
118 if (!done) {
119 logger.warning(" Programming error - a thread invoked anneal when it was already running.");
120 return;
121 }
122 done = false;
123 logger.info(" Beginning simulated annealing");
124 begin();
125 }
126
127
128
129
130
131
132 public double getKineticEnergy() {
133 return molecularDynamics.getKineticEnergy();
134 }
135
136
137
138
139
140
141 public double getPotentialEnergy() {
142 return molecularDynamics.getPotentialEnergy();
143 }
144
145
146
147
148
149
150 public double getTemperature() {
151 return molecularDynamics.getTemperature();
152 }
153
154
155
156
157
158
159 public double getTotalEnergy() {
160 return molecularDynamics.getTotalEnergy();
161 }
162
163
164
165
166
167
168 @Override
169 public void run() {
170 done = false;
171 terminate = false;
172
173 int minMdSteps = (int) (mdSteps * schedule.minWindowLength());
174 if (minMdSteps < trajSteps) {
175 logger.warning(String.format(
176 " Minimum number of MD steps per annealing cycle %d was less than steps per OpenMM MD cycle %d! Setting steps per MD cycle to %d",
177 minMdSteps, trajSteps, minMdSteps));
178 setTrajectorySteps(minMdSteps);
179 }
180
181 int nWindows = schedule.getNumWindows();
182 boolean forceFirstReinit = (dynFile == null);
183
184 for (int i = 0; i < nWindows; i++) {
185 double temperature = schedule.getTemperature(i);
186 int nSteps = (int) (schedule.windowLength(i) * mdSteps);
187 logger.info(
188 String.format(" Annealing window %d: %d steps at %9.4g K", (i + 1), nSteps, temperature));
189 molecularDynamics.dynamic(nSteps, timeStep, printInterval, saveFrequency, temperature,
190 (reInitVelocity || forceFirstReinit), dynFile);
191 if (dynFile == null) {
192 dynFile = molecularDynamics.getDynFile();
193 }
194 forceFirstReinit = false;
195 if (terminate) {
196 logger.info(String.format("\n Terminating at temperature %8.3f.\n", temperature));
197 break;
198 }
199 }
200 if (!terminate) {
201 logger.info(String.format(" Completed %8d annealing steps\n", nWindows));
202 }
203
204 done = true;
205 terminate = false;
206 }
207
208
209
210
211
212
213 public void setPrintInterval(double printInterval) {
214 this.printInterval = printInterval;
215 }
216
217
218
219
220
221
222
223 public void setRestartFrequency(double restart) throws IllegalArgumentException {
224 if (Double.isFinite(restart) && restart > 0) {
225 molecularDynamics.setRestartFrequency(restart);
226 } else {
227 throw new IllegalArgumentException(
228 String.format(" Restart frequency must be positive finite, was %10.4g", restart));
229 }
230 }
231
232
233
234
235
236
237 public void setSaveFrequency(double save) {
238 this.saveFrequency = save;
239 }
240
241
242
243
244
245
246 public void setTrajectorySteps(int trajectorySteps) {
247 molecularDynamics.setIntervalSteps(trajectorySteps);
248 }
249
250
251 @Override
252 public void terminate() {
253 terminate = true;
254 while (!done) {
255 synchronized (this) {
256 try {
257 wait(1);
258 } catch (Exception e) {
259 logger.log(Level.WARNING, "Exception terminating annealing.\n", e);
260 }
261 }
262 }
263 }
264
265 private void begin() {
266 logger.info(String.format(" Initial temperature: %8.3f (Kelvin)", schedule.getHighTemp()));
267 logger.info(String.format(" Final temperature: %8.3f (Kelvin)", schedule.getLowTemp()));
268 logger.info(String.format(" Annealing steps: %8d", schedule.getNumWindows()));
269 logger.info(String.format(" MD steps/temperature: %8d", mdSteps));
270 logger.info(String.format(" MD time step: %8.3f (fs)", timeStep));
271
272 Thread annealingThread = new Thread(this);
273 annealingThread.start();
274 synchronized (this) {
275 try {
276 while (annealingThread.isAlive()) {
277 wait(100);
278 }
279 } catch (Exception e) {
280 String message = "Simulated annealing interrupted.";
281 logger.log(Level.WARNING, message, e);
282 }
283 }
284 }
285
286
287 public enum Schedules {
288 EXP(ExpAnnealSchedule::new, "EXP", "EXPONENTIAL"), LINEAR(LinearAnnealSchedule::new, "LINEAR");
289
290 private final ScheduleConstructor sc;
291 private final Set<String> aliases;
292
293 Schedules(ScheduleConstructor sc, String... names) {
294 this.sc = sc;
295 aliases = Set.of(names);
296 }
297
298
299
300
301
302
303
304 public static Schedules parse(String name) {
305 name = name.toUpperCase();
306 for (Schedules s : values()) {
307 if (s.aliases.contains(name)) {
308 return s;
309 }
310 }
311 return valueOf(name);
312 }
313
314
315
316
317
318
319
320
321
322 public AnnealingSchedule generate(int nWindows, double tLow, double tHigh) {
323 return sc.asConstruct(nWindows, tLow, tHigh);
324 }
325 }
326
327
328 private interface ScheduleConstructor {
329
330 AnnealingSchedule asConstruct(int nWindows, double tLow, double tHigh);
331 }
332 }