1 /* 2 * This file is part of the GROMACS molecular simulation package. 3 * 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,2021, 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. 11 * 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. 16 * 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. 21 * 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. 26 * 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. 34 * 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. 37 */ 38 #ifndef GMX_MDTYPES_INPUTREC_H 39 #define GMX_MDTYPES_INPUTREC_H 40 41 #include <cstdio> 42 43 #include <memory> 44 #include <vector> 45 46 #include "gromacs/math/vectypes.h" 47 #include "gromacs/mdtypes/md_enums.h" 48 #include "gromacs/utility/basedefinitions.h" 49 #include "gromacs/utility/real.h" 50 51 #define EGP_EXCL (1 << 0) 52 #define EGP_TABLE (1 << 1) 53 54 struct gmx_enfrot; 55 struct gmx_enfrotgrp; 56 struct pull_params_t; 57 58 namespace gmx 59 { 60 class Awh; 61 struct AwhParams; 62 class KeyValueTreeObject; 63 struct MtsLevel; 64 } // namespace gmx 65 66 enum class PbcType; 67 68 struct t_grpopts 69 { 70 //! Number of T-Coupl groups 71 int ngtc; 72 //! Number of of Nose-Hoover chains per group 73 int nhchainlength; 74 //! Number of Accelerate groups 75 int ngacc; 76 //! Number of Freeze groups 77 int ngfrz; 78 //! Number of Energy groups 79 int ngener; 80 //! Number of degrees of freedom in a group 81 real* nrdf; 82 //! Coupling temperature per group 83 real* ref_t; 84 //! No/simple/periodic simulated annealing for each group 85 int* annealing; 86 //! Number of annealing time points per group 87 int* anneal_npoints; 88 //! For each group: Time points 89 real** anneal_time; 90 //! For each group: Temperature at these times. Final temp after all intervals is ref_t 91 real** anneal_temp; 92 //! Tau coupling time 93 real* tau_t; 94 //! Acceleration per group 95 rvec* acc; 96 //! Whether the group will be frozen in each direction 97 ivec* nFreeze; 98 //! Exclusions/tables of energy group pairs 99 int* egp_flags; 100 101 /* QMMM stuff */ 102 //! Number of QM groups 103 int ngQM; 104 }; 105 106 struct t_simtemp 107 { 108 //! Simulated temperature scaling; linear or exponential 109 int eSimTempScale; 110 //! The low temperature for simulated tempering 111 real simtemp_low; 112 //! The high temperature for simulated tempering 113 real simtemp_high; 114 //! The range of temperatures used for simulated tempering 115 real* temperatures; 116 }; 117 118 struct t_lambda 119 { 120 //! The frequency for calculating dhdl 121 int nstdhdl; 122 //! 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) 123 double init_lambda; 124 //! The initial number of the state 125 int init_fep_state; 126 //! Change of lambda per time step (fraction of (0.1) 127 double delta_lambda; 128 //! Print no, total or potential energies in dhdl 129 int edHdLPrintEnergy; 130 //! The number of foreign lambda points 131 int n_lambda; 132 //! The array of all lambda values 133 double** all_lambda; 134 //! The number of neighboring lambda states to calculate the energy for in up and down directions (-1 for all) 135 int lambda_neighbors; 136 //! The first lambda to calculate energies for 137 int lambda_start_n; 138 //! The last lambda +1 to calculate energies for 139 int lambda_stop_n; 140 //! Free energy soft-core parameter 141 real sc_alpha; 142 //! Lambda power for soft-core interactions 143 int sc_power; 144 //! R power for soft-core interactions 145 real sc_r_power; 146 //! Free energy soft-core sigma when c6 or c12=0 147 real sc_sigma; 148 //! Free energy soft-core sigma for ????? 149 real sc_sigma_min; 150 //! Use softcore for the coulomb portion as well (default FALSE) 151 gmx_bool bScCoul; 152 //! Whether to print the dvdl term associated with this term; if it is not specified as separate, it is lumped with the FEP term 153 gmx_bool separate_dvdl[efptNR]; 154 //! Whether to write a separate dhdl.xvg file note: NOT a gmx_bool, but an enum 155 int separate_dhdl_file; 156 //! Whether to calculate+write dhdl derivatives note: NOT a gmx_bool, but an enum 157 int dhdl_derivatives; 158 //! The maximum table size for the dH histogram 159 int dh_hist_size; 160 //! The spacing for the dH histogram 161 double dh_hist_spacing; 162 }; 163 164 struct t_expanded 165 { 166 //! The frequency of expanded ensemble state changes 167 int nstexpanded; 168 //! Which type of move updating do we use for lambda monte carlo (or no for none) 169 int elamstats; 170 //! What move set will be we using for state space moves 171 int elmcmove; 172 //! The method we use to decide of we have equilibrated the weights 173 int elmceq; 174 //! The minumum number of samples at each lambda for deciding whether we have reached a minimum 175 int equil_n_at_lam; 176 //! Wang-Landau delta at which we stop equilibrating weights 177 real equil_wl_delta; 178 //! Use the ratio of weights (ratio of minimum to maximum) to decide when to stop equilibrating 179 real equil_ratio; 180 //! After equil_steps steps we stop equilibrating the weights 181 int equil_steps; 182 //! After equil_samples total samples (steps/nstfep), we stop equilibrating the weights 183 int equil_samples; 184 //! Random number seed for lambda mc switches 185 int lmc_seed; 186 //! Whether to use minumum variance weighting 187 gmx_bool minvar; 188 //! The number of samples needed before kicking into minvar routine 189 int minvarmin; 190 //! The offset for the variance in MinVar 191 real minvar_const; 192 //! Range of cvalues used for BAR 193 int c_range; 194 //! Whether to print symmetrized matrices 195 gmx_bool bSymmetrizedTMatrix; 196 //! How frequently to print the transition matrices 197 int nstTij; 198 //! Number of repetitions in the MC lambda jumps MRS -- VERIFY THIS 199 int lmc_repeats; 200 //! Minimum number of samples for each state before free sampling MRS -- VERIFY THIS! 201 int lmc_forced_nstart; 202 //! Distance in lambda space for the gibbs interval 203 int gibbsdeltalam; 204 //! Scaling factor for Wang-Landau 205 real wl_scale; 206 //! Ratio between largest and smallest number for freezing the weights 207 real wl_ratio; 208 //! Starting delta for Wang-Landau 209 real init_wl_delta; 210 //! Use one over t convergence for Wang-Landau when the delta get sufficiently small 211 gmx_bool bWLoneovert; 212 //! Did we initialize the weights? TODO: REMOVE FOR 5.0, no longer needed with new logic 213 gmx_bool bInit_weights; 214 //! To override the main temperature, or define it if it's not defined 215 real mc_temp; 216 //! User-specified initial weights to start with 217 real* init_lambda_weights; 218 }; 219 220 struct t_rotgrp 221 { 222 //! Rotation type for this group 223 int eType; 224 //! Use mass-weighed positions? 225 int bMassW; 226 //! Number of atoms in the group 227 int nat; 228 //! The global atoms numbers 229 int* ind; 230 //! The reference positions 231 rvec* x_ref; 232 //! The normalized rotation vector 233 rvec inputVec; 234 //! Rate of rotation (degree/ps) 235 real rate; 236 //! Force constant (kJ/(mol nm^2) 237 real k; 238 //! Pivot point of rotation axis (nm) 239 rvec pivot; 240 //! Type of fit to determine actual group angle 241 int eFittype; 242 //! Number of angles around the reference angle for which the rotation potential is also evaluated (for fit type 'potential' only) 243 int PotAngle_nstep; 244 //! Distance between two angles in degrees (for fit type 'potential' only) 245 real PotAngle_step; 246 //! Slab distance (nm) 247 real slab_dist; 248 //! Minimum value the gaussian must have so that the force is actually evaluated 249 real min_gaussian; 250 //! Additive constant for radial motion2 and flexible2 potentials (nm^2) 251 real eps; 252 }; 253 254 struct t_rot 255 { 256 //! Number of rotation groups 257 int ngrp; 258 //! Output frequency for main rotation outfile 259 int nstrout; 260 //! Output frequency for per-slab data 261 int nstsout; 262 //! Groups to rotate 263 t_rotgrp* grp; 264 }; 265 266 struct t_IMD 267 { 268 //! Number of interactive atoms 269 int nat; 270 //! The global indices of the interactive atoms 271 int* ind; 272 }; 273 274 struct t_swapGroup 275 { 276 //! Name of the swap group, e.g. NA, CL, SOL 277 char* molname; 278 //! Number of atoms in this group 279 int nat; 280 //! The global ion group atoms numbers 281 int* ind; 282 //! Requested number of molecules of this type per compartment 283 int nmolReq[eCompNR]; 284 }; 285 286 struct t_swapcoords 287 { 288 //! Period between when a swap is attempted 289 int nstswap; 290 //! Use mass-weighted positions in split group 291 gmx_bool massw_split[2]; 292 /*! \brief Split cylinders defined by radius, upper and lower 293 * extension. The split cylinders define the channels and are 294 * each anchored in the center of the split group */ 295 /**@{*/ 296 real cyl0r, cyl1r; 297 real cyl0u, cyl1u; 298 real cyl0l, cyl1l; 299 /**@}*/ 300 //! Coupling constant (number of swap attempt steps) 301 int nAverage; 302 //! Ion counts may deviate from the requested values by +-threshold before a swap is done 303 real threshold; 304 //! Offset of the swap layer (='bulk') with respect to the compartment-defining layers 305 real bulkOffset[eCompNR]; 306 //! Number of groups to be controlled 307 int ngrp; 308 //! All swap groups, including split and solvent 309 t_swapGroup* grp; 310 }; 311 312 struct t_inputrec // NOLINT (clang-analyzer-optin.performance.Padding) 313 { 314 t_inputrec(); 315 explicit t_inputrec(const t_inputrec&) = delete; 316 t_inputrec& operator=(const t_inputrec&) = delete; 317 ~t_inputrec(); 318 319 //! Integration method 320 int eI; 321 //! Number of steps to be taken 322 int64_t nsteps; 323 //! Used in checkpointing to separate chunks 324 int simulation_part; 325 //! Start at a stepcount >0 (used w. convert-tpr) 326 int64_t init_step; 327 //! Frequency of energy calc. and T/P coupl. upd. 328 int nstcalcenergy; 329 //! Group or verlet cutoffs 330 int cutoff_scheme; 331 //! Number of steps before pairlist is generated 332 int nstlist; 333 //! Number of steps after which center of mass motion is removed 334 int nstcomm; 335 //! Center of mass motion removal algorithm 336 int comm_mode; 337 //! Number of steps after which print to logfile 338 int nstlog; 339 //! Number of steps after which X is output 340 int nstxout; 341 //! Number of steps after which V is output 342 int nstvout; 343 //! Number of steps after which F is output 344 int nstfout; 345 //! Number of steps after which energies printed 346 int nstenergy; 347 //! Number of steps after which compressed trj (.xtc,.tng) is output 348 int nstxout_compressed; 349 //! Initial time (ps) 350 double init_t; 351 //! Time step (ps) 352 double delta_t; 353 //! Whether we use multiple time stepping 354 bool useMts; 355 //! The multiple time stepping levels 356 std::vector<gmx::MtsLevel> mtsLevels; 357 //! Precision of x in compressed trajectory file 358 real x_compression_precision; 359 //! Requested fourier_spacing, when nk? not set 360 real fourier_spacing; 361 //! Number of k vectors in x dimension for fourier methods for long range electrost. 362 int nkx; 363 //! Number of k vectors in y dimension for fourier methods for long range electrost. 364 int nky; 365 //! Number of k vectors in z dimension for fourier methods for long range electrost. 366 int nkz; 367 //! Interpolation order for PME 368 int pme_order; 369 //! Real space tolerance for Ewald, determines the real/reciprocal space relative weight 370 real ewald_rtol; 371 //! Real space tolerance for LJ-Ewald 372 real ewald_rtol_lj; 373 //! Normal/3D ewald, or pseudo-2D LR corrections 374 int ewald_geometry; 375 //! Epsilon for PME dipole correction 376 real epsilon_surface; 377 //! Type of combination rule in LJ-PME 378 int ljpme_combination_rule; 379 //! Type of periodic boundary conditions 380 PbcType pbcType; 381 //! Periodic molecules 382 bool bPeriodicMols; 383 //! Continuation run: starting state is correct (ie. constrained) 384 gmx_bool bContinuation; 385 //! Temperature coupling 386 int etc; 387 //! Interval in steps for temperature coupling 388 int nsttcouple; 389 //! Whether to print nose-hoover chains 390 gmx_bool bPrintNHChains; 391 //! Pressure coupling 392 int epc; 393 //! Pressure coupling type 394 int epct; 395 //! Interval in steps for pressure coupling 396 int nstpcouple; 397 //! Pressure coupling time (ps) 398 real tau_p; 399 //! Reference pressure (kJ/(mol nm^3)) 400 tensor ref_p; 401 //! Compressibility ((mol nm^3)/kJ) 402 tensor compress; 403 //! How to scale absolute reference coordinates 404 int refcoord_scaling; 405 //! The COM of the posres atoms 406 rvec posres_com; 407 //! The B-state COM of the posres atoms 408 rvec posres_comB; 409 //! Random seed for Andersen thermostat (obsolete) 410 int andersen_seed; 411 //! Per atom pair energy drift tolerance (kJ/mol/ps/atom) for list buffer 412 real verletbuf_tol; 413 //! Short range pairlist cut-off (nm) 414 real rlist; 415 //! Radius for test particle insertion 416 real rtpi; 417 //! Type of electrostatics treatment 418 int coulombtype; 419 //! Modify the Coulomb interaction 420 int coulomb_modifier; 421 //! Coulomb switch range start (nm) 422 real rcoulomb_switch; 423 //! Coulomb cutoff (nm) 424 real rcoulomb; 425 //! Relative dielectric constant 426 real epsilon_r; 427 //! Relative dielectric constant of the RF 428 real epsilon_rf; 429 //! Always false (no longer supported) 430 bool implicit_solvent; 431 //! Type of Van der Waals treatment 432 int vdwtype; 433 //! Modify the Van der Waals interaction 434 int vdw_modifier; 435 //! Van der Waals switch range start (nm) 436 real rvdw_switch; 437 //! Van der Waals cutoff (nm) 438 real rvdw; 439 //! Perform Long range dispersion corrections 440 int eDispCorr; 441 //! Extension of the table beyond the cut-off, as well as the table length for 1-4 interac. 442 real tabext; 443 //! Tolerance for shake 444 real shake_tol; 445 //! Free energy calculations 446 int efep; 447 //! Data for the FEP state 448 t_lambda* fepvals; 449 //! Whether to do simulated tempering 450 gmx_bool bSimTemp; 451 //! Variables for simulated tempering 452 t_simtemp* simtempvals; 453 //! Whether expanded ensembles are used 454 gmx_bool bExpanded; 455 //! Expanded ensemble parameters 456 t_expanded* expandedvals; 457 //! Type of distance restraining 458 int eDisre; 459 //! Force constant for time averaged distance restraints 460 real dr_fc; 461 //! Type of weighting of pairs in one restraints 462 int eDisreWeighting; 463 //! Use combination of time averaged and instantaneous violations 464 gmx_bool bDisreMixed; 465 //! Frequency of writing pair distances to enx 466 int nstdisreout; 467 //! Time constant for memory function in disres 468 real dr_tau; 469 //! Force constant for orientational restraints 470 real orires_fc; 471 //! Time constant for memory function in orires 472 real orires_tau; 473 //! Frequency of writing tr(SD) to energy output 474 int nstorireout; 475 //! The stepsize for updating 476 real em_stepsize; 477 //! The tolerance 478 real em_tol; 479 //! Number of iterations for convergence of steepest descent in relax_shells 480 int niter; 481 //! Stepsize for directional minimization in relax_shells 482 real fc_stepsize; 483 //! Number of steps after which a steepest descents step is done while doing cg 484 int nstcgsteep; 485 //! Number of corrections to the Hessian to keep 486 int nbfgscorr; 487 //! Type of constraint algorithm 488 int eConstrAlg; 489 //! Order of the LINCS Projection Algorithm 490 int nProjOrder; 491 //! Warn if any bond rotates more than this many degrees 492 real LincsWarnAngle; 493 //! Number of iterations in the final LINCS step 494 int nLincsIter; 495 //! Use successive overrelaxation for shake 496 gmx_bool bShakeSOR; 497 //! Friction coefficient for BD (amu/ps) 498 real bd_fric; 499 //! Random seed for SD and BD 500 int64_t ld_seed; 501 //! The number of walls 502 int nwall; 503 //! The type of walls 504 int wall_type; 505 //! The potentail is linear for r<=wall_r_linpot 506 real wall_r_linpot; 507 //! The atom type for walls 508 int wall_atomtype[2]; 509 //! Number density for walls 510 real wall_density[2]; 511 //! Scaling factor for the box for Ewald 512 real wall_ewald_zfac; 513 514 /* COM pulling data */ 515 //! Do we do COM pulling? 516 gmx_bool bPull; 517 //! The data for center of mass pulling 518 std::unique_ptr<pull_params_t> pull; 519 520 /* AWH bias data */ 521 //! Whether to use AWH biasing for PMF calculations 522 gmx_bool bDoAwh; 523 //! AWH biasing parameters 524 gmx::AwhParams* awhParams; 525 526 /* Enforced rotation data */ 527 //! Whether to calculate enforced rotation potential(s) 528 gmx_bool bRot; 529 //! The data for enforced rotation potentials 530 t_rot* rot; 531 532 //! Whether to do ion/water position exchanges (CompEL) 533 int eSwapCoords; 534 //! Swap data structure. 535 t_swapcoords* swap; 536 537 //! Whether the tpr makes an interactive MD session possible. 538 gmx_bool bIMD; 539 //! Interactive molecular dynamics 540 t_IMD* imd; 541 542 //! Acceleration for viscosity calculation 543 real cos_accel; 544 //! Triclinic deformation velocities (nm/ps) 545 tensor deform; 546 /*! \brief User determined parameters */ 547 /**@{*/ 548 int userint1; 549 int userint2; 550 int userint3; 551 int userint4; 552 real userreal1; 553 real userreal2; 554 real userreal3; 555 real userreal4; 556 /**@}*/ 557 //! Group options 558 t_grpopts opts; 559 //! QM/MM calculation 560 gmx_bool bQMMM; 561 562 /* Fields for removed features go here (better caching) */ 563 //! Whether AdResS is enabled - always false if a valid .tpr was read 564 gmx_bool bAdress; 565 //! Whether twin-range scheme is active - always false if a valid .tpr was read 566 gmx_bool useTwinRange; 567 568 //! KVT object that contains input parameters converted to the new style. 569 gmx::KeyValueTreeObject* params; 570 571 //! KVT for storing simulation parameters that are not part of the mdp file. 572 std::unique_ptr<gmx::KeyValueTreeObject> internalParameters; 573 }; 574 575 int ir_optimal_nstcalcenergy(const t_inputrec* ir); 576 577 int tcouple_min_integration_steps(int etc); 578 579 int ir_optimal_nsttcouple(const t_inputrec* ir); 580 581 int pcouple_min_integration_steps(int epc); 582 583 int ir_optimal_nstpcouple(const t_inputrec* ir); 584 585 /* Returns if the Coulomb force or potential is switched to zero */ 586 gmx_bool ir_coulomb_switched(const t_inputrec* ir); 587 588 /* Returns if the Coulomb interactions are zero beyond the rcoulomb. 589 * Note: always returns TRUE for the Verlet cut-off scheme. 590 */ 591 gmx_bool ir_coulomb_is_zero_at_cutoff(const t_inputrec* ir); 592 593 /* As ir_coulomb_is_zero_at_cutoff, but also returns TRUE for user tabulated 594 * interactions, since these might be zero beyond rcoulomb. 595 */ 596 gmx_bool ir_coulomb_might_be_zero_at_cutoff(const t_inputrec* ir); 597 598 /* Returns if the Van der Waals force or potential is switched to zero */ 599 gmx_bool ir_vdw_switched(const t_inputrec* ir); 600 601 /* Returns if the Van der Waals interactions are zero beyond the rvdw. 602 * Note: always returns TRUE for the Verlet cut-off scheme. 603 */ 604 gmx_bool ir_vdw_is_zero_at_cutoff(const t_inputrec* ir); 605 606 /* As ir_vdw_is_zero_at_cutoff, but also returns TRUE for user tabulated 607 * interactions, since these might be zero beyond rvdw. 608 */ 609 gmx_bool ir_vdw_might_be_zero_at_cutoff(const t_inputrec* ir); 610 611 /*! \brief Free memory from input record. 612 * 613 * All arrays and pointers will be freed. 614 * 615 * \param[in] ir The data structure 616 */ 617 void done_inputrec(t_inputrec* ir); 618 619 void pr_inputrec(FILE* fp, int indent, const char* title, const t_inputrec* ir, gmx_bool bMDPformat); 620 621 void cmp_inputrec(FILE* fp, const t_inputrec* ir1, const t_inputrec* ir2, real ftol, real abstol); 622 623 void comp_pull_AB(FILE* fp, const pull_params_t& pull, real ftol, real abstol); 624 625 626 gmx_bool inputrecDeform(const t_inputrec* ir); 627 628 gmx_bool inputrecDynamicBox(const t_inputrec* ir); 629 630 gmx_bool inputrecPreserveShape(const t_inputrec* ir); 631 632 gmx_bool inputrecNeedMutot(const t_inputrec* ir); 633 634 gmx_bool inputrecTwinRange(const t_inputrec* ir); 635 636 gmx_bool inputrecExclForces(const t_inputrec* ir); 637 638 gmx_bool inputrecNptTrotter(const t_inputrec* ir); 639 640 gmx_bool inputrecNvtTrotter(const t_inputrec* ir); 641 642 gmx_bool inputrecNphTrotter(const t_inputrec* ir); 643 644 /*! \brief Return true if the simulation is 2D periodic with two walls. */ 645 bool inputrecPbcXY2Walls(const t_inputrec* ir); 646 647 //! \brief Return true if the simulation has frozen atoms (non-trivial freeze groups). 648 bool inputrecFrozenAtoms(const t_inputrec* ir); 649 650 /*! \brief Returns true for MD integator with T and/or P-coupling that supports 651 * calculating a conserved energy quantity. 652 * 653 * Note that dynamical integrators without T and P coupling (ie NVE) 654 * return false, i.e. the return value refers to whether there 655 * is a conserved quantity different than the total energy. 656 */ 657 bool integratorHasConservedEnergyQuantity(const t_inputrec* ir); 658 659 /*! \brief Returns true when temperature is coupled or constant. */ 660 bool integratorHasReferenceTemperature(const t_inputrec* ir); 661 662 /*! \brief Return the number of bounded dimensions 663 * 664 * \param[in] ir The input record with MD parameters 665 * \return the number of dimensions in which 666 * the coordinates of the particles are bounded, starting at X. 667 */ 668 int inputrec2nboundeddim(const t_inputrec* ir); 669 670 /*! \brief Returns the number of degrees of freedom in center of mass motion 671 * 672 * \param[in] ir The inputrec structure 673 * \return the number of degrees of freedom of the center of mass 674 */ 675 int ndof_com(const t_inputrec* ir); 676 677 /*! \brief Returns the maximum reference temperature over all coupled groups 678 * 679 * Returns 0 for energy minimization and normal mode computation. 680 * Returns -1 for MD without temperature coupling. 681 * 682 * \param[in] ir The inputrec structure 683 */ 684 real maxReferenceTemperature(const t_inputrec& ir); 685 686 /*! \brief Returns whether there is an Ewald surface contribution 687 */ 688 bool haveEwaldSurfaceContribution(const t_inputrec& ir); 689 690 /*! \brief Check if calculation of the specific FEP type was requested. 691 * 692 * \param[in] ir Input record. 693 * \param[in] fepType Free-energy perturbation type to check for. 694 * 695 * \returns If the \p fepType is perturbed in this run. 696 */ 697 bool haveFreeEnergyType(const t_inputrec& ir, int fepType); 698 699 #endif /* GMX_MDTYPES_INPUTREC_H */ 700