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.cli;
39
40 import static java.lang.String.format;
41
42 import ffx.algorithms.AlgorithmFunctions;
43 import ffx.algorithms.AlgorithmListener;
44 import ffx.algorithms.optimize.RotamerOptimization;
45 import ffx.crystal.CrystalPotential;
46 import ffx.numerics.Potential;
47 import ffx.potential.DualTopologyEnergy;
48 import ffx.potential.MolecularAssembly;
49 import ffx.potential.QuadTopologyEnergy;
50 import ffx.potential.bonded.LambdaInterface;
51 import ffx.potential.bonded.Polymer;
52 import ffx.potential.bonded.Residue;
53 import ffx.potential.bonded.RotamerLibrary;
54 import ffx.potential.cli.AlchemicalOptions;
55 import ffx.potential.cli.TopologyOptions;
56 import java.io.File;
57 import java.util.ArrayList;
58 import java.util.List;
59 import java.util.logging.Logger;
60 import java.util.regex.Matcher;
61 import java.util.regex.Pattern;
62 import org.apache.commons.configuration2.CompositeConfiguration;
63 import picocli.CommandLine.ArgGroup;
64 import picocli.CommandLine.Option;
65
66
67
68
69
70
71
72
73
74
75 public class MultiDynamicsOptions {
76
77 private static final Logger logger = Logger.getLogger(MultiDynamicsOptions.class.getName());
78
79
80
81
82 @ArgGroup(heading = "%n Multiple Walker Options for MPI Simulations%n", validate = false)
83 private final MultiDynamicsOptionGroup group = new MultiDynamicsOptionGroup();
84
85
86
87
88
89
90
91
92
93
94
95
96 public void distribute(MolecularAssembly[] molecularAssemblies, Potential[] potentials,
97 CrystalPotential crystalPotential, AlgorithmFunctions algorithmFunctions, int rank,
98 int worldSize) {
99 if (!group.distributeWalkersString.equalsIgnoreCase("AUTO")
100 && !group.distributeWalkersString.equalsIgnoreCase("OFF")) {
101 logger.info(" Distributing walker conformations.");
102 int nSys = molecularAssemblies.length;
103 assert nSys == potentials.length;
104 switch (nSys) {
105 case 1 -> optStructure(molecularAssemblies[0], crystalPotential, algorithmFunctions, rank,
106 worldSize);
107 case 2 -> {
108 DualTopologyEnergy dte = (DualTopologyEnergy) crystalPotential;
109 if (dte.getNumSharedVariables() == dte.getNumberOfVariables()) {
110 logger.info(" Generating starting structures based on dual-topology:");
111 optStructure(molecularAssemblies[0], dte, algorithmFunctions, rank, worldSize);
112 } else {
113 logger.info(
114 " Generating separate starting structures for each topology of the dual topology:");
115 optStructure(molecularAssemblies[0], potentials[0], algorithmFunctions, rank, worldSize);
116 optStructure(molecularAssemblies[1], potentials[1], algorithmFunctions, rank, worldSize);
117 }
118 }
119 case 4 -> {
120 QuadTopologyEnergy qte = (QuadTopologyEnergy) crystalPotential;
121 optStructure(molecularAssemblies[0], qte.getDualTopA(), algorithmFunctions, rank,
122 worldSize);
123 optStructure(molecularAssemblies[3], qte.getDualTopB(), algorithmFunctions, rank,
124 worldSize);
125 }
126
127 default -> logger.severe(" First: must have 1, 2, or 4 topologies.");
128 }
129 } else {
130 logger.finer(" Skipping RO-based distribution of initial configurations.");
131 }
132 }
133
134
135
136
137
138
139
140
141
142
143
144 public void distribute(MolecularAssembly[] molecularAssemblies, CrystalPotential crystalPotential,
145 AlgorithmFunctions algorithmFunctions, int rank, int worldSize) {
146 int ntops = molecularAssemblies.length;
147 Potential[] energies = new Potential[ntops];
148 for (int i = 0; i < ntops; i++) {
149 energies[i] = molecularAssemblies[i].getPotentialEnergy();
150 }
151 distribute(molecularAssemblies, energies, crystalPotential, algorithmFunctions, rank, worldSize);
152 }
153
154
155
156
157
158
159
160
161 public boolean isSynchronous() {
162 return group.synchronous;
163 }
164
165 public void setSynchronous(boolean synchronous) {
166 group.synchronous = synchronous;
167 }
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183 public MolecularAssembly openFile(AlgorithmFunctions algorithmFunctions,
184 TopologyOptions topologyOptions, int threadsPer, String toOpen, int topNum,
185 AlchemicalOptions alchemicalOptions, File structureFile, int rank) {
186 boolean autoDist = group.distributeWalkersString.equalsIgnoreCase("AUTO");
187
188 if (autoDist) {
189 String openName = format("%s_%d", toOpen, rank + 1);
190 File testFile = new File(openName);
191 if (testFile.exists()) {
192 toOpen = openName;
193 } else {
194 logger.warning(format(" File %s does not exist; using default %s", openName, toOpen));
195 }
196 }
197 MolecularAssembly assembly = alchemicalOptions.openFile(algorithmFunctions, topologyOptions,
198 threadsPer, toOpen, topNum);
199 assembly.setFile(structureFile);
200 return assembly;
201 }
202
203
204
205
206
207
208 private String[] parseDistributed() {
209 String distributeWalkersString = group.distributeWalkersString;
210 if (distributeWalkersString.equalsIgnoreCase("OFF")
211 || distributeWalkersString.equalsIgnoreCase("AUTO") || distributeWalkersString.isEmpty()) {
212 return null;
213 }
214 return distributeWalkersString.split("\\.");
215 }
216
217
218
219
220
221
222
223
224 public String getDistributeWalkersString() {
225 return group.distributeWalkersString;
226 }
227
228
229
230
231
232
233
234 private void optStructure(MolecularAssembly molecularAssembly, Potential potential,
235 AlgorithmFunctions algorithmFunctions, int rank, int worldSize) {
236 RotamerLibrary rLib = new RotamerLibrary(false);
237 String[] distribRes = parseDistributed();
238
239 if (distribRes == null || distribRes.length == 0) {
240 throw new IllegalArgumentException(
241 " Programming error: Must have list of residues to split on!");
242 }
243
244 LambdaInterface lambdaInterface =
245 (potential instanceof LambdaInterface) ? (LambdaInterface) potential : null;
246 double initLam = -1.0;
247 if (lambdaInterface != null) {
248 initLam = lambdaInterface.getLambda();
249 lambdaInterface.setLambda(0.5);
250 }
251
252 Pattern chainMatcher = Pattern.compile("^([a-zA-Z])?([0-9]+)$");
253
254 List<Residue> residueList = new ArrayList<>(distribRes.length);
255
256 for (String ts : distribRes) {
257 Matcher m = chainMatcher.matcher(ts);
258 Character chainID;
259 int resNum;
260 if (m.find()) {
261 if (m.groupCount() == 2) {
262 chainID = m.group(1).charAt(0);
263 resNum = Integer.parseInt(m.group(2));
264 } else {
265 chainID = ' ';
266 resNum = Integer.parseInt(m.group(1));
267 }
268 } else {
269 logger.warning(format(" Could not parse %s as a valid residue!", ts));
270 continue;
271 }
272 logger.info(format(" Looking for chain %c residue %d", chainID, resNum));
273
274 for (Polymer p : molecularAssembly.getChains()) {
275 if (p.getChainID() == chainID) {
276 for (Residue r : p.getResidues()) {
277 if (r.getResidueNumber() == resNum && r.setRotamers(rLib) != null) {
278 residueList.add(r);
279 }
280 }
281 }
282 }
283 }
284
285 if (residueList.isEmpty()) {
286 throw new IllegalArgumentException(" No valid entries for distWalkers!");
287 }
288
289 AlgorithmListener alist = algorithmFunctions.getDefaultListener();
290 RotamerOptimization ropt = new RotamerOptimization(molecularAssembly, potential, alist);
291 ropt.setRotamerLibrary(rLib);
292
293 ropt.setThreeBodyEnergy(false);
294
295 CompositeConfiguration properties = molecularAssembly.getProperties();
296 if (!properties.containsKey("ro-ensembleNumber") && !properties.containsKey(
297 "ro-ensembleEnergy")) {
298 logger.info(format(" Setting ensemble to default of number of walkers %d", worldSize));
299 ropt.setEnsemble(worldSize);
300 }
301
302 ropt.setPrintFiles(false);
303 ropt.setResiduesIgnoreNull(residueList);
304
305 RotamerLibrary.measureRotamers(residueList, false);
306
307 boolean lazyMat = properties.getBoolean("ro-lazyMatrix", false);
308 if (!lazyMat) {
309 properties.clearProperty("ro-lazyMatrix");
310 properties.addProperty("ro-lazyMatrix", true);
311 }
312
313 ropt.optimize(RotamerOptimization.Algorithm.ALL);
314 ropt.setCoordinatesToEnsemble(rank);
315
316
317
318 double[] xyz = new double[potential.getNumberOfVariables()];
319 potential.getCoordinates(xyz);
320 logger.info(" Final Optimized Energy:");
321 potential.energy(xyz, true);
322
323 if (lambdaInterface != null) {
324 lambdaInterface.setLambda(initLam);
325 }
326
327 if (!lazyMat) {
328 properties.clearProperty("ro-lazyMatrix");
329 }
330 }
331
332 public void setDistributeWalkersString(String distributeWalkersString) {
333 group.distributeWalkersString = distributeWalkersString;
334 }
335
336 public int getFirstDir() {
337 return group.firstDir;
338 }
339
340 public void setFirstDir(int firstDir) {
341 group.firstDir = firstDir;
342 }
343
344
345
346
347 private static class MultiDynamicsOptionGroup {
348
349
350 @Option(names = {"--firstDir"}, defaultValue = "0", paramLabel = "0",
351 description = "The first directory to use for multiple walker jobs.")
352 private int firstDir = 0;
353
354
355 @Option(names = {"-y", "--synchronous"}, defaultValue = "false",
356 description = "Walker communication is synchronous")
357 private boolean synchronous = false;
358
359
360
361
362
363
364 @Option(names = {"--dw", "--distributeWalkers"}, paramLabel = "OFF", defaultValue = "OFF",
365 description = "AUTO: Pick up per-walker configurations as [filename.pdb]_[num], or specify a residue to distribute on.")
366 private String distributeWalkersString = "OFF";
367 }
368 }