View Javadoc
1   // ******************************************************************************
2   //
3   // Title:       Force Field X.
4   // Description: Force Field X - Software for Molecular Biophysics.
5   // Copyright:   Copyright (c) Michael J. Schnieders 2001-2024.
6   //
7   // This file is part of Force Field X.
8   //
9   // Force Field X is free software; you can redistribute it and/or modify it
10  // under the terms of the GNU General Public License version 3 as published by
11  // the Free Software Foundation.
12  //
13  // Force Field X is distributed in the hope that it will be useful, but WITHOUT
14  // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
15  // FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
16  // details.
17  //
18  // You should have received a copy of the GNU General Public License along with
19  // Force Field X; if not, write to the Free Software Foundation, Inc., 59 Temple
20  // Place, Suite 330, Boston, MA 02111-1307 USA
21  //
22  // Linking this library statically or dynamically with other modules is making a
23  // combined work based on this library. Thus, the terms and conditions of the
24  // GNU General Public License cover the whole combination.
25  //
26  // As a special exception, the copyright holders of this library give you
27  // permission to link this library with independent modules to produce an
28  // executable, regardless of the license terms of these independent modules, and
29  // to copy and distribute the resulting executable under terms of your choice,
30  // provided that you also meet, for each linked independent module, the terms
31  // and conditions of the license of that module. An independent module is a
32  // module which is not derived from or based on this library. If you modify this
33  // library, you may extend this exception to your version of the library, but
34  // you are not obligated to do so. If you do not wish to do so, delete this
35  // exception statement from your version.
36  //
37  // ******************************************************************************
38  package ffx.algorithms.cli;
39  
40  import static java.lang.Double.parseDouble;
41  import static java.lang.Integer.parseInt;
42  import static java.lang.String.format;
43  
44  import ffx.algorithms.dynamics.Barostat;
45  import ffx.crystal.CrystalPotential;
46  import ffx.potential.MolecularAssembly;
47  
48  import java.util.logging.Logger;
49  
50  import picocli.CommandLine.ArgGroup;
51  import picocli.CommandLine.Option;
52  
53  /**
54   * Represents command line options for scripts that use a barostat/NPT.
55   *
56   * @author Michael J. Schnieders
57   * @author Hernan V. Bernabe
58   * @since 1.0
59   */
60  public class BarostatOptions {
61  
62    /**
63     * Default maximum density constraint on the barostat that prevents reduction in unit cell
64     * (particularly at or near vapor states).
65     */
66    public static final String DEFAULT_MAX_DENSITY = "1.6";
67    /**
68     * Default "tin box" constraint on the barostat that prevents expansion of the unit cell
69     * (particularly at or near vapor states).
70     */
71    public static final String DEFAULT_MIN_DENSITY = "0.75";
72    /**
73     * Default width of proposed crystal angle moves (uniformly distributed) in degrees.
74     */
75    public static final String DEFAULT_MAX_ANGLE_MOVE = "0.5";
76    /**
77     * Default size of proposed unit cell volume moves (uniformly distributed) in Angstroms^3.
78     */
79    public static final String DEFAULT_MAX_VOLUME_MOVE = "1.0";
80    /**
81     * Default mean number of MD steps (Poisson distribution) between barostat move proposals.
82     */
83    public static final String DEFAULT_BAROSTAT_INTERVAL = "10";
84    /**
85     * Default number of Barostat moves between print statements.
86     */
87    public static final String DEFAULT_BAROSTAT_PRINT_INTERVAL = "1000";
88  
89    private static final Logger logger = Logger.getLogger(BarostatOptions.class.getName());
90  
91    /**
92     * The ArgGroup keeps the BarostatOptions together when printing help.
93     */
94    @ArgGroup(heading = "%n Monte Carlo Pressure Options%n", validate = false)
95    private final BarostatOptionGroup group = new BarostatOptionGroup();
96  
97    /**
98     * If pressure has been set > 0, creates a Barostat around a CrystalPotential, else returns the
99     * original, unmodified CrystalPotential.
100    *
101    * @param molecularAssembly Primary molecularAssembly of the CrystalPotential.
102    * @param crystalPotential  A CrystalPotential.
103    * @return Either a Barostat (NPT enabled) or the CrystalPotential.
104    */
105   public CrystalPotential checkNPT(MolecularAssembly molecularAssembly,
106                                    CrystalPotential crystalPotential) {
107     if (group.pressure > 0) {
108       return createBarostat(molecularAssembly, crystalPotential);
109     } else {
110       return crystalPotential;
111     }
112   }
113 
114   /**
115    * Creates a Barostat around a CrystalPotential.
116    *
117    * @param assembly         Primary assembly of the CrystalPotential.
118    * @param crystalPotential A CrystalPotential.
119    * @return An NPT potential.
120    * @throws IllegalArgumentException If this BarostatOptions has pressure <= 0.
121    */
122   public Barostat createBarostat(MolecularAssembly assembly, CrystalPotential crystalPotential)
123       throws IllegalArgumentException {
124     if (group.pressure > 0) {
125       Barostat barostat = new Barostat(assembly, crystalPotential);
126       barostat.setPressure(group.pressure);
127       barostat.setIsotropic(group.isotropic);
128       barostat.setMaxDensity(group.maxD);
129       barostat.setMinDensity(group.minD);
130       barostat.setBarostatPrintFrequency(group.barPrint);
131       double dens = barostat.density();
132       if (dens < group.minD) {
133         logger.info(format(
134             " Barostat: initial density %9.4g < minimum density %9.4g, resetting to minimum density",
135             dens, group.minD));
136         barostat.setDensity(group.minD);
137       } else if (dens > group.maxD) {
138         logger.info(format(
139             " Barostat: initial density %9.4g > maximum density %9.4g, resetting to maximum density",
140             dens, group.maxD));
141         barostat.setDensity(group.maxD);
142       }
143       barostat.setMaxAngleMove(group.maxAM);
144       barostat.setMaxVolumeMove(group.maxV);
145       barostat.setMeanBarostatInterval(group.barInt);
146       return barostat;
147     } else {
148       throw new IllegalArgumentException(" Pressure is <= 0; cannot create a Barostat!");
149     }
150   }
151 
152   /**
153    * -p or --npt Specify use of a Monte Carlo Barostat at the given pressure (default 0 = constant
154    * volume).
155    *
156    * @return Returns the pressure (atm).
157    */
158   public double getPressure() {
159     return group.pressure;
160   }
161 
162   public void setPressure(double pressure) {
163     group.pressure = pressure;
164   }
165 
166   /**
167    * Restrict the MC Barostat to isotropic moves. The lattice angles are * held fixed, and lattice
168    * lengths are scaled equally.
169    *
170    * @return Returns true if the Barostat is isotropic.
171    */
172   public boolean isIsotropic() {
173     return group.isotropic;
174   }
175 
176   public void setIsotropic(boolean isotropic) {
177     group.isotropic = isotropic;
178   }
179 
180   /**
181    * The maximum density accepted by the MC Barostat (g/cc).
182    *
183    * @return Returns the maximum density.
184    */
185   public double getMaxD() {
186     return group.maxD;
187   }
188 
189   public void setMaxD(double maxD) {
190     group.maxD = maxD;
191   }
192 
193   /**
194    * The minimum density accepted by the MC Barostat (g/cc).
195    *
196    * @return Returns the minimum density.
197    */
198   public double getMinD() {
199     return group.minD;
200   }
201 
202   public void setMinD(double minD) {
203     group.minD = minD;
204   }
205 
206   /**
207    * The volume of a proposed unit cell volume moves (uniformly distributed) in Angstroms^3.
208    *
209    * @return Returns the maximum volume move.
210    */
211   public double getMaxV() {
212     return group.maxV;
213   }
214 
215   /**
216    * The width of proposed crystal angle moves (uniformly distributed) in degrees.
217    *
218    * @return Returns the width of angle moves.
219    */
220   public double getMaxAM() {
221     return group.maxAM;
222   }
223 
224   public void setMaxAM(double maxAM) {
225     group.maxAM = maxAM;
226   }
227 
228   public void setMaxV(double maxV) {
229     group.maxV = maxV;
230   }
231 
232   /**
233    * The mean number of MD steps (Poisson distribution) between barostat move proposals.
234    *
235    * @return Returns the mean number of MD steps between barostat MC trials.
236    */
237   public int getBarInt() {
238     return group.barInt;
239   }
240 
241   public void setBarInt(int barInt) {
242     group.barInt = barInt;
243   }
244 
245   /**
246    * --bpi or --barostatPrintInterval Sets the number of Barostat MC cycles between print
247    * statements.
248    */
249   public int getPrintInt() {
250     return group.barPrint;
251   }
252 
253   public void setPrintInt(int printInterval) {
254     group.barPrint = printInterval;
255   }
256 
257   /**
258    * Collection of Barostat Options.
259    */
260   private static class BarostatOptionGroup {
261 
262     /**
263      * -p or --npt Specify use of a Monte Carlo Barostat at the given pressure (default 0 = constant
264      * volume).
265      */
266     @Option(names = {"-p", "--npt"}, paramLabel = "0", defaultValue = "0",
267         description = "Specify use of a MC Barostat at the given pressure; the default 0 disables NPT (atm).")
268     private double pressure = 0;
269     /**
270      * -iso or --isotropic Restrict the MC Barostat to isotropic moves. The lattice angles are held
271      * fixed, and lattice lengths are scaled equally.
272      */
273     @Option(names = {"--iso", "--isotropic"},
274         defaultValue = "false", description = "Restrict the MC Barostat to isotropic moves.")
275     private boolean isotropic = false;
276     /**
277      * --maxD or --maxDensity Specify the maximum density accepted by the MC Barostat (g/cc).
278      */
279     @Option(names = {"--maxD", " --maxDensity"},
280         paramLabel = DEFAULT_MAX_DENSITY, defaultValue = DEFAULT_MAX_DENSITY,
281         description = "Specify the maximum density accepted by the MC Barostat (g/cc).")
282     private double maxD = parseDouble(DEFAULT_MAX_DENSITY);
283     /**
284      * --minD or --minDensity Specify the minimum density accepted by the MC Barostat (g/cc).
285      */
286     @Option(names = {"--minD", " --minDensity"},
287         paramLabel = DEFAULT_MIN_DENSITY, defaultValue = DEFAULT_MIN_DENSITY,
288         description = "Specify the minimum density accepted by the MC Barostat (g/cc).")
289     private double minD = parseDouble(DEFAULT_MIN_DENSITY);
290     @Option(names = {"--maxAM", "--maxAngleMove"},
291         paramLabel = DEFAULT_MAX_ANGLE_MOVE, defaultValue = DEFAULT_MAX_ANGLE_MOVE,
292         description = "Sets the width of proposed crystal angle moves (uniformly distributed) in degrees.")
293     private double maxAM = parseDouble(DEFAULT_MAX_ANGLE_MOVE);
294     /**
295      * --maxV or --maxVolumeMove Sets the volume move size (uniformly distributed) in Angstroms^3.
296      */
297     @Option(names = {"--maxV", "--maxVolumeMove"},
298         paramLabel = DEFAULT_MAX_VOLUME_MOVE, defaultValue = DEFAULT_MAX_VOLUME_MOVE,
299         description = "Default width of proposed unit cell side length moves (uniformly distributed) in Angstroms.")
300     private double maxV = parseDouble(DEFAULT_MAX_VOLUME_MOVE);
301 
302     /**
303      * --barInt or --meanBarostatInterval Sets the mean number of MD steps (Poisson distribution)
304      * between barostat move proposals.
305      */
306     @Option(names = {"--barInt", "--meanBarostatInterval"},
307         paramLabel = DEFAULT_BAROSTAT_INTERVAL, defaultValue = DEFAULT_BAROSTAT_INTERVAL,
308         description = "Sets the mean number of MD steps between barostat move proposals.")
309     private int barInt = parseInt(DEFAULT_BAROSTAT_INTERVAL);
310     /**
311      * --bpi or --barostatPrintInterval Sets the number of Barostat MC cycles between print
312      * statements.
313      */
314     @Option(names = {"--bpi", "--barostatPrintInterval"},
315         paramLabel = DEFAULT_BAROSTAT_PRINT_INTERVAL, defaultValue = DEFAULT_BAROSTAT_PRINT_INTERVAL,
316         description = "Sets the number of Barostat MC cycles between print statements.")
317     private int barPrint = parseInt(DEFAULT_BAROSTAT_PRINT_INTERVAL);
318   }
319 }