2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
7 * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
8 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
9 * and including many others, as listed in the AUTHORS file in the
10 * top-level source directory and at http://www.gromacs.org.
12 * GROMACS is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public License
14 * as published by the Free Software Foundation; either version 2.1
15 * of the License, or (at your option) any later version.
17 * GROMACS is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 * Lesser General Public License for more details.
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with GROMACS; if not, see
24 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
25 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
27 * If you want to redistribute modifications to GROMACS, please
28 * consider that scientific software is very special. Version
29 * control is crucial - bugs must be traceable. We will be happy to
30 * consider code for inclusion in the official distribution, but
31 * derived work must not be called official GROMACS. Details are found
32 * in the README & COPYING files - if they are missing, get the
33 * official version at http://www.gromacs.org.
35 * To help us fund GROMACS development, we humbly ask that you cite
36 * the research papers on the package. Check out http://www.gromacs.org.
38 #ifndef GMX_MDTYPES_INPUTREC_H
39 #define GMX_MDTYPES_INPUTREC_H
45 #include "gromacs/math/vectypes.h"
46 #include "gromacs/mdtypes/md_enums.h"
47 #include "gromacs/utility/basedefinitions.h"
48 #include "gromacs/utility/real.h"
50 #define EGP_EXCL (1 << 0)
51 #define EGP_TABLE (1 << 1)
61 class KeyValueTreeObject
;
68 //! Number of T-Coupl groups
70 //! Number of of Nose-Hoover chains per group
72 //! Number of Accelerate groups
74 //! Number of Freeze groups
76 //! Number of Energy groups
78 //! Number of degrees of freedom in a group
80 //! Coupling temperature per group
82 //! No/simple/periodic simulated annealing for each group
84 //! Number of annealing time points per group
86 //! For each group: Time points
88 //! For each group: Temperature at these times. Final temp after all intervals is ref_t
92 //! Acceleration per group
94 //! Whether the group will be frozen in each direction
96 //! Exclusions/tables of energy group pairs
100 //! Number of QM groups
102 //! Level of theory in the QM calculation
104 //! Basisset in the QM calculation
106 //! Total charge in the QM region
108 //! Spin multiplicicty in the QM region
110 //! Surface hopping (diabatic hop only)
112 //! Number of orbiatls in the active space
114 //! Number of electrons in the active space
116 //! At which gap (A.U.) the SA is switched on
118 //! At which gap (A.U.) the SA is switched off
120 //! In how many steps SA goes from 1-1 to 0.5-0.5
126 //! Simulated temperature scaling; linear or exponential
128 //! The low temperature for simulated tempering
130 //! The high temperature for simulated tempering
132 //! The range of temperatures used for simulated tempering
138 //! The frequency for calculating dhdl
140 //! Fractional value of lambda (usually will use init_fep_state, this will only be for slow growth, and for legacy free energy code. Only has a valid value if positive)
142 //! The initial number of the state
144 //! Change of lambda per time step (fraction of (0.1)
146 //! Print no, total or potential energies in dhdl
147 int edHdLPrintEnergy
;
148 //! The number of foreign lambda points
150 //! The array of all lambda values
152 //! The number of neighboring lambda states to calculate the energy for in up and down directions (-1 for all)
153 int lambda_neighbors
;
154 //! The first lambda to calculate energies for
156 //! The last lambda +1 to calculate energies for
158 //! Free energy soft-core parameter
160 //! Lambda power for soft-core interactions
162 //! R power for soft-core interactions
164 //! Free energy soft-core sigma when c6 or c12=0
166 //! Free energy soft-core sigma for ?????
168 //! Use softcore for the coulomb portion as well (default FALSE)
170 //! Whether to print the dvdl term associated with this term; if it is not specified as separate, it is lumped with the FEP term
171 gmx_bool separate_dvdl
[efptNR
];
172 //! Whether to write a separate dhdl.xvg file note: NOT a gmx_bool, but an enum
173 int separate_dhdl_file
;
174 //! Whether to calculate+write dhdl derivatives note: NOT a gmx_bool, but an enum
175 int dhdl_derivatives
;
176 //! The maximum table size for the dH histogram
178 //! The spacing for the dH histogram
179 double dh_hist_spacing
;
184 //! The frequency of expanded ensemble state changes
186 //! Which type of move updating do we use for lambda monte carlo (or no for none)
188 //! What move set will be we using for state space moves
190 //! The method we use to decide of we have equilibrated the weights
192 //! The minumum number of samples at each lambda for deciding whether we have reached a minimum
194 //! Wang-Landau delta at which we stop equilibrating weights
196 //! Use the ratio of weights (ratio of minimum to maximum) to decide when to stop equilibrating
198 //! After equil_steps steps we stop equilibrating the weights
200 //! After equil_samples total samples (steps/nstfep), we stop equilibrating the weights
202 //! Random number seed for lambda mc switches
204 //! Whether to use minumum variance weighting
206 //! The number of samples needed before kicking into minvar routine
208 //! The offset for the variance in MinVar
210 //! Range of cvalues used for BAR
212 //! Whether to print symmetrized matrices
213 gmx_bool bSymmetrizedTMatrix
;
214 //! How frequently to print the transition matrices
216 //! Number of repetitions in the MC lambda jumps MRS -- VERIFY THIS
218 //! Minimum number of samples for each state before free sampling MRS -- VERIFY THIS!
219 int lmc_forced_nstart
;
220 //! Distance in lambda space for the gibbs interval
222 //! Scaling factor for Wang-Landau
224 //! Ratio between largest and smallest number for freezing the weights
226 //! Starting delta for Wang-Landau
228 //! Use one over t convergence for Wang-Landau when the delta get sufficiently small
229 gmx_bool bWLoneovert
;
230 //! Did we initialize the weights? TODO: REMOVE FOR 5.0, no longer needed with new logic
231 gmx_bool bInit_weights
;
232 //! To override the main temperature, or define it if it's not defined
234 //! User-specified initial weights to start with
235 real
* init_lambda_weights
;
240 //! Rotation type for this group
242 //! Use mass-weighed positions?
244 //! Number of atoms in the group
246 //! The global atoms numbers
248 //! The reference positions
250 //! The normalized rotation vector
252 //! Rate of rotation (degree/ps)
254 //! Force constant (kJ/(mol nm^2)
256 //! Pivot point of rotation axis (nm)
258 //! Type of fit to determine actual group angle
260 //! Number of angles around the reference angle for which the rotation potential is also evaluated (for fit type 'potential' only)
262 //! Distance between two angles in degrees (for fit type 'potential' only)
264 //! Slab distance (nm)
266 //! Minimum value the gaussian must have so that the force is actually evaluated
268 //! Additive constant for radial motion2 and flexible2 potentials (nm^2)
274 //! Number of rotation groups
276 //! Output frequency for main rotation outfile
278 //! Output frequency for per-slab data
286 //! Number of interactive atoms
288 //! The global indices of the interactive atoms
294 //! Name of the swap group, e.g. NA, CL, SOL
296 //! Number of atoms in this group
298 //! The global ion group atoms numbers
300 //! Requested number of molecules of this type per compartment
301 int nmolReq
[eCompNR
];
306 //! Period between when a swap is attempted
308 //! Use mass-weighted positions in split group
309 gmx_bool massw_split
[2];
310 /*! \brief Split cylinders defined by radius, upper and lower
311 * extension. The split cylinders define the channels and are
312 * each anchored in the center of the split group */
318 //! Coupling constant (number of swap attempt steps)
320 //! Ion counts may deviate from the requested values by +-threshold before a swap is done
322 //! Offset of the swap layer (='bulk') with respect to the compartment-defining layers
323 real bulkOffset
[eCompNR
];
324 //! Number of groups to be controlled
326 //! All swap groups, including split and solvent
330 struct t_inputrec
// NOLINT (clang-analyzer-optin.performance.Padding)
333 explicit t_inputrec(const t_inputrec
&) = delete;
334 t_inputrec
& operator=(const t_inputrec
&) = delete;
337 //! Integration method
339 //! Number of steps to be taken
341 //! Used in checkpointing to separate chunks
343 //! Start at a stepcount >0 (used w. convert-tpr)
345 //! Frequency of energy calc. and T/P coupl. upd.
347 //! Group or verlet cutoffs
349 //! Number of steps before pairlist is generated
351 //! Number of steps after which center of mass motion is removed
353 //! Center of mass motion removal algorithm
355 //! Number of steps after which print to logfile
357 //! Number of steps after which X is output
359 //! Number of steps after which V is output
361 //! Number of steps after which F is output
363 //! Number of steps after which energies printed
365 //! Number of steps after which compressed trj (.xtc,.tng) is output
366 int nstxout_compressed
;
367 //! Initial time (ps)
371 //! Precision of x in compressed trajectory file
372 real x_compression_precision
;
373 //! Requested fourier_spacing, when nk? not set
374 real fourier_spacing
;
375 //! Number of k vectors in x dimension for fourier methods for long range electrost.
377 //! Number of k vectors in y dimension for fourier methods for long range electrost.
379 //! Number of k vectors in z dimension for fourier methods for long range electrost.
381 //! Interpolation order for PME
383 //! Real space tolerance for Ewald, determines the real/reciprocal space relative weight
385 //! Real space tolerance for LJ-Ewald
387 //! Normal/3D ewald, or pseudo-2D LR corrections
389 //! Epsilon for PME dipole correction
390 real epsilon_surface
;
391 //! Type of combination rule in LJ-PME
392 int ljpme_combination_rule
;
393 //! Type of periodic boundary conditions
395 //! Periodic molecules
397 //! Continuation run: starting state is correct (ie. constrained)
398 gmx_bool bContinuation
;
399 //! Temperature coupling
401 //! Interval in steps for temperature coupling
403 //! Whether to print nose-hoover chains
404 gmx_bool bPrintNHChains
;
405 //! Pressure coupling
407 //! Pressure coupling type
409 //! Interval in steps for pressure coupling
411 //! Pressure coupling time (ps)
413 //! Reference pressure (kJ/(mol nm^3))
415 //! Compressibility ((mol nm^3)/kJ)
417 //! How to scale absolute reference coordinates
418 int refcoord_scaling
;
419 //! The COM of the posres atoms
421 //! The B-state COM of the posres atoms
423 //! Random seed for Andersen thermostat (obsolete)
425 //! Per atom pair energy drift tolerance (kJ/mol/ps/atom) for list buffer
427 //! Short range pairlist cut-off (nm)
429 //! Radius for test particle insertion
431 //! Type of electrostatics treatment
433 //! Modify the Coulomb interaction
434 int coulomb_modifier
;
435 //! Coulomb switch range start (nm)
436 real rcoulomb_switch
;
437 //! Coulomb cutoff (nm)
439 //! Relative dielectric constant
441 //! Relative dielectric constant of the RF
443 //! Always false (no longer supported)
444 bool implicit_solvent
;
445 //! Type of Van der Waals treatment
447 //! Modify the Van der Waals interaction
449 //! Van der Waals switch range start (nm)
451 //! Van der Waals cutoff (nm)
453 //! Perform Long range dispersion corrections
455 //! Extension of the table beyond the cut-off, as well as the table length for 1-4 interac.
457 //! Tolerance for shake
459 //! Free energy calculations
461 //! Data for the FEP state
463 //! Whether to do simulated tempering
465 //! Variables for simulated tempering
466 t_simtemp
* simtempvals
;
467 //! Whether expanded ensembles are used
469 //! Expanded ensemble parameters
470 t_expanded
* expandedvals
;
471 //! Type of distance restraining
473 //! Force constant for time averaged distance restraints
475 //! Type of weighting of pairs in one restraints
477 //! Use combination of time averaged and instantaneous violations
478 gmx_bool bDisreMixed
;
479 //! Frequency of writing pair distances to enx
481 //! Time constant for memory function in disres
483 //! Force constant for orientational restraints
485 //! Time constant for memory function in orires
487 //! Frequency of writing tr(SD) to energy output
489 //! The stepsize for updating
493 //! Number of iterations for convergence of steepest descent in relax_shells
495 //! Stepsize for directional minimization in relax_shells
497 //! Number of steps after which a steepest descents step is done while doing cg
499 //! Number of corrections to the Hessian to keep
501 //! Type of constraint algorithm
503 //! Order of the LINCS Projection Algorithm
505 //! Warn if any bond rotates more than this many degrees
507 //! Number of iterations in the final LINCS step
509 //! Use successive overrelaxation for shake
511 //! Friction coefficient for BD (amu/ps)
513 //! Random seed for SD and BD
515 //! The number of walls
517 //! The type of walls
519 //! The potentail is linear for r<=wall_r_linpot
521 //! The atom type for walls
522 int wall_atomtype
[2];
523 //! Number density for walls
524 real wall_density
[2];
525 //! Scaling factor for the box for Ewald
526 real wall_ewald_zfac
;
528 /* COM pulling data */
529 //! Do we do COM pulling?
531 //! The data for center of mass pulling
535 //! Whether to use AWH biasing for PMF calculations
537 //! AWH biasing parameters
538 gmx::AwhParams
* awhParams
;
540 /* Enforced rotation data */
541 //! Whether to calculate enforced rotation potential(s)
543 //! The data for enforced rotation potentials
546 //! Whether to do ion/water position exchanges (CompEL)
548 //! Swap data structure.
551 //! Whether the tpr makes an interactive MD session possible.
553 //! Interactive molecular dynamics
556 //! Acceleration for viscosity calculation
558 //! Triclinic deformation velocities (nm/ps)
560 /*! \brief User determined parameters */
573 //! QM/MM calculation
575 //! Constraints on QM bonds
577 //! Scheme: ONIOM or normal
579 //! Factor for scaling the MM charges in QM calc.
582 /* Fields for removed features go here (better caching) */
583 //! Whether AdResS is enabled - always false if a valid .tpr was read
585 //! Whether twin-range scheme is active - always false if a valid .tpr was read
586 gmx_bool useTwinRange
;
588 //! KVT object that contains input parameters converted to the new style.
589 gmx::KeyValueTreeObject
* params
;
591 //! KVT for storing simulation parameters that are not part of the mdp file.
592 std::unique_ptr
<gmx::KeyValueTreeObject
> internalParameters
;
595 int ir_optimal_nstcalcenergy(const t_inputrec
* ir
);
597 int tcouple_min_integration_steps(int etc
);
599 int ir_optimal_nsttcouple(const t_inputrec
* ir
);
601 int pcouple_min_integration_steps(int epc
);
603 int ir_optimal_nstpcouple(const t_inputrec
* ir
);
605 /* Returns if the Coulomb force or potential is switched to zero */
606 gmx_bool
ir_coulomb_switched(const t_inputrec
* ir
);
608 /* Returns if the Coulomb interactions are zero beyond the rcoulomb.
609 * Note: always returns TRUE for the Verlet cut-off scheme.
611 gmx_bool
ir_coulomb_is_zero_at_cutoff(const t_inputrec
* ir
);
613 /* As ir_coulomb_is_zero_at_cutoff, but also returns TRUE for user tabulated
614 * interactions, since these might be zero beyond rcoulomb.
616 gmx_bool
ir_coulomb_might_be_zero_at_cutoff(const t_inputrec
* ir
);
618 /* Returns if the Van der Waals force or potential is switched to zero */
619 gmx_bool
ir_vdw_switched(const t_inputrec
* ir
);
621 /* Returns if the Van der Waals interactions are zero beyond the rvdw.
622 * Note: always returns TRUE for the Verlet cut-off scheme.
624 gmx_bool
ir_vdw_is_zero_at_cutoff(const t_inputrec
* ir
);
626 /* As ir_vdw_is_zero_at_cutoff, but also returns TRUE for user tabulated
627 * interactions, since these might be zero beyond rvdw.
629 gmx_bool
ir_vdw_might_be_zero_at_cutoff(const t_inputrec
* ir
);
631 /*! \brief Free memory from input record.
633 * All arrays and pointers will be freed.
635 * \param[in] ir The data structure
637 void done_inputrec(t_inputrec
* ir
);
639 void pr_inputrec(FILE* fp
, int indent
, const char* title
, const t_inputrec
* ir
, gmx_bool bMDPformat
);
641 void cmp_inputrec(FILE* fp
, const t_inputrec
* ir1
, const t_inputrec
* ir2
, real ftol
, real abstol
);
643 void comp_pull_AB(FILE* fp
, pull_params_t
* pull
, real ftol
, real abstol
);
646 gmx_bool
inputrecDeform(const t_inputrec
* ir
);
648 gmx_bool
inputrecDynamicBox(const t_inputrec
* ir
);
650 gmx_bool
inputrecPreserveShape(const t_inputrec
* ir
);
652 gmx_bool
inputrecNeedMutot(const t_inputrec
* ir
);
654 gmx_bool
inputrecTwinRange(const t_inputrec
* ir
);
656 gmx_bool
inputrecExclForces(const t_inputrec
* ir
);
658 gmx_bool
inputrecNptTrotter(const t_inputrec
* ir
);
660 gmx_bool
inputrecNvtTrotter(const t_inputrec
* ir
);
662 gmx_bool
inputrecNphTrotter(const t_inputrec
* ir
);
664 /*! \brief Return true if the simulation is 2D periodic with two walls. */
665 bool inputrecPbcXY2Walls(const t_inputrec
* ir
);
667 /*! \brief Returns true for MD integator with T and/or P-coupling that supports
668 * calculating the conserved energy quantity.
670 bool integratorHasConservedEnergyQuantity(const t_inputrec
* ir
);
672 /*! \brief Returns true when temperature is coupled or constant. */
673 bool integratorHasReferenceTemperature(const t_inputrec
* ir
);
675 /*! \brief Return the number of bounded dimensions
677 * \param[in] ir The input record with MD parameters
678 * \return the number of dimensions in which
679 * the coordinates of the particles are bounded, starting at X.
681 int inputrec2nboundeddim(const t_inputrec
* ir
);
683 /*! \brief Returns the number of degrees of freedom in center of mass motion
685 * \param[in] ir The inputrec structure
686 * \return the number of degrees of freedom of the center of mass
688 int ndof_com(const t_inputrec
* ir
);
690 /*! \brief Returns the maximum reference temperature over all coupled groups
692 * Returns 0 for energy minimization and normal mode computation.
693 * Returns -1 for MD without temperature coupling.
695 * \param[in] ir The inputrec structure
697 real
maxReferenceTemperature(const t_inputrec
& ir
);
699 /*! \brief Returns whether there is an Ewald surface contribution
701 bool haveEwaldSurfaceContribution(const t_inputrec
& ir
);
703 #endif /* GMX_MDTYPES_INPUTREC_H */