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,2018,2019, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
37 #ifndef GMX_MDTYPES_INPUTREC_H
38 #define GMX_MDTYPES_INPUTREC_H
44 #include "gromacs/math/vectypes.h"
45 #include "gromacs/mdtypes/md_enums.h"
46 #include "gromacs/utility/basedefinitions.h"
47 #include "gromacs/utility/real.h"
49 #define EGP_EXCL (1<<0)
50 #define EGP_TABLE (1<<1)
56 typedef struct t_swap *gmx_swapcoords_t;
62 class KeyValueTreeObject;
67 //! Number of T-Coupl groups
69 //! Number of of Nose-Hoover chains per group
71 //! Number of Accelerate groups
73 //! Number of Freeze groups
75 //! Number of Energy groups
77 //! Number of degrees of freedom in a group
79 //! Coupling temperature per group
81 //! No/simple/periodic simulated annealing for each group
83 //! Number of annealing time points per group
85 //! For each group: Time points
87 //! For each group: Temperature at these times. Final temp after all intervals is ref_t
91 //! Acceleration per group
93 //! Whether the group will be frozen in each direction
95 //! Exclusions/tables of energy group pairs
99 //! Number of QM groups
101 //! Level of theory in the QM calculation
103 //! Basisset in the QM calculation
105 //! Total charge in the QM region
107 //! Spin multiplicicty in the QM region
109 //! Surface hopping (diabatic hop only)
111 //! Number of orbiatls in the active space
113 //! Number of electrons in the active space
115 //! At which gap (A.U.) the SA is switched on
117 //! At which gap (A.U.) the SA is switched off
119 //! In how many steps SA goes from 1-1 to 0.5-0.5
125 //! Simulated temperature scaling; linear or exponential
127 //! The low temperature for simulated tempering
129 //! The high temperature for simulated tempering
131 //! The range of temperatures used for simulated tempering
137 //! The frequency for calculating dhdl
139 //! 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)
141 //! The initial number of the state
143 //! Change of lambda per time step (fraction of (0.1)
145 //! Print no, total or potential energies in dhdl
146 int edHdLPrintEnergy;
147 //! The number of foreign lambda points
149 //! The array of all lambda values
151 //! The number of neighboring lambda states to calculate the energy for in up and down directions (-1 for all)
152 int lambda_neighbors;
153 //! The first lambda to calculate energies for
155 //! The last lambda +1 to calculate energies for
157 //! Free energy soft-core parameter
159 //! Lambda power for soft-core interactions
161 //! R power for soft-core interactions
163 //! Free energy soft-core sigma when c6 or c12=0
165 //! Free energy soft-core sigma for ?????
167 //! Use softcore for the coulomb portion as well (default FALSE)
169 //! Whether to print the dvdl term associated with this term; if it is not specified as separate, it is lumped with the FEP term
170 gmx_bool separate_dvdl[efptNR];
171 //! Whether to write a separate dhdl.xvg file note: NOT a gmx_bool, but an enum
172 int separate_dhdl_file;
173 //! Whether to calculate+write dhdl derivatives note: NOT a gmx_bool, but an enum
174 int dhdl_derivatives;
175 //! The maximum table size for the dH histogram
177 //! The spacing for the dH histogram
178 double dh_hist_spacing;
182 //! The frequency of expanded ensemble state changes
184 //! Which type of move updating do we use for lambda monte carlo (or no for none)
186 //! What move set will be we using for state space moves
188 //! The method we use to decide of we have equilibrated the weights
190 //! The minumum number of samples at each lambda for deciding whether we have reached a minimum
192 //! Wang-Landau delta at which we stop equilibrating weights
194 //! Use the ratio of weights (ratio of minimum to maximum) to decide when to stop equilibrating
196 //! After equil_steps steps we stop equilibrating the weights
198 //! After equil_samples total samples (steps/nstfep), we stop equilibrating the weights
200 //! Random number seed for lambda mc switches
202 //! Whether to use minumum variance weighting
204 //! The number of samples needed before kicking into minvar routine
206 //! The offset for the variance in MinVar
208 //! Range of cvalues used for BAR
210 //! Whether to print symmetrized matrices
211 gmx_bool bSymmetrizedTMatrix;
212 //! How frequently to print the transition matrices
214 //! Number of repetitions in the MC lambda jumps MRS -- VERIFY THIS
216 //! Minimum number of samples for each state before free sampling MRS -- VERIFY THIS!
217 int lmc_forced_nstart;
218 //! Distance in lambda space for the gibbs interval
220 //! Scaling factor for Wang-Landau
222 //! Ratio between largest and smallest number for freezing the weights
224 //! Starting delta for Wang-Landau
226 //! Use one over t convergence for Wang-Landau when the delta get sufficiently small
227 gmx_bool bWLoneovert;
228 //! Did we initialize the weights? TODO: REMOVE FOR 5.0, no longer needed with new logic
229 gmx_bool bInit_weights;
230 //! To override the main temperature, or define it if it's not defined
232 //! User-specified initial weights to start with
233 real *init_lambda_weights;
238 //! Rotation type for this group
240 //! Use mass-weighed positions?
242 //! Number of atoms in the group
244 //! The global atoms numbers
246 //! The reference positions
248 //! The normalized rotation vector
250 //! Rate of rotation (degree/ps)
252 //! Force constant (kJ/(mol nm^2)
254 //! Pivot point of rotation axis (nm)
256 //! Type of fit to determine actual group angle
258 //! Number of angles around the reference angle for which the rotation potential is also evaluated (for fit type 'potential' only)
260 //! Distance between two angles in degrees (for fit type 'potential' only)
262 //! Slab distance (nm)
264 //! Minimum value the gaussian must have so that the force is actually evaluated
266 //! Additive constant for radial motion2 and flexible2 potentials (nm^2)
272 //! Number of rotation groups
274 //! Output frequency for main rotation outfile
276 //! Output frequency for per-slab data
284 //! Number of interactive atoms
286 //! The global indices of the interactive atoms
292 //! Name of the swap group, e.g. NA, CL, SOL
294 //! Number of atoms in this group
296 //! The global ion group atoms numbers
298 //! Requested number of molecules of this type per compartment
299 int nmolReq[eCompNR];
304 //! Period between when a swap is attempted
306 //! Use mass-weighted positions in split group
307 gmx_bool massw_split[2];
308 /*! \brief Split cylinders defined by radius, upper and lower
309 * extension. The split cylinders define the channels and are
310 * each anchored in the center of the split group */
316 //! Coupling constant (number of swap attempt steps)
318 //! Ion counts may deviate from the requested values by +-threshold before a swap is done
320 //! Offset of the swap layer (='bulk') with respect to the compartment-defining layers
321 real bulkOffset[eCompNR];
322 //! Number of groups to be controlled
324 //! All swap groups, including split and solvent
326 //! Swap private data accessible in swapcoords.cpp
327 gmx_swapcoords_t si_priv;
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 //! Which neighbor search method should we use?
351 //! Number of steps before pairlist is generated
353 //! Number of cells per rlong
355 //! Number of steps after which center of mass motion is removed
357 //! Center of mass motion removal algorithm
359 //! Number of steps after which print to logfile
361 //! Number of steps after which X is output
363 //! Number of steps after which V is output
365 //! Number of steps after which F is output
367 //! Number of steps after which energies printed
369 //! Number of steps after which compressed trj (.xtc,.tng) is output
370 int nstxout_compressed;
371 //! Initial time (ps)
375 //! Precision of x in compressed trajectory file
376 real x_compression_precision;
377 //! Requested fourier_spacing, when nk? not set
378 real fourier_spacing;
379 //! Number of k vectors in x dimension for fourier methods for long range electrost.
381 //! Number of k vectors in y dimension for fourier methods for long range electrost.
383 //! Number of k vectors in z dimension for fourier methods for long range electrost.
385 //! Interpolation order for PME
387 //! Real space tolerance for Ewald, determines the real/reciprocal space relative weight
389 //! Real space tolerance for LJ-Ewald
391 //! Normal/3D ewald, or pseudo-2D LR corrections
393 //! Epsilon for PME dipole correction
394 real epsilon_surface;
395 //! Type of combination rule in LJ-PME
396 int ljpme_combination_rule;
397 //! Type of periodic boundary conditions
399 //! Periodic molecules
401 //! Continuation run: starting state is correct (ie. constrained)
402 gmx_bool bContinuation;
403 //! Temperature coupling
405 //! Interval in steps for temperature coupling
407 //! Whether to print nose-hoover chains
408 gmx_bool bPrintNHChains;
409 //! Pressure coupling
411 //! Pressure coupling type
413 //! Interval in steps for pressure coupling
415 //! Pressure coupling time (ps)
417 //! Reference pressure (kJ/(mol nm^3))
419 //! Compressibility ((mol nm^3)/kJ)
421 //! How to scale absolute reference coordinates
422 int refcoord_scaling;
423 //! The COM of the posres atoms
425 //! The B-state COM of the posres atoms
427 //! Random seed for Andersen thermostat (obsolete)
429 //! Per atom pair energy drift tolerance (kJ/mol/ps/atom) for list buffer
431 //! Short range pairlist cut-off (nm)
433 //! Radius for test particle insertion
435 //! Type of electrostatics treatment
437 //! Modify the Coulomb interaction
438 int coulomb_modifier;
439 //! Coulomb switch range start (nm)
440 real rcoulomb_switch;
441 //! Coulomb cutoff (nm)
443 //! Relative dielectric constant
445 //! Relative dielectric constant of the RF
447 //! Always false (no longer supported)
448 bool implicit_solvent;
449 //! Type of Van der Waals treatment
451 //! Modify the Van der Waals interaction
453 //! Van der Waals switch range start (nm)
455 //! Van der Waals cutoff (nm)
457 //! Perform Long range dispersion corrections
459 //! Extension of the table beyond the cut-off, as well as the table length for 1-4 interac.
461 //! Tolerance for shake
463 //! Free energy calculations
465 //! Data for the FEP state
467 //! Whether to do simulated tempering
469 //! Variables for simulated tempering
470 t_simtemp *simtempvals;
471 //! Whether expanded ensembles are used
473 //! Expanded ensemble parameters
474 t_expanded *expandedvals;
475 //! Type of distance restraining
477 //! Force constant for time averaged distance restraints
479 //! Type of weighting of pairs in one restraints
481 //! Use combination of time averaged and instantaneous violations
482 gmx_bool bDisreMixed;
483 //! Frequency of writing pair distances to enx
485 //! Time constant for memory function in disres
487 //! Force constant for orientational restraints
489 //! Time constant for memory function in orires
491 //! Frequency of writing tr(SD) to energy output
493 //! The stepsize for updating
497 //! Number of iterations for convergence of steepest descent in relax_shells
499 //! Stepsize for directional minimization in relax_shells
501 //! Number of steps after which a steepest descents step is done while doing cg
503 //! Number of corrections to the Hessian to keep
505 //! Type of constraint algorithm
507 //! Order of the LINCS Projection Algorithm
509 //! Warn if any bond rotates more than this many degrees
511 //! Number of iterations in the final LINCS step
513 //! Use successive overrelaxation for shake
515 //! Friction coefficient for BD (amu/ps)
517 //! Random seed for SD and BD
519 //! The number of walls
521 //! The type of walls
523 //! The potentail is linear for r<=wall_r_linpot
525 //! The atom type for walls
526 int wall_atomtype[2];
527 //! Number density for walls
528 real wall_density[2];
529 //! Scaling factor for the box for Ewald
530 real wall_ewald_zfac;
532 /* COM pulling data */
533 //! Do we do COM pulling?
535 //! The data for center of mass pulling
537 // TODO: Remove this by converting pull into a ForceProvider
538 //! The COM pull force calculation data structure
542 //! Whether to use AWH biasing for PMF calculations
544 //! AWH biasing parameters
545 gmx::AwhParams *awhParams;
547 /* Enforced rotation data */
548 //! Whether to calculate enforced rotation potential(s)
550 //! The data for enforced rotation potentials
553 //! Whether to do ion/water position exchanges (CompEL)
555 //! Swap data structure.
558 //! Whether the tpr makes an interactive MD session possible.
560 //! Interactive molecular dynamics
563 //! Acceleration for viscosity calculation
565 //! Triclinic deformation velocities (nm/ps)
567 /*! \brief User determined parameters */
580 //! QM/MM calculation
582 //! Constraints on QM bonds
584 //! Scheme: ONIOM or normal
586 //! Factor for scaling the MM charges in QM calc.
589 /* Fields for removed features go here (better caching) */
590 //! Whether AdResS is enabled - always false if a valid .tpr was read
592 //! Whether twin-range scheme is active - always false if a valid .tpr was read
593 gmx_bool useTwinRange;
595 //! KVT object that contains input parameters converted to the new style.
596 gmx::KeyValueTreeObject *params;
599 int ir_optimal_nstcalcenergy(const t_inputrec *ir);
601 int tcouple_min_integration_steps(int etc);
603 int ir_optimal_nsttcouple(const t_inputrec *ir);
605 int pcouple_min_integration_steps(int epc);
607 int ir_optimal_nstpcouple(const t_inputrec *ir);
609 /* Returns if the Coulomb force or potential is switched to zero */
610 gmx_bool ir_coulomb_switched(const t_inputrec *ir);
612 /* Returns if the Coulomb interactions are zero beyond the rcoulomb.
613 * Note: always returns TRUE for the Verlet cut-off scheme.
615 gmx_bool ir_coulomb_is_zero_at_cutoff(const t_inputrec *ir);
617 /* As ir_coulomb_is_zero_at_cutoff, but also returns TRUE for user tabulated
618 * interactions, since these might be zero beyond rcoulomb.
620 gmx_bool ir_coulomb_might_be_zero_at_cutoff(const t_inputrec *ir);
622 /* Returns if the Van der Waals force or potential is switched to zero */
623 gmx_bool ir_vdw_switched(const t_inputrec *ir);
625 /* Returns if the Van der Waals interactions are zero beyond the rvdw.
626 * Note: always returns TRUE for the Verlet cut-off scheme.
628 gmx_bool ir_vdw_is_zero_at_cutoff(const t_inputrec *ir);
630 /* As ir_vdw_is_zero_at_cutoff, but also returns TRUE for user tabulated
631 * interactions, since these might be zero beyond rvdw.
633 gmx_bool ir_vdw_might_be_zero_at_cutoff(const t_inputrec *ir);
635 /*! \brief Free memory from input record.
637 * All arrays and pointers will be freed.
639 * \param[in] ir The data structure
641 void done_inputrec(t_inputrec *ir);
643 void pr_inputrec(FILE *fp, int indent, const char *title, const t_inputrec *ir,
644 gmx_bool bMDPformat);
646 void cmp_inputrec(FILE *fp, const t_inputrec *ir1, const t_inputrec *ir2, real ftol, real abstol);
648 void comp_pull_AB(FILE *fp, pull_params_t *pull, real ftol, real abstol);
651 gmx_bool inputrecDeform(const t_inputrec *ir);
653 gmx_bool inputrecDynamicBox(const t_inputrec *ir);
655 gmx_bool inputrecPreserveShape(const t_inputrec *ir);
657 gmx_bool inputrecNeedMutot(const t_inputrec *ir);
659 gmx_bool inputrecTwinRange(const t_inputrec *ir);
661 gmx_bool inputrecExclForces(const t_inputrec *ir);
663 gmx_bool inputrecNptTrotter(const t_inputrec *ir);
665 gmx_bool inputrecNvtTrotter(const t_inputrec *ir);
667 gmx_bool inputrecNphTrotter(const t_inputrec *ir);
669 /*! \brief Return true if the simulation is 2D periodic with two walls. */
670 bool inputrecPbcXY2Walls(const t_inputrec *ir);
672 /*! \brief Returns true for MD integator with T and/or P-coupling that supports
673 * calculating the conserved energy quantity.
675 bool integratorHasConservedEnergyQuantity(const t_inputrec *ir);
677 /*! \brief Returns true when temperature is coupled or constant. */
678 bool integratorHasReferenceTemperature(const t_inputrec *ir);
680 /*! \brief Return the number of bounded dimensions
682 * \param[in] ir The input record with MD parameters
683 * \return the number of dimensions in which
684 * the coordinates of the particles are bounded, starting at X.
686 int inputrec2nboundeddim(const t_inputrec *ir);
688 /*! \brief Returns the number of degrees of freedom in center of mass motion
690 * \param[in] ir The inputrec structure
691 * \return the number of degrees of freedom of the center of mass
693 int ndof_com(const t_inputrec *ir);
695 /*! \brief Returns the maximum reference temperature over all coupled groups
697 * Returns 0 for energy minimization and normal mode computation.
698 * Returns -1 for MD without temperature coupling.
700 * \param[in] ir The inputrec structure
702 real maxReferenceTemperature(const t_inputrec &ir);
704 #endif /* GMX_MDTYPES_INPUTREC_H */