#ifndef GMX_MATH_UNITS_H
#define GMX_MATH_UNITS_H
+#include <cmath>
+
/*
* Physical constants to be used in Gromacs.
* No constants (apart from 0, 1 or 2) should
* be anywhere else in the code.
*/
-#include "gromacs/math/utilities.h"
-
-#define ANGSTROM (1e-10) /* Old... */
-#define KILO (1e3) /* Thousand */
-#define NANO (1e-9) /* A Number */
-#define PICO (1e-12) /* A Number */
-#define A2NM (ANGSTROM / NANO) /* NANO */
-#define NM2A (NANO / ANGSTROM) /* 10.0 */
-#define RAD2DEG (180.0 / M_PI) /* Conversion */
-#define DEG2RAD (M_PI / 180.0) /* id */
-#define CAL2JOULE (4.184) /* Exact definition of the calorie */
-#define E_CHARGE (1.602176634e-19) /* Exact definition, Coulomb NIST 2018 CODATA */
-
-#define AMU (1.66053906660e-27) /* kg, NIST 2018 CODATA */
-#define BOLTZMANN (1.380649e-23) /* (J/K, Exact definition, NIST 2018 CODATA */
-#define AVOGADRO (6.02214076e23) /* 1/mol, Exact definition, NIST 2018 CODATA */
-#define RGAS (BOLTZMANN * AVOGADRO) /* (J/(mol K)) */
-#define BOLTZ (RGAS / KILO) /* (kJ/(mol K)) */
-#define FARADAY (E_CHARGE * AVOGADRO) /* (C/mol) */
-#define ELECTRONVOLT (E_CHARGE * AVOGADRO / KILO) /* (kJ/mol) */
-#define PLANCK1 (6.62607015e-34) /* J/Hz, Exact definition, NIST 2018 CODATA */
-#define PLANCK (PLANCK1 * AVOGADRO / (PICO * KILO)) /* (kJ/mol) ps */
-
-#define EPSILON0_SI (8.8541878128e-12) /* F/m, NIST 2018 CODATA */
+namespace gmx
+{
+
+constexpr double c_angstrom = 1e-10;
+constexpr double c_kilo = 1e3;
+constexpr double c_nano = 1e-9;
+constexpr double c_pico = 1e-12;
+constexpr double c_nm2A = c_nano / c_angstrom;
+constexpr double c_rad2Deg = 180.0 / M_PI;
+constexpr double c_deg2Rad = M_PI / 180.0;
+constexpr double c_cal2Joule = 4.184; /* Exact definition of the calorie */
+constexpr double c_electronCharge = 1.602176634e-19; /* Exact definition, Coulomb NIST 2018 CODATA */
+
+constexpr double c_amu = 1.66053906660e-27; /* kg, NIST 2018 CODATA */
+constexpr double c_boltzmann = 1.380649e-23; /* (J/K, Exact definition, NIST 2018 CODATA */
+constexpr double c_avogadro = 6.02214076e23; /* 1/mol, Exact definition, NIST 2018 CODATA */
+constexpr double c_universalGasConstant = c_boltzmann * c_avogadro; /* (J/(mol K)) */
+constexpr double c_boltz = c_universalGasConstant / c_kilo; /* (kJ/(mol K)) */
+constexpr double c_faraday = c_electronCharge * c_avogadro; /* (C/mol) */
+constexpr double c_planck1 = 6.62607015e-34; /* J/Hz, Exact definition, NIST 2018 CODATA */
+constexpr double c_planck = (c_planck1 * c_avogadro / (c_pico * c_kilo)); /* (kJ/mol) ps */
+
+constexpr double c_epsilon0Si = 8.8541878128e-12; /* F/m, NIST 2018 CODATA */
/* Epsilon in our MD units: (e^2 / Na (kJ nm)) == (e^2 mol/(kJ nm)) */
-#define EPSILON0 ((EPSILON0_SI * NANO * KILO) / (E_CHARGE * E_CHARGE * AVOGADRO))
+constexpr double c_epsilon0 =
+ ((c_epsilon0Si * c_nano * c_kilo) / (c_electronCharge * c_electronCharge * c_avogadro));
-#define SPEED_OF_LIGHT (2.99792458e05) /* units of nm/ps, Exact definition, NIST 2018 CODATA */
-#define ATOMICMASS_keV (931494.10242) /* Atomic mass in keV, NIST 2018 CODATA */
-#define ELECTRONMASS_keV (510.99895000) /* Electron mas in keV, NIST 2018 CODATA */
+constexpr double c_speedOfLight =
+ 2.99792458e05; /* units of nm/ps, Exact definition, NIST 2018 CODATA */
-#define RYDBERG (1.0973731568160e-02) /* nm^-1, NIST 2018 CODATA */
+constexpr double c_rydberg = 1.0973731568160e-02; /* nm^-1, NIST 2018 CODATA */
-#define ONE_4PI_EPS0 (1.0 / (4.0 * M_PI * EPSILON0))
-#define FACEL (10.0 * ONE_4PI_EPS0)
+constexpr double c_one4PiEps0 = (1.0 / (4.0 * M_PI * c_epsilon0));
/* Pressure in MD units is:
- * 1 bar = 1e5 Pa = 1e5 kg m^-1 s^-2 = 1e-28 kg nm^-1 ps^-2 = 1e-28 / AMU amu nm^1 ps ^2
+ * 1 bar = 1e5 Pa = 1e5 kg m^-1 s^-2 = 1e-28 kg nm^-1 ps^-2 = 1e-28 / amu amu nm^1 ps ^2
*/
-#define BAR_MDUNITS (1e5 * NANO * PICO * PICO / AMU)
-#define PRESFAC (1.0 / BAR_MDUNITS)
+constexpr double c_barMdunits = (1e5 * c_nano * c_pico * c_pico / c_amu);
+constexpr double c_presfac = 1.0 / c_barMdunits;
-/* DEBYE2ENM should be (1e-21*PICO)/(SPEED_OF_LIGHT*E_CHARGE*NANO*NANO),
+/* c_debye2Enm should be (1e-21*c_pico)/(c_speedOfLight*c_electronCharge*c_nano*c_nano),
* but we need to factor out some of the exponents to avoid single-precision overflows.
*/
-#define DEBYE2ENM (1e-15 / (SPEED_OF_LIGHT * E_CHARGE))
-#define ENM2DEBYE (1.0 / DEBYE2ENM)
+constexpr double c_debye2Enm = (1e-15 / (c_speedOfLight * c_electronCharge));
+constexpr double c_enm2Debye = 1.0 / c_debye2Enm;
/* to convert from a acceleration in (e V)/(amu nm) */
-/* FIELDFAC is also Faraday's constant and E_CHARGE/(1e6 AMU) */
-#define FIELDFAC (FARADAY / KILO)
+/* c_fieldfac is also Faraday's constant and c_electronCharge/(1e6 amu) */
+constexpr double c_fieldfac = c_faraday / c_kilo;
/* to convert AU to MD units: */
-#define HARTREE2KJ ((2.0 * RYDBERG * PLANCK * SPEED_OF_LIGHT) / AVOGADRO)
-#define BOHR2NM (0.0529177210903) /* nm^-1, NIST 2018 CODATA */
-#define HARTREE_BOHR2MD (HARTREE2KJ * AVOGADRO / BOHR2NM)
+constexpr double c_hartree2Kj = ((2.0 * c_rydberg * c_planck * c_speedOfLight) / c_avogadro);
+constexpr double c_bohr2Nm = 0.0529177210903; /* nm^-1, NIST 2018 CODATA */
+constexpr double c_hartreeBohr2Md = (c_hartree2Kj * c_avogadro / c_bohr2Nm);
+} // namespace gmx
/* The four basic units */
#define unit_length "nm"
#define unit_mass "u"
#define unit_energy "kJ/mol"
-/* Temperature unit, T in this unit times BOLTZ give energy in unit_energy */
+/* Temperature unit, T in this unit times c_boltz give energy in unit_energy */
#define unit_temp_K "K"
-/* Charge unit, electron charge, involves ONE_4PI_EPS0 */
+/* Charge unit, electron charge, involves c_one4PiEps0 */
#define unit_charge_e "e"
-/* Pressure unit, pressure in basic units times PRESFAC gives this unit */
+/* Pressure unit, pressure in basic units times c_presfac gives this unit */
#define unit_pres_bar "bar"
-/* Dipole unit, debye, conversion from the unit_charge_e involves ENM2DEBYE */
+/* Dipole unit, debye, conversion from the unit_charge_e involves c_enm2Debye */
#define unit_dipole_D "D"
/* Derived units from basic units only */
/* Other derived units */
#define unit_surft_bar unit_pres_bar " " unit_length
-/* SI units, conversion from basic units involves NANO, PICO and AMU */
+/* SI units, conversion from basic units involves c_nano, c_pico and amu */
#define unit_length_SI "m"
#define unit_time_SI "s"
#define unit_mass_SI "kg"
}
/* Initialize all the biases */
- const double beta = 1 / (BOLTZ * inputRecord.opts.ref_t[0]);
+ const double beta = 1 / (gmx::c_boltz * inputRecord.opts.ref_t[0]);
for (int k = 0; k < awhParams.numBias; k++)
{
const AwhBiasParams& awhBiasParams = awhParams.awhBiasParams[k];
GMX_THROW(InvalidInputError(
"Pull geometry 'direction-periodic' is not supported by AWH"));
}
- double conversionFactor = pull_coordinate_is_angletype(&pullCoord) ? DEG2RAD : 1;
+ double conversionFactor = pull_coordinate_is_angletype(&pullCoord) ? gmx::c_deg2Rad : 1;
pullCoordIndex.push_back(awhDimParams.coordIndex);
dimParams.push_back(DimParams::pullDimParams(
conversionFactor, awhDimParams.forceConstant, beta));
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2015,2016,2017,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2015,2016,2017,2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
* Builder function for pull dimension parameters.
*
* \param[in] conversionFactor Conversion factor from user coordinate units to bias internal
- * units (=DEG2RAD for angles).
+ * units (=c_deg2Rad for angles).
* \param[in] forceConstant The harmonic force constant.
* \param[in] beta 1/(k_B T).
*/
for (int m = 0; (m < DIM); m++)
{
- const real fieldStrength = FIELDFAC * field(m, t);
+ const real fieldStrength = gmx::c_fieldfac * field(m, t);
if (fieldStrength != 0)
{
}
/* 1/(Vol*e0) */
- real scaleRecip = 4.0 * M_PI / (boxDiag[XX] * boxDiag[YY] * boxDiag[ZZ]) * ONE_4PI_EPS0 / ir.epsilon_r;
+ real scaleRecip =
+ 4.0 * M_PI / (boxDiag[XX] * boxDiag[YY] * boxDiag[ZZ]) * gmx::c_one4PiEps0 / ir.epsilon_r;
snew(eir, et->kmax);
for (n = 0; n < et->kmax; n++)
/* Apply charge correction */
vol = box[XX][XX] * box[YY][YY] * box[ZZ][ZZ];
- fac = M_PI * ONE_4PI_EPS0
+ fac = M_PI * gmx::c_one4PiEps0
/ (fr->ic->epsilon_r * 2.0 * vol * vol * gmx::square(fr->ic->ewaldcoeff_q));
qs2A = fr->qsum[0] * fr->qsum[0];
EwaldBoxZScaler boxScaler(ir);
boxScaler.scaleBox(box, scaledBox);
- one_4pi_eps = ONE_4PI_EPS0 / fr.ic->epsilon_r;
+ one_4pi_eps = gmx::c_one4PiEps0 / fr.ic->epsilon_r;
Vexcl_q = 0;
dvdl_excl_q = 0;
Vdipole[0] = 0;
*/
for (i = 0; (i < DIM); i++)
{
- mutot[0][i] = mu_tot[0][i] * DEBYE2ENM;
- mutot[1][i] = mu_tot[1][i] * DEBYE2ENM;
+ mutot[0][i] = mu_tot[0][i] * gmx::c_debye2Enm;
+ mutot[1][i] = mu_tot[1][i] * gmx::c_debye2Enm;
dipcorrA[i] = 0;
dipcorrB[i] = 0;
}
case EwaldGeometry::ThreeD:
if (ir.epsilon_surface != 0)
{
- dipole_coeff = 2 * M_PI * ONE_4PI_EPS0
+ dipole_coeff = 2 * M_PI * gmx::c_one4PiEps0
/ ((2 * ir.epsilon_surface + fr.ic->epsilon_r) * boxVolume);
for (i = 0; (i < DIM); i++)
{
GMX_ASSERT(pmeGpu->common->epsilon_r != 0.0F, "PME GPU: bad electrostatic coefficient");
auto* kernelParamsPtr = pme_gpu_get_kernel_params_base_ptr(pmeGpu);
- kernelParamsPtr->constants.elFactor = ONE_4PI_EPS0 / pmeGpu->common->epsilon_r;
+ kernelParamsPtr->constants.elFactor = gmx::c_one4PiEps0 / pmeGpu->common->epsilon_r;
}
void pme_gpu_get_real_grid_sizes(const PmeGpu* pmeGpu, gmx::IVec* gridSize, gmx::IVec* paddedGridSize)
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2016,2017,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2016,2017,2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
*/
struct PmeGpuConstParams
{
- /*! \brief Electrostatics coefficient = ONE_4PI_EPS0 / pme->epsilon_r */
+ /*! \brief Electrostatics coefficient = c_one4PiEps0 / pme->epsilon_r */
float elFactor;
/*! \brief Virial and energy GPU array. Size is c_virialAndEnergyCount (7) floats.
* The element order is virxx, viryy, virzz, virxy, virxz, viryz, energy. */
ivec local_ndata, local_offset, local_size;
real elfac;
- elfac = ONE_4PI_EPS0 / pme->epsilon_r;
+ elfac = gmx::c_one4PiEps0 / pme->epsilon_r;
nx = pme->nkx;
ny = pme->nky;
if (norm2(box[YY]) * norm2(box[ZZ]) != 0)
{
- alpha = RAD2DEG * gmx_angle(box[YY], box[ZZ]);
+ alpha = gmx::c_rad2Deg * gmx_angle(box[YY], box[ZZ]);
}
else
{
}
if (norm2(box[XX]) * norm2(box[ZZ]) != 0)
{
- beta = RAD2DEG * gmx_angle(box[XX], box[ZZ]);
+ beta = gmx::c_rad2Deg * gmx_angle(box[XX], box[ZZ]);
}
else
{
}
if (norm2(box[XX]) * norm2(box[YY]) != 0)
{
- gamma = RAD2DEG * gmx_angle(box[XX], box[YY]);
+ gamma = gmx::c_rad2Deg * gmx_angle(box[XX], box[YY]);
}
else
{
{
if (alpha != 90.0)
{
- cosa = std::cos(alpha * DEG2RAD);
+ cosa = std::cos(alpha * gmx::c_deg2Rad);
}
else
{
}
if (beta != 90.0)
{
- cosb = std::cos(beta * DEG2RAD);
+ cosb = std::cos(beta * gmx::c_deg2Rad);
}
else
{
}
if (gamma != 90.0)
{
- cosg = std::cos(gamma * DEG2RAD);
- sing = std::sin(gamma * DEG2RAD);
+ cosg = std::cos(gamma * gmx::c_deg2Rad);
+ sing = std::sin(gamma * gmx::c_deg2Rad);
}
else
{
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2013,2014,2015,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#if GMX_USE_TNG
gmx_tng_add_mtop(gmx_tng, mtop);
set_writing_intervals(gmx_tng, FALSE, ir);
- tng_time_per_frame_set(gmx_tng->tng, ir->delta_t * PICO);
+ tng_time_per_frame_set(gmx_tng->tng, ir->delta_t * gmx::c_pico);
gmx_tng->timePerFrameIsSet = true;
#else
GMX_UNUSED_VALUE(gmx_tng);
gmx_tng_add_mtop(gmx_tng, mtop);
add_selection_groups(gmx_tng, mtop);
set_writing_intervals(gmx_tng, TRUE, ir);
- tng_time_per_frame_set(gmx_tng->tng, ir->delta_t * PICO);
+ tng_time_per_frame_set(gmx_tng->tng, ir->delta_t * gmx::c_pico);
gmx_tng->timePerFrameIsSet = true;
gmx_tng_set_compression_precision(gmx_tng, ir->x_compression_precision);
#else
# else
static write_data_func_pointer write_data = tng_util_generic_with_time_write;
# endif
- double elapsedSeconds = elapsedPicoSeconds * PICO;
+ double elapsedSeconds = elapsedPicoSeconds * gmx::c_pico;
int64_t nParticles;
char compression;
tng_num_frames_get(tng, &nFrames);
tng_util_time_of_frame_get(tng, nFrames - 1, &time);
- fTime = time / PICO;
+ fTime = time / gmx::c_pico;
return fTime;
#else
GMX_UNUSED_VALUE(gmx_tng);
// GROMACS expects distances in nm
switch (exp)
{
- case 9: distanceScaleFactor = NANO / NANO; break;
- case 10: distanceScaleFactor = NANO / ANGSTROM; break;
+ case 9: distanceScaleFactor = gmx::c_nano / gmx::c_nano; break;
+ case 10: distanceScaleFactor = gmx::c_nano / gmx::c_angstrom; break;
default: distanceScaleFactor = pow(10.0, exp + 9.0);
}
fr->bStep = TRUE;
// Convert the time to ps
- fr->time = frameTime / PICO;
+ fr->time = frameTime / gmx::c_pico;
fr->bTime = (frameTime > 0);
// TODO This does not leak, but is not exception safe.
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
static int calc_Nbin(real phi, int multiplicity, real core_frac)
{
- static const real r360 = 360 * DEG2RAD;
+ static const real r360 = 360 * gmx::c_deg2Rad;
real rot_width, core_width, core_offset, low, hi;
int bin;
/* with multiplicity 3 and core_frac 0.5
{
low = ((bin - 1) * rot_width) + core_offset;
hi = ((bin - 1) * rot_width) + core_offset + core_width;
- low *= DEG2RAD;
- hi *= DEG2RAD;
+ low *= gmx::c_deg2Rad;
+ hi *= gmx::c_deg2Rad;
if ((phi > low) && (phi < hi))
{
return bin;
for (i = 0; i < nangles; i++)
{
- angle = angles[i] * RAD2DEG;
+ angle = angles[i] * gmx::c_rad2Deg;
if (angle > 135 && angle < 225)
{
#ifndef GMX_GMXANA_ANGLE_CORRECTION_H
#define GMX_GMXANA_ANGLE_CORRECTION_H
-#include "gromacs/math/units.h"
#include "gromacs/utility/real.h"
/*! \brief
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
out = xvgropen(opt2fn("-ov", NFILE, fnm), title, "Time (ps)", "Angle (degrees)", oenv);
for (i = 0; (i < nframes); i++)
{
- fprintf(out, "%10.5f %8.3f", time[i], aver_angle[i] * RAD2DEG);
+ fprintf(out, "%10.5f %8.3f", time[i], aver_angle[i] * gmx::c_rad2Deg);
if (bALL)
{
for (j = 0; (j < nangles); j++)
if (bPBC)
{
real dd = dih[j][i];
- fprintf(out, " %8.3f", std::atan2(std::sin(dd), std::cos(dd)) * RAD2DEG);
+ fprintf(out, " %8.3f", std::atan2(std::sin(dd), std::cos(dd)) * gmx::c_rad2Deg);
}
else
{
- fprintf(out, " %8.3f", dih[j][i] * RAD2DEG);
+ fprintf(out, " %8.3f", dih[j][i] * gmx::c_rad2Deg);
}
}
}
if (bChandler)
{
- real dval, sixty = DEG2RAD * 60;
+ real dval, sixty = gmx::c_deg2Rad * 60;
gmx_bool bTest;
for (i = 0; (i < nangles); i++)
}
aver /= nframes;
double aversig = correctRadianAngleRange(aver);
- aversig *= RAD2DEG;
- aver *= RAD2DEG;
+ aversig *= gmx::c_rad2Deg;
+ aver *= gmx::c_rad2Deg;
printf(" < angle > = %g\n", aversig);
if (mult == 3)
* This file is part of the GROMACS molecular simulation package.
*
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
AwhGraphSelection awhGraphSelection =
(moreGraphs ? AwhGraphSelection::All : AwhGraphSelection::Pmf);
EnergyUnit energyUnit = (kTUnit ? EnergyUnit::KT : EnergyUnit::KJPerMol);
- awhReader = std::make_unique<AwhReader>(
- ir.awhParams, nfile, fnm, awhGraphSelection, energyUnit, BOLTZ * ir.opts.ref_t[0], block);
+ awhReader = std::make_unique<AwhReader>(ir.awhParams,
+ nfile,
+ fnm,
+ awhGraphSelection,
+ energyUnit,
+ gmx::c_boltz * ir.opts.ref_t[0],
+ block);
}
awhReader->processAwhFrame(*block, frame->t, oenv);
double DG0, DG1, DG2, dDG1;
double n1, n2; /* numbers of samples as doubles */
- kT = BOLTZ * temp;
+ kT = gmx::c_boltz * temp;
beta = 1 / kT;
/* count the numbers of samples */
double Wfac1, Wfac2;
double n1, n2;
- kT = BOLTZ * temp;
+ kT = gmx::c_boltz * temp;
beta = 1 / kT;
/* count the numbers of samples */
double Wfac1, Wfac2;
double n1, n2;
- kT = BOLTZ * temp;
+ kT = gmx::c_boltz * temp;
beta = 1 / kT;
/* count the numbers of samples */
}
/* print results in kT */
- kT = BOLTZ * temp;
+ kT = gmx::c_boltz * temp;
printf("\nTemperature: %g K\n", temp);
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
fprintf(flen, " %6g", bun.len[i]);
fprintf(fdist, " %6g", norm(bun.mid[i]));
fprintf(fz, " %6g", bun.mid[i][ZZ]);
- fprintf(ftilt, " %6g", RAD2DEG * acos(bun.dir[i][ZZ]));
+ fprintf(ftilt, " %6g", gmx::c_rad2Deg * acos(bun.dir[i][ZZ]));
comp = bun.mid[i][XX] * bun.dir[i][XX] + bun.mid[i][YY] * bun.dir[i][YY];
- fprintf(ftiltr, " %6g", RAD2DEG * std::asin(comp / std::hypot(comp, bun.dir[i][ZZ])));
+ fprintf(ftiltr, " %6g", gmx::c_rad2Deg * std::asin(comp / std::hypot(comp, bun.dir[i][ZZ])));
comp = bun.mid[i][YY] * bun.dir[i][XX] - bun.mid[i][XX] * bun.dir[i][YY];
- fprintf(ftiltl, " %6g", RAD2DEG * std::asin(comp / std::hypot(comp, bun.dir[i][ZZ])));
+ fprintf(ftiltl, " %6g", gmx::c_rad2Deg * std::asin(comp / std::hypot(comp, bun.dir[i][ZZ])));
if (bKink)
{
rvec_sub(bun.end[0][i], bun.end[2][i], va);
rvec_sub(bun.end[2][i], bun.end[1][i], vb);
unitv(va, va);
unitv(vb, vb);
- fprintf(fkink, " %6g", RAD2DEG * acos(iprod(va, vb)));
+ fprintf(fkink, " %6g", gmx::c_rad2Deg * acos(iprod(va, vb)));
cprod(va, vb, vc);
copy_rvec(bun.mid[i], vr);
vr[ZZ] = 0;
unitv(vr, vr);
- fprintf(fkinkr, " %6g", RAD2DEG * std::asin(iprod(vc, vr)));
+ fprintf(fkinkr, " %6g", gmx::c_rad2Deg * std::asin(iprod(vc, vr)));
vl[XX] = vr[YY];
vl[YY] = -vr[XX];
vl[ZZ] = 0;
- fprintf(fkinkl, " %6g", RAD2DEG * std::asin(iprod(vc, vl)));
+ fprintf(fkinkl, " %6g", gmx::c_rad2Deg * std::asin(iprod(vc, vl)));
}
}
fprintf(flen, "\n");
#define NPP asize(map)
int x, y;
-#define INDEX(ppp) (((static_cast<int>(360 + (ppp)*RAD2DEG)) % 360) / 6)
+#define INDEX(ppp) (((static_cast<int>(360 + (ppp)*gmx::c_rad2Deg)) % 360) / 6)
x = INDEX(phi);
y = INDEX(psi);
#undef INDEX
Psi = dlist[i].j0[edPsi];
for (j = 0; (j < nf); j++)
{
- phi = RAD2DEG * dih[Phi][j];
- psi = RAD2DEG * dih[Psi][j];
+ phi = gmx::c_rad2Deg * dih[Phi][j];
+ psi = gmx::c_rad2Deg * dih[Psi][j];
fprintf(fp, "%10g %10g\n", phi, psi);
if (bViol)
{
- fprintf(gp, "%d\n", static_cast<int>(!bAllowed(dih[Phi][j], RAD2DEG * dih[Psi][j])));
+ fprintf(gp,
+ "%d\n",
+ static_cast<int>(!bAllowed(dih[Phi][j], gmx::c_rad2Deg * dih[Psi][j])));
}
if (bOm)
{
- omega = RAD2DEG * dih[Om][j];
+ omega = gmx::c_rad2Deg * dih[Om][j];
mat[static_cast<int>(((phi * NMAT) / 360) + gmx::exactDiv(NMAT, 2))]
[static_cast<int>(((psi * NMAT) / 360) + gmx::exactDiv(NMAT, 2))] += omega;
}
Xi2 = dlist[i].j0[edChi2];
for (j = 0; (j < nf); j++)
{
- fprintf(fp, "%10g %10g\n", RAD2DEG * dih[Xi1][j], RAD2DEG * dih[Xi2][j]);
+ fprintf(fp, "%10g %10g\n", gmx::c_rad2Deg * dih[Xi1][j], gmx::c_rad2Deg * dih[Xi2][j]);
}
xvgrclose(fp);
}
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2007, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
ekin += 0.5 * m * iprod(v[ai], v[ai]);
}
}
- temp = (ekin * 2.0) / (3.0 * tfac * max_clust_size * BOLTZ);
+ temp = (ekin * 2.0) / (3.0 * tfac * max_clust_size * gmx::c_boltz);
fprintf(tp, "%10.3f %10.3f\n", frameTime, temp);
}
}
*
* Copyright (c) 2008,2009,2010,2011,2012 by the GROMACS development team.
* Copyright (c) 2013,2014,2015,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/smalloc.h"
-#define EPSI0 (EPSILON0 * E_CHARGE * E_CHARGE * AVOGADRO / (KILO * NANO)) /* EPSILON0 in SI units */
+constexpr double EPSI0 = (gmx::c_epsilon0 * gmx::c_electronCharge * gmx::c_electronCharge * gmx::c_avogadro
+ / (gmx::c_kilo * gmx::c_nano)); /* c_epsilon0 in SI units */
static void index_atom2mol(int* n, int* index, t_block* mols)
{
volume_av /= refr;
prefactor = 1.0;
- prefactor /= 3.0 * EPSILON0 * volume_av * BOLTZ * temp;
+ prefactor /= 3.0 * gmx::c_epsilon0 * volume_av * gmx::c_boltz * temp;
- prefactorav = E_CHARGE * E_CHARGE;
- prefactorav /= volume_av * BOLTZMANN * temp * NANO * 6.0;
+ prefactorav = gmx::c_electronCharge * gmx::c_electronCharge;
+ prefactorav /= volume_av * gmx::c_boltzmann * temp * gmx::c_nano * 6.0;
fprintf(stderr, "Prefactor fit E-H: 1 / 6.0*V*k_B*T: %g\n", prefactorav);
{
printf("\nCalculating current autocorrelation ... \n");
- sgk = calc_cacf(outf, prefactorav / PICO, cacf, time, nvfr, vfr, ie, nshift);
+ sgk = calc_cacf(outf, prefactorav / gmx::c_pico, cacf, time, nvfr, vfr, ie, nshift);
if (ie > ii)
{
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
}
if (dens_opt[0][0] == 'm')
{
- fprintf(den, " %12g", ddd * AMU / (NANO * NANO * NANO));
+ fprintf(den, " %12g", ddd * gmx::c_amu / (gmx::c_nano * gmx::c_nano * gmx::c_nano));
}
else
{
* This file is part of the GROMACS molecular simulation package.
*
* Copyright (c) 2010-2018, The GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
(*Densdevel)[*tblock] = Densslice;
}
- dscale = (*xslices) * (*yslices) * (*zslices) * AMU
- / (box[ax1][ax1] * box[ax2][ax2] * box[axis][axis] * nsttblock * (NANO * NANO * NANO));
+ dscale = (*xslices) * (*yslices) * (*zslices) * gmx::c_amu
+ / (box[ax1][ax1] * box[ax2][ax2] * box[axis][axis] * nsttblock
+ * (gmx::c_nano * gmx::c_nano * gmx::c_nano));
if (bCenter)
{
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/smalloc.h"
-#define e2d(x) ENM2DEBYE*(x)
-#define EANG2CM (E_CHARGE * 1.0e-10) /* e Angstrom to Coulomb meter */
-#define CM2D (SPEED_OF_LIGHT * 1.0e+24) /* Coulomb meter to Debye */
+#define e2d(x) gmx::c_enm2Debye*(x)
+constexpr double EANG2CM = gmx::c_electronCharge * 1.0e-10; /* e Angstrom to Coulomb meter */
+constexpr double CM2D = gmx::c_speedOfLight * 1.0e+24; /* Coulomb meter to Debye */
typedef struct
{
{
cosav = 0;
}
- ener = -0.5 * cosav * ONE_4PI_EPS0 / (x1 * x1 * x1);
+ ener = -0.5 * cosav * gmx::c_one4PiEps0 / (x1 * x1 * x1);
fprintf(fp, "%10.5e %12.5e %12.5e %12.5e %12.5e %12.5e\n", x1, Gkr, cosav, hOO, gOO, ener);
double eps_0 = 8.854187817e-12; /* epsilon_0 in C^2 J^-1 m^-1 */
double fac = 1.112650021e-59; /* converts Debye^2 to C^2 m^2 */
- A = M_diff * fac / (3 * eps_0 * volume * NANO * NANO * NANO * BOLTZMANN * temp);
+ A = M_diff * fac
+ / (3 * eps_0 * volume * gmx::c_nano * gmx::c_nano * gmx::c_nano * gmx::c_boltzmann * temp);
if (epsRF == 0.0)
{
* This file is part of the GROMACS molecular simulation package.
*
* Copyright (c) 2011-2018, The GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
{
double fy = f * y;
- return BOLTZ * (std::log(calc_compress(fy)) + fy * (3 * fy - 4) / gmx::square(1 - fy));
+ return gmx::c_boltz * (std::log(calc_compress(fy)) + fy * (3 * fy - 4) / gmx::square(1 - fy));
}
static real wCsolid(real nu, real beta)
{
- real bhn = beta * PLANCK * nu;
+ real bhn = beta * gmx::c_planck * nu;
real ebn, koko;
if (bhn == 0)
static real wSsolid(real nu, real beta)
{
- real bhn = beta * PLANCK * nu;
+ real bhn = beta * gmx::c_planck * nu;
if (bhn == 0)
{
static real wAsolid(real nu, real beta)
{
- real bhn = beta * PLANCK * nu;
+ real bhn = beta * gmx::c_planck * nu;
if (bhn == 0)
{
static real wEsolid(real nu, real beta)
{
- real bhn = beta * PLANCK * nu;
+ real bhn = beta * gmx::c_planck * nu;
if (bhn == 0)
{
return 0;
}
- beta = 1 / (Temp * BOLTZ);
+ beta = 1 / (Temp * gmx::c_boltz);
fplog = gmx_fio_fopen(ftp2fn(efLOG, NFILE, fnm), "w");
fprintf(fplog, "Doing density of states analysis based on trajectory.\n");
DoS0 = dos[DOS][0];
/* Note this eqn. is incorrect in Pascal2011a! */
- Delta = ((2 * DoS0 / (9 * Natom)) * std::sqrt(M_PI * BOLTZ * Temp * Natom / tmass)
+ Delta = ((2 * DoS0 / (9 * Natom)) * std::sqrt(M_PI * gmx::c_boltz * Temp * Natom / tmass)
* std::pow((Natom / V), 1.0 / 3.0) * std::pow(6.0 / M_PI, 2.0 / 3.0));
f = calc_fluidicity(Delta, toler);
y = calc_y(f, Delta, toler);
z = calc_compress(y);
- Sig = BOLTZ
- * (5.0 / 2.0 + std::log(2 * M_PI * BOLTZ * Temp / (gmx::square(PLANCK)) * V / (f * Natom)));
+ Sig = gmx::c_boltz
+ * (5.0 / 2.0
+ + std::log(2 * M_PI * gmx::c_boltz * Temp / (gmx::square(gmx::c_planck)) * V / (f * Natom)));
Shs = Sig + calc_Shs(f, y);
- rho = (tmass * AMU) / (V * NANO * NANO * NANO);
+ rho = (tmass * gmx::c_amu) / (V * gmx::c_nano * gmx::c_nano * gmx::c_nano);
sigHS = std::cbrt(6 * y * V / (M_PI * Natom));
fprintf(fplog, "System = \"%s\"\n", *top.name);
"\\f{4}S(\\f{12}n\\f{4})",
oenv);
xvgr_legend(fp, asize(DoSlegend), DoSlegend, oenv);
- recip_fac = bRecip ? (1e7 / SPEED_OF_LIGHT) : 1.0;
+ recip_fac = bRecip ? (1e7 / gmx::c_speedOfLight) : 1.0;
for (j = 0; (j < nframes / 4); j++)
{
dos[DOS_DIFF][j] = DoS0 / (1 + gmx::square(DoS0 * M_PI * nu[j] / (6 * f * Natom)));
/* Finally analyze the results! */
wCdiff = 0.5;
- wSdiff = Shs / (3 * BOLTZ); /* Is this correct? */
+ wSdiff = Shs / (3 * gmx::c_boltz); /* Is this correct? */
wEdiff = 0.5;
wAdiff = wEdiff - wSdiff;
for (j = 0; (j < nframes / 4); j++)
fprintf(fplog, "Diffusion coefficient from VACF %g 10^-5 cm^2/s\n", DiffCoeff);
fprintf(fplog, "Diffusion coefficient from DoS %g 10^-5 cm^2/s\n", 1000 * DoS0 / (12 * tmass * beta));
- cP = BOLTZ * evaluate_integral(nframes / 4, nu, dos[DOS_CP], nullptr, int{ nframes / 4 }, &stddev);
+ cP = gmx::c_boltz
+ * evaluate_integral(nframes / 4, nu, dos[DOS_CP], nullptr, int{ nframes / 4 }, &stddev);
fprintf(fplog, "Heat capacity %g J/mol K\n", 1000 * cP / Nmol);
fprintf(fplog, "\nArrivederci!\n");
gmx_fio_fclose(fplog);
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
}
eaver[i] /= nenergy;
}
- beta = 1.0 / (BOLTZ * reftemp);
+ beta = 1.0 / (gmx::c_boltz * reftemp);
snew(efree, ngroups);
snew(edif, ngroups);
for (i = 0; (i < ngroups); i++)
}
}
/* Convert to SI for the viscosity */
- fac = (V * NANO * NANO * NANO * PICO * 1e10) / (2 * BOLTZMANN * T) / (nint - i);
+ fac = (V * gmx::c_nano * gmx::c_nano * gmx::c_nano * gmx::c_pico * 1e10)
+ / (2 * gmx::c_boltzmann * T) / (nint - i);
fprintf(fp0, "%10g", i * dt);
for (m = 0; (m <= nsets); m++)
{
const char* my_ener[] = { "Volume", "Enthalpy", "Temperature", "Total Energy" };
int ii[eNR];
- NANO3 = NANO * NANO * NANO;
+ NANO3 = gmx::c_nano * gmx::c_nano * gmx::c_nano;
if (!bDriftCorr)
{
fprintf(fp,
{
vv = edat->s[ii[eVol]].av * NANO3;
varv = gmx::square(edat->s[ii[eVol]].rmsd * NANO3);
- kappa = (varv / vv) / (BOLTZMANN * tt);
+ kappa = (varv / vv) / (gmx::c_boltzmann * tt);
}
/* Enthalpy */
hh = varh = NOTSET;
if ((ii[eEnth] < nset) && (ii[eTemp] < nset))
{
- hh = KILO * edat->s[ii[eEnth]].av / AVOGADRO;
- varh = gmx::square(KILO * edat->s[ii[eEnth]].rmsd / AVOGADRO);
- cp = AVOGADRO * ((varh / nmol) / (BOLTZMANN * tt * tt));
+ hh = gmx::c_kilo * edat->s[ii[eEnth]].av / gmx::c_avogadro;
+ varh = gmx::square(gmx::c_kilo * edat->s[ii[eEnth]].rmsd / gmx::c_avogadro);
+ cp = gmx::c_avogadro * ((varh / nmol) / (gmx::c_boltzmann * tt * tt));
}
/* Total energy */
if ((ii[eEtot] < nset) && (hh == NOTSET) && (tt != NOTSET))
by checking whether the enthalpy was computed.
*/
varet = gmx::square(edat->s[ii[eEtot]].rmsd);
- cv = KILO * ((varet / nmol) / (BOLTZ * tt * tt));
+ cv = gmx::c_kilo * ((varet / nmol) / (gmx::c_boltz * tt * tt));
}
/* Alpha, dcp */
if ((ii[eVol] < nset) && (ii[eEnth] < nset) && (ii[eTemp] < nset))
for (j = 0; (j < edat->nframes); j++)
{
v = edat->s[ii[eVol]].ener[j] * NANO3;
- h = KILO * edat->s[ii[eEnth]].ener[j] / AVOGADRO;
+ h = gmx::c_kilo * edat->s[ii[eEnth]].ener[j] / gmx::c_avogadro;
v_sum += v;
h_sum += h;
vh_sum += (v * h);
vh_aver = vh_sum / edat->nframes;
v_aver = v_sum / edat->nframes;
h_aver = h_sum / edat->nframes;
- alpha = (vh_aver - v_aver * h_aver) / (v_aver * BOLTZMANN * tt * tt);
- dcp = (v_aver * AVOGADRO / nmol) * tt * gmx::square(alpha) / (kappa);
+ alpha = (vh_aver - v_aver * h_aver) / (v_aver * gmx::c_boltzmann * tt * tt);
+ dcp = (v_aver * gmx::c_avogadro / nmol) * tt * gmx::square(alpha) / (kappa);
}
if (tt != NOTSET)
{
if (varv != NOTSET)
{
- fprintf(fp, "varv = %10g (m^6)\n", varv * AVOGADRO / nmol);
+ fprintf(fp, "varv = %10g (m^6)\n", varv * gmx::c_avogadro / nmol);
}
}
if (vv != NOTSET)
{
- fprintf(fp, "Volume = %10g m^3/mol\n", vv * AVOGADRO / nmol);
+ fprintf(fp, "Volume = %10g m^3/mol\n", vv * gmx::c_avogadro / nmol);
}
if (varh != NOTSET)
{
- fprintf(fp, "Enthalpy = %10g kJ/mol\n", hh * AVOGADRO / (KILO * nmol));
+ fprintf(fp,
+ "Enthalpy = %10g kJ/mol\n",
+ hh * gmx::c_avogadro / (gmx::c_kilo * nmol));
}
if (alpha != NOTSET)
{
expEtot = 0;
if (bFee)
{
- beta = 1.0 / (BOLTZ * reftemp);
+ beta = 1.0 / (gmx::c_boltz * reftemp);
snew(fee, nset);
}
for (i = 0; (i < nset); i++)
0.0,
0);
- factor = (Vaver * 1e-26 / (BOLTZMANN * Temp)) * Dt;
+ factor = (Vaver * 1e-26 / (gmx::c_boltzmann * Temp)) * Dt;
fp = xvgropen(visfn, buf, "Time (ps)", "\\8h\\4 (cp)", oenv);
xvgr_legend(fp, asize(leg), leg, oenv);
}
fprintf(stdout, "\n%-24s %10s\n", "Energy", "dF = -kT ln < exp(-(EB-EA)/kT) >A");
sum = 0;
- beta = 1.0 / (BOLTZ * reftemp);
+ beta = 1.0 / (gmx::c_boltz * reftemp);
for (i = 0; i < nset; i++)
{
if (gmx_strcasecmp(leg[i], enm[set[i]].name) != 0)
sum += std::exp(-dE * beta);
if (fp)
{
- fprintf(fp, "%10g %10g %10g\n", time[j], dE, -BOLTZ * reftemp * std::log(sum / (j + 1)));
+ fprintf(fp, "%10g %10g %10g\n", time[j], dE, -gmx::c_boltz * reftemp * std::log(sum / (j + 1)));
}
}
- aver = -BOLTZ * reftemp * std::log(sum / nenergy);
+ aver = -gmx::c_boltz * reftemp * std::log(sum / nenergy);
fprintf(stdout, "%-24s %10g\n", leg[i], aver);
}
if (fp)
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2008, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
{
real kbt;
- kbt = BOLTZ * temp;
+ kbt = gmx::c_boltz * temp;
if (tau <= 0)
{
return -666;
}
else
{
- return kbt * std::log(kbt * tau / PLANCK);
+ return kbt * std::log(kbt * tau / gmx::c_planck);
}
}
compute_weighted_rates(
n, t, ct, nt, kt, sigma_ct, sigma_nt, sigma_kt, &k, &kp, &sigma_k, &sigma_kp, fit_start);
Q = 0; /* quality_of_fit(chi2, 2);*/
- ddg = BOLTZ * temp * sigma_k / k;
+ ddg = gmx::c_boltz * temp * sigma_k / k;
printf("Fitting paramaters chi^2 = %10g, Quality of fit = %10g\n", chi2, Q);
printf("The Rate and Delta G are followed by an error estimate\n");
printf("----------------------------------------------------------\n"
1 / k,
calc_dg(1 / k, temp),
ddg);
- ddg = BOLTZ * temp * sigma_kp / kp;
+ ddg = gmx::c_boltz * temp * sigma_kp / kp;
printf("Backward %10.3f %6.2f %8.3f %10.3f %6.2f\n",
kp,
sigma_kp,
/* process input */
bSelected = FALSE;
- ccut = std::cos(acut * DEG2RAD);
+ ccut = std::cos(acut * gmx::c_deg2Rad);
if (bContact)
{
"for an hbond: %f",
dist);
}
- ang *= RAD2DEG;
+ ang *= gmx::c_rad2Deg;
__ADIST[static_cast<int>(ang / abin)]++;
__RDIST[static_cast<int>(dist / rbin)]++;
if (!bTwo)
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
static double cv_corr(double nu, double T)
{
- double x = PLANCK * nu / (BOLTZ * T);
+ double x = gmx::c_planck * nu / (gmx::c_boltz * T);
double ex = std::exp(x);
if (nu <= 0)
{
- return BOLTZ * KILO;
+ return gmx::c_boltz * gmx::c_kilo;
}
else
{
- return BOLTZ * KILO * (ex * gmx::square(x) / gmx::square(ex - 1) - 1);
+ return gmx::c_boltz * gmx::c_kilo * (ex * gmx::square(x) / gmx::square(ex - 1) - 1);
}
}
static double u_corr(double nu, double T)
{
- double x = PLANCK * nu / (BOLTZ * T);
+ double x = gmx::c_planck * nu / (gmx::c_boltz * T);
double ex = std::exp(x);
if (nu <= 0)
{
- return BOLTZ * T;
+ return gmx::c_boltz * T;
}
else
{
- return BOLTZ * T * (0.5 * x - 1 + x / (ex - 1));
+ return gmx::c_boltz * T * (0.5 * x - 1 + x / (ex - 1));
}
}
*/
static double calcTranslationalHeatCapacity()
{
- return RGAS * 1.5;
+ return gmx::c_universalGasConstant * 1.5;
}
/*! \brief Compute internal energy due to translational motion
*/
static double calcTranslationalInternalEnergy(double T)
{
- return BOLTZ * T * 1.5;
+ return gmx::c_boltz * T * 1.5;
}
/*! \brief Compute heat capacity due to rotational motion
{
if (linear)
{
- return BOLTZ * T;
+ return gmx::c_boltz * T;
}
else
{
- return BOLTZ * T * 1.5;
+ return gmx::c_boltz * T * 1.5;
}
}
{
if (linear)
{
- return RGAS;
+ return gmx::c_universalGasConstant;
}
else
{
- return RGAS * 1.5;
+ return gmx::c_universalGasConstant * 1.5;
}
}
principal_comp(index.size(), index.data(), top.atoms.atom, as_rvec_array(x_com.data()), trans, inertia);
bool linear = (inertia[XX] / inertia[YY] < linear_toler && inertia[XX] / inertia[ZZ] < linear_toler);
// (kJ/mol ps)^2/(Dalton nm^2 kJ/mol K) =
- // KILO kg m^2 ps^2/(s^2 mol g/mol nm^2 K) =
- // KILO^2 10^18 / 10^24 K = 1/K
- double rot_const = gmx::square(PLANCK) / (8 * gmx::square(M_PI) * BOLTZ);
+ // c_kilo kg m^2 ps^2/(s^2 mol g/mol nm^2 K) =
+ // c_kilo^2 10^18 / 10^24 K = 1/K
+ double rot_const = gmx::square(gmx::c_planck) / (8 * gmx::square(M_PI) * gmx::c_boltz);
// Rotational temperature (1/K)
rvec theta = { 0, 0, 0 };
if (linear)
* light. Do this by first converting to omega^2 (units 1/s), take the square
* root, and finally divide by the speed of light (nm/ps in gromacs).
*/
- factor_gmx_to_omega2 = 1.0E21 / (AVOGADRO * AMU);
- factor_omega_to_wavenumber = 1.0E-5 / (2.0 * M_PI * SPEED_OF_LIGHT);
+ factor_gmx_to_omega2 = 1.0E21 / (gmx::c_avogadro * gmx::c_amu);
+ factor_omega_to_wavenumber = 1.0E-5 / (2.0 * M_PI * gmx::c_speedOfLight);
value = 0;
for (i = begin; (i <= end); i++)
qu = u_corr(nu, T);
if (i > end - nharm)
{
- qcv += BOLTZ * KILO;
- qu += BOLTZ * T;
+ qcv += gmx::c_boltz * gmx::c_kilo;
+ qu += gmx::c_boltz * T;
}
fprintf(qc, "%6d %15g %15g\n", i, qcv, qu);
qcvtot += qcv;
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
/* (r-0.5) n times: var_n = n * var_1 = n/12
n=4: var_n = 1/3, so multiply with 3 */
- rfac = std::sqrt(3.0 * BOLTZ * temp / eigval[iout[j]]);
+ rfac = std::sqrt(3.0 * gmx::c_boltz * temp / eigval[iout[j]]);
rhalf = 2.0 * rfac;
rfac = rfac / im;
*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2017,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2017,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
/* Derive amplitude from temperature and eigenvalue if we can */
/* Convert eigenvalue to angular frequency, in units s^(-1) */
- omega = std::sqrt(eigval[kmode] * 1.0E21 / (AVOGADRO * AMU));
+ omega = std::sqrt(eigval[kmode] * 1.0E21 / (gmx::c_avogadro * gmx::c_amu));
/* Harmonic motion will be x=x0 + A*sin(omega*t)*eigenvec.
* The velocity is thus:
*
/* Convert Ekin from amu*(nm/s)^2 to J, i.e., kg*(m/s)^2
* This will also be proportional to A^2
*/
- Ekin *= AMU * 1E-18;
+ Ekin *= gmx::c_amu * 1E-18;
/* Set the amplitude so the energy is kT/2 */
- amplitude[i] = std::sqrt(0.5 * BOLTZMANN * temp / Ekin);
+ amplitude[i] = std::sqrt(0.5 * gmx::c_boltzmann * temp / Ekin);
}
else
{
*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2017,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2017,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
for (i = 0; (i < xr->npp); i++)
{
- phi = xr->dih[xr->pp[i].iphi].ang * RAD2DEG;
- psi = xr->dih[xr->pp[i].ipsi].ang * RAD2DEG;
+ phi = xr->dih[xr->pp[i].iphi].ang * gmx::c_rad2Deg;
+ psi = xr->dih[xr->pp[i].ipsi].ang * gmx::c_rad2Deg;
fprintf(out, "%g %g %s\n", phi, psi, xr->pp[i].label);
}
}
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#include "gromacs/gmxana/gstat.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/topology/topology.h"
#include "gromacs/utility/arraysize.h"
bfac[i] = ibox[i] / (max_eig[i] - min_eig[i]);
}
/* Do the binning */
- bref = 1 / (BOLTZ * Tref);
+ bref = 1 / (gmx::c_boltz * Tref);
snew(bE, n);
if (bGE || nenerT == 2)
{
}
else
{
- bE[j] = (bref - 1 / (BOLTZ * enerT[1][j])) * enerT[0][j];
+ bE[j] = (bref - 1 / (gmx::c_boltz * enerT[1][j])) * enerT[0][j];
}
Emin = std::min(Emin, static_cast<double>(bE[j]));
}
}
else if (idim[i] == -1)
{
- efac /= std::sin(DEG2RAD * eig[i][j]);
+ efac /= std::sin(gmx::c_deg2Rad * eig[i][j]);
}
}
/* Update the probability */
if (P[i] != 0)
{
Pmax = std::max(P[i], Pmax);
- W[i] = -BOLTZ * Tref * std::log(P[i]);
+ W[i] = -gmx::c_boltz * Tref * std::log(P[i]);
if (W[i] < Wmin)
{
Wmin = W[i];
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2012,2013,2014,2015,2017 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
static real pot(real x, real qq, real c6, real cn, int npow)
{
- return cn * pow(x, -npow) - c6 / gmx::power6(x) + qq * ONE_4PI_EPS0 / x;
+ return cn * pow(x, -npow) - c6 / gmx::power6(x) + qq * gmx::c_one4PiEps0 / x;
}
static real bhpot(real x, real A, real B, real C)
static real dpot(real x, real qq, real c6, real cn, int npow)
{
return -(npow * cn * std::pow(x, -npow - 1) - 6 * c6 / (x * gmx::power6(x))
- + qq * ONE_4PI_EPS0 / gmx::square(x));
+ + qq * gmx::c_one4PiEps0 / gmx::square(x));
}
int gmx_sigeps(int argc, char* argv[])
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
}
unitv(dir, dir);
- svmul(ENM2DEBYE, dip, dip);
+ svmul(gmx::c_enm2Debye, dip, dip);
dip2 = norm2(dip);
sdip += std::sqrt(dip2);
sdip2 += dip2;
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
fitparms[0] = 1;
fitparms[1] = 1;
do_lmfit(ncorr, tcaf[k], sig, dt, nullptr, 0, ncorr * dt, oenv, bDebugMode(), effnVAC, fitparms, 0, nullptr);
- eta = 1000 * fitparms[1] * rho / (4 * fitparms[0] * PICO * norm2(kfac[k]) / (NANO * NANO));
+ eta = 1000 * fitparms[1] * rho
+ / (4 * fitparms[0] * gmx::c_pico * norm2(kfac[k]) / (gmx::c_nano * gmx::c_nano));
fprintf(stdout, "k %6.3f tau %6.3f eta %8.5f 10^-3 kg/(m s)\n", norm(kfac[k]), fitparms[0], eta);
fprintf(fp_vk, "%6.3f %g\n", norm(kfac[k]), eta);
for (i = 0; i < ncorr; i++)
fitparms[1] = 1;
do_lmfit(ncorr, tcafc[k], sig, dt, nullptr, 0, ncorr * dt, oenv, bDebugMode(), effnVAC, fitparms, 0, nullptr);
eta = 1000 * fitparms[1] * rho
- / (4 * fitparms[0] * PICO * norm2(kfac[kset_c[k]]) / (NANO * NANO));
+ / (4 * fitparms[0] * gmx::c_pico * norm2(kfac[kset_c[k]]) / (gmx::c_nano * gmx::c_nano));
fprintf(stdout,
"k %6.3f tau %6.3f Omega %6.3f eta %8.5f 10^-3 kg/(m s)\n",
norm(kfac[kset_c[k]]),
dt = (t1 - t0) / (nframes - 1);
- rho *= sysmass / nframes * AMU / (NANO * NANO * NANO);
+ rho *= sysmass / nframes * gmx::c_amu / (gmx::c_nano * gmx::c_nano * gmx::c_nano);
fprintf(stdout, "Density = %g (kg/m^3)\n", rho);
process_tcaf(nframes,
dt,
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
ekin2 += mass[j] * norm2(v[j]);
}
- return ekin2 / (3 * isize * BOLTZ);
+ return ekin2 / (3 * isize * gmx::c_boltz);
}
static void remove_jump(matrix box, int natoms, rvec xp[], rvec x[])
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
* The timestep between saving the trajectory is
* 1e7 is to convert nanometer to cm
*/
- recip_fac = bRecip ? (1e7 / SPEED_OF_LIGHT) : 1.0;
+ recip_fac = bRecip ? (1e7 / gmx::c_speedOfLight) : 1.0;
for (i = 0; (i < n); i += 2)
{
nu = i / (2 * dt);
}
}
/* Note: there are two contributions to bin k in the wham equations:
- i) N[j]*exp(- U/(BOLTZ*opt->Temperature) + window[i].z[j])
- ii) exp(- U/(BOLTZ*opt->Temperature))
+ i) N[j]*exp(- U/(c_boltz*opt->Temperature) + window[i].z[j])
+ ii) exp(- U/(c_boltz*opt->Temperature))
where U is the umbrella potential
If any of these number is larger wham_contrib_lim, I set contrib=TRUE
*/
{
U = tabulated_pot(distance, opt); /* Use tabulated potential */
}
- contrib1 = profile[k] * std::exp(-U / (BOLTZ * opt->Temperature));
- contrib2 = window[i].N[j] * std::exp(-U / (BOLTZ * opt->Temperature) + window[i].z[j]);
+ contrib1 = profile[k] * std::exp(-U / (gmx::c_boltz * opt->Temperature));
+ contrib2 = window[i].N[j]
+ * std::exp(-U / (gmx::c_boltz * opt->Temperature) + window[i].z[j]);
window[i].bContrib[j][k] = (contrib1 > wham_contrib_lim || contrib2 > wham_contrib_lim);
bAnyContrib = bAnyContrib || window[i].bContrib[j][k];
if (window[i].bContrib[j][k])
U = tabulated_pot(distance, opt); /* Use tabulated potential */
}
denom += invg * window[j].N[k]
- * std::exp(-U / (BOLTZ * opt->Temperature) + window[j].z[k]);
+ * std::exp(-U / (gmx::c_boltz * opt->Temperature) + window[j].z[k]);
}
}
profile[i] = num / denom;
{
U = tabulated_pot(distance, opt); /* Use tabulated potential */
}
- total += profile[k] * std::exp(-U / (BOLTZ * opt->Temperature));
+ total += profile[k] * std::exp(-U / (gmx::c_boltz * opt->Temperature));
}
/* Avoid floating point exception if window is far outside min and max */
if (total != 0.0)
}
else if (opt->unit == en_kJ)
{
- unit_factor = BOLTZ * opt->Temperature;
+ unit_factor = gmx::c_boltz * opt->Temperature;
}
else if (opt->unit == en_kCal)
{
- unit_factor = (BOLTZ / CAL2JOULE) * opt->Temperature;
+ unit_factor = (gmx::c_boltz / gmx::c_cal2Joule) * opt->Temperature;
}
else
{
*/
for (j = 0; j < opt->bins; ++j)
{
- pot[j] = std::exp(-pot[j] / (BOLTZ * opt->Temperature));
+ pot[j] = std::exp(-pot[j] / (gmx::c_boltz * opt->Temperature));
}
calc_z(pot, window, nWindows, opt, TRUE);
xx = cp * x2[XX] + sp * x2[YY];
yy = -sp * x2[XX] + cp * x2[YY];
- dphi = RAD2DEG * std::atan2(yy, xx);
+ dphi = gmx::c_rad2Deg * std::atan2(yy, xx);
return dphi;
}
aj = index[i + 1];
ak = index[i + 2];
al = index[i + 3];
- phi = RAD2DEG
+ phi = gmx::c_rad2Deg
* dih_angle(x[ai], x[aj], x[ak], x[al], nullptr, r_ij, r_kj, r_kl, m, n, &t1, &t2, &t3);
phitot += phi;
}
bb[i].d5 = norm(dx);
}
- bb[i].phi = RAD2DEG
+ bb[i].phi = gmx::c_rad2Deg
* dih_angle(x[bb[i].Cprev],
x[bb[i].N],
x[bb[i].CA],
&t1,
&t2,
&t3);
- bb[i].psi = RAD2DEG
+ bb[i].psi = gmx::c_rad2Deg
* dih_angle(x[bb[i].N],
x[bb[i].CA],
x[bb[i].C],
&t3);
bb[i].pprms2 = gmx::square(bb[i].phi - PHI_AHX) + gmx::square(bb[i].psi - PSI_AHX);
- bb[i].jcaha += 1.4 * std::sin((bb[i].psi + 138.0) * DEG2RAD)
- - 4.1 * std::cos(2.0 * DEG2RAD * (bb[i].psi + 138.0))
- + 2.0 * std::cos(2.0 * DEG2RAD * (bb[i].phi + 30.0));
+ bb[i].jcaha += 1.4 * std::sin((bb[i].psi + 138.0) * gmx::c_deg2Rad)
+ - 4.1 * std::cos(2.0 * gmx::c_deg2Rad * (bb[i].psi + 138.0))
+ + 2.0 * std::cos(2.0 * gmx::c_deg2Rad * (bb[i].phi + 30.0));
}
}
*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
{
snew(newdata[i], nny);
phi = i * 2 * M_PI / (nnx - 1);
- x_phi[i] = phi * RAD2DEG - 180;
+ x_phi[i] = phi * gmx::c_rad2Deg - 180;
for (j = 0; (j < nny); j++)
{
psi = j * 2 * M_PI / (nny - 1);
if (i == 0)
{
- y_psi[j] = psi * RAD2DEG - 180;
+ y_psi[j] = psi * gmx::c_rad2Deg - 180;
}
/*if (((i % nfac) == 0) && ((j % nfac) == 0))
newdata[i][j] = sd->data[i/nfac][j/nfac];
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2017,2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
static double eigval_to_frequency(double eigval)
{
- double factor_gmx_to_omega2 = 1.0E21 / (AVOGADRO * AMU);
+ double factor_gmx_to_omega2 = 1.0E21 / (gmx::c_avogadro * gmx::c_amu);
return std::sqrt(std::max(0.0, eigval) * factor_gmx_to_omega2);
}
double calcZeroPointEnergy(gmx::ArrayRef<const real> eigval, real scale_factor)
{
// Convert frequency (ps^-1) to energy (kJ/mol)
- double factor = PLANCK * PICO / (2.0 * M_PI);
+ double factor = gmx::c_planck * gmx::c_pico / (2.0 * M_PI);
double zpe = 0;
for (auto& r : eigval)
{
{
size_t nskip = linear ? 5 : 6;
double Evib = 0;
- double hbar = PLANCK1 / (2 * M_PI);
+ double hbar = gmx::c_planck1 / (2 * M_PI);
for (gmx::index i = nskip; i < eigval.ssize(); i++)
{
if (eigval[i] > 0)
{
double omega = scale_factor * eigval_to_frequency(eigval[i]);
- double hwkT = (hbar * omega) / (BOLTZMANN * temperature);
+ double hwkT = (hbar * omega) / (gmx::c_boltzmann * temperature);
// Prevent overflow by checking for unreasonably large numbers.
if (hwkT < 100)
{
}
}
}
- return temperature * BOLTZ * Evib;
+ return temperature * gmx::c_boltz * Evib;
}
double calcVibrationalHeatCapacity(gmx::ArrayRef<const real> eigval, real temperature, gmx_bool linear, real scale_factor)
{
size_t nskip = linear ? 5 : 6;
double cv = 0;
- double hbar = PLANCK1 / (2 * M_PI);
+ double hbar = gmx::c_planck1 / (2 * M_PI);
for (gmx::index i = nskip; i < eigval.ssize(); i++)
{
if (eigval[i] > 0)
{
double omega = scale_factor * eigval_to_frequency(eigval[i]);
- double hwkT = (hbar * omega) / (BOLTZMANN * temperature);
+ double hwkT = (hbar * omega) / (gmx::c_boltzmann * temperature);
// Prevent overflow by checking for unreasonably large numbers.
if (hwkT < 100)
{
}
}
}
- return RGAS * cv;
+ return gmx::c_universalGasConstant * cv;
}
double calcTranslationalEntropy(real mass, real temperature, real pressure)
{
- double kT = BOLTZ * temperature;
+ double kT = gmx::c_boltz * temperature;
GMX_RELEASE_ASSERT(mass > 0, "Molecular mass should be larger than zero");
GMX_RELEASE_ASSERT(pressure > 0, "Pressure should be larger than zero");
GMX_RELEASE_ASSERT(temperature > 0, "Temperature should be larger than zero");
// Convert bar to Pascal
- double P = pressure * 1e5;
- double qT = (std::pow(2 * M_PI * mass * kT / gmx::square(PLANCK), 1.5) * (kT / P) * (1e30 / AVOGADRO));
- return RGAS * (std::log(qT) + 2.5);
+ double P = pressure * 1e5;
+ double qT = (std::pow(2 * M_PI * mass * kT / gmx::square(gmx::c_planck), 1.5) * (kT / P)
+ * (1e30 / gmx::c_avogadro));
+ return gmx::c_universalGasConstant * (std::log(qT) + 2.5);
}
double calcRotationalEntropy(real temperature, int natom, gmx_bool linear, const rvec theta, real sigma_r)
{
GMX_RELEASE_ASSERT(theta[0] > 0, "Theta should be larger than zero");
double qR = temperature / (sigma_r * theta[0]);
- sR = RGAS * (std::log(qR) + 1);
+ sR = gmx::c_universalGasConstant * (std::log(qR) + 1);
}
else
{
double Q = theta[XX] * theta[YY] * theta[ZZ];
GMX_RELEASE_ASSERT(Q > 0, "Q should be larger than zero");
double qR = std::sqrt(M_PI * std::pow(temperature, 3) / Q) / sigma_r;
- sR = RGAS * (std::log(qR) + 1.5);
+ sR = gmx::c_universalGasConstant * (std::log(qR) + 1.5);
}
}
return sR;
{
size_t nskip = bLinear ? 5 : 6;
double S = 0;
- double hbar = PLANCK1 / (2 * M_PI);
+ double hbar = gmx::c_planck1 / (2 * M_PI);
for (gmx::index i = nskip; (i < eigval.ssize()); i++)
{
if (eigval[i] > 0)
{
double omega = scale_factor * eigval_to_frequency(eigval[i]);
- double hwkT = (hbar * omega) / (BOLTZMANN * temperature);
+ double hwkT = (hbar * omega) / (gmx::c_boltzmann * temperature);
double dS = (hwkT / std::expm1(hwkT) - std::log1p(-std::exp(-hwkT)));
S += dS;
if (debug)
fprintf(debug, "eigval[%d] = %g\n", static_cast<int>(i + 1), static_cast<double>(eigval[i]));
}
}
- return S * RGAS;
+ return S * gmx::c_universalGasConstant;
}
double calcSchlitterEntropy(gmx::ArrayRef<const real> eigval, real temperature, gmx_bool bLinear)
{
size_t nskip = bLinear ? 5 : 6;
- double hbar = PLANCK1 / (2 * M_PI); // J s
- double kt = BOLTZMANN * temperature; // J
+ double hbar = gmx::c_planck1 / (2 * M_PI); // J s
+ double kt = gmx::c_boltzmann * temperature; // J
double kteh = kt * std::exp(2.0) / (hbar * hbar); // 1/(J s^2) = 1/(kg m^2)
- double evcorr = NANO * NANO * AMU;
+ double evcorr = gmx::c_nano * gmx::c_nano * gmx::c_amu;
if (debug)
{
fprintf(debug, "n = %td, kteh = %g evcorr = %g\n", ssize(eigval), kteh, evcorr);
double dd = 1 + kteh * eigval[i] * evcorr;
deter += std::log(dd);
}
- return 0.5 * RGAS * deter;
+ return 0.5 * gmx::c_universalGasConstant * deter;
}
*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2010,2014,2015,2016,2019, by the GROMACS development team, led by
+ * Copyright (c) 2010,2014,2015,2016,2019,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#define alfaHpl (2 * M_PI / 3) /* 120 degrees */
#define distH 0.1
-#define alfaCOM (DEG2RAD * 117)
-#define alfaCO (DEG2RAD * 121)
-#define alfaCOA (DEG2RAD * 115)
+#define alfaCOM (gmx::c_deg2Rad * 117)
+#define alfaCO (gmx::c_deg2Rad * 121)
+#define alfaCOA (gmx::c_deg2Rad * 115)
#define distO 0.123
#define distOA 0.125
{
case F_G96ANGLES:
/* Post processing of input data: store cosine iso angle itself */
- newparam->harmonic.rA = cos(old[0] * DEG2RAD);
+ newparam->harmonic.rA = cos(old[0] * gmx::c_deg2Rad);
newparam->harmonic.krA = old[1];
- newparam->harmonic.rB = cos(old[2] * DEG2RAD);
+ newparam->harmonic.rB = cos(old[2] * gmx::c_deg2Rad);
newparam->harmonic.krB = old[3];
break;
case F_G96BONDS:
newparam->vsite.f = old[5];
break;
case F_VSITE3FAD:
- newparam->vsite.a = old[1] * cos(DEG2RAD * old[0]);
- newparam->vsite.b = old[1] * sin(DEG2RAD * old[0]);
+ newparam->vsite.a = old[1] * cos(gmx::c_deg2Rad * old[0]);
+ newparam->vsite.b = old[1] * sin(gmx::c_deg2Rad * old[0]);
newparam->vsite.c = old[2];
newparam->vsite.d = old[3];
newparam->vsite.e = old[4];
printf(" center :%7.3f%7.3f%7.3f (nm)\n", gc[XX], gc[YY], gc[ZZ]);
printf(" box vectors :%7.3f%7.3f%7.3f (nm)\n", norm(box[XX]), norm(box[YY]), norm(box[ZZ]));
printf(" box angles :%7.2f%7.2f%7.2f (degrees)\n",
- norm2(box[ZZ]) == 0 ? 0 : RAD2DEG * gmx_angle(box[YY], box[ZZ]),
- norm2(box[ZZ]) == 0 ? 0 : RAD2DEG * gmx_angle(box[XX], box[ZZ]),
- norm2(box[YY]) == 0 ? 0 : RAD2DEG * gmx_angle(box[XX], box[YY]));
+ norm2(box[ZZ]) == 0 ? 0 : gmx::c_rad2Deg * gmx_angle(box[YY], box[ZZ]),
+ norm2(box[ZZ]) == 0 ? 0 : gmx::c_rad2Deg * gmx_angle(box[XX], box[ZZ]),
+ norm2(box[YY]) == 0 ? 0 : gmx::c_rad2Deg * gmx_angle(box[XX], box[YY]));
printf(" box volume :%7.2f (nm^3)\n", det(box));
}
real vol, dens;
vol = det(box);
- dens = (mass * AMU) / (vol * NANO * NANO * NANO);
+ dens = (mass * gmx::c_amu) / (vol * gmx::c_nano * gmx::c_nano * gmx::c_nano);
fprintf(stderr, "Volume of input %g (nm^3)\n", vol);
fprintf(stderr, "Mass of input %g (a.m.u.)\n", mass);
fprintf(stderr, "Density of input %g (g/l)\n", dens);
rotangles[ZZ]);
for (i = 0; i < DIM; i++)
{
- rotangles[i] *= DEG2RAD;
+ rotangles[i] *= gmx::c_deg2Rad;
}
rotate_conf(natom, x, v, rotangles[XX], rotangles[YY], rotangles[ZZ]);
}
{
printf("new box vectors :%7.3f%7.3f%7.3f (nm)\n", norm(box[XX]), norm(box[YY]), norm(box[ZZ]));
printf("new box angles :%7.2f%7.2f%7.2f (degrees)\n",
- norm2(box[ZZ]) == 0 ? 0 : RAD2DEG * gmx_angle(box[YY], box[ZZ]),
- norm2(box[ZZ]) == 0 ? 0 : RAD2DEG * gmx_angle(box[XX], box[ZZ]),
- norm2(box[YY]) == 0 ? 0 : RAD2DEG * gmx_angle(box[XX], box[YY]));
+ norm2(box[ZZ]) == 0 ? 0 : gmx::c_rad2Deg * gmx_angle(box[YY], box[ZZ]),
+ norm2(box[ZZ]) == 0 ? 0 : gmx::c_rad2Deg * gmx_angle(box[XX], box[ZZ]),
+ norm2(box[YY]) == 0 ? 0 : gmx::c_rad2Deg * gmx_angle(box[XX], box[YY]));
printf("new box volume :%7.2f (nm^3)\n", det(box));
}
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
static void low_mspeed(real tempi, gmx_mtop_t* mtop, rvec v[], gmx::ThreeFry2x64<>* rng, const gmx::MDLogger& logger)
{
int nrdf;
- real boltz;
real ekin, temp;
gmx::TabulatedNormalDistribution<real> normalDist;
- boltz = BOLTZ * tempi;
- ekin = 0.0;
- nrdf = 0;
+ ekin = 0.0;
+ nrdf = 0;
for (const AtomProxy atomP : AtomRange(*mtop))
{
const t_atom& local = atomP.atom();
if (mass > 0)
{
rng->restart(i, 0);
- real sd = std::sqrt(boltz / mass);
+ real sd = std::sqrt(gmx::c_boltz * tempi / mass);
for (int m = 0; (m < DIM); m++)
{
v[i][m] = sd * normalDist(*rng);
nrdf += DIM;
}
}
- temp = (2.0 * ekin) / (nrdf * BOLTZ);
+ temp = (2.0 * ekin) / (nrdf * gmx::c_boltz);
if (temp > 0)
{
real scal = std::sqrt(tempi / temp);
} /* for i */
}
-#define ANGLE_6RING (DEG2RAD * 120)
+#define ANGLE_6RING (gmx::c_deg2Rad * 120)
/* cosine rule: a^2 = b^2 + c^2 - 2 b c cos(alpha) */
/* get a^2 when a, b and alpha are given: */
b_CE3_HE3 = get_ddb_bond(vsitetop, "TRP", "CE3", "HE3");
b_CZ3_HZ3 = get_ddb_bond(vsitetop, "TRP", "CZ3", "HZ3");
- a_NE1_CE2_CD2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "NE1", "CE2", "CD2");
- a_CE2_CD2_CG = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CE2", "CD2", "CG");
- a_CB_CG_CD2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CB", "CG", "CD2");
- a_CD2_CG_CD1 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CD2", "CG", "CD1");
-
- a_CE2_CD2_CE3 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CE2", "CD2", "CE3");
- a_CD2_CE2_CZ2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CD2", "CE2", "CZ2");
- a_CD2_CE3_CZ3 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CD2", "CE3", "CZ3");
- a_CE3_CZ3_HZ3 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CE3", "CZ3", "HZ3");
- a_CZ2_CH2_HH2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CZ2", "CH2", "HH2");
- a_CE2_CZ2_HZ2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CE2", "CZ2", "HZ2");
- a_CE2_CZ2_CH2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CE2", "CZ2", "CH2");
- a_CG_CD1_HD1 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CG", "CD1", "HD1");
- a_HE1_NE1_CE2 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "HE1", "NE1", "CE2");
- a_CD2_CE3_HE3 = DEG2RAD * get_ddb_angle(vsitetop, "TRP", "CD2", "CE3", "HE3");
+ a_NE1_CE2_CD2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "NE1", "CE2", "CD2");
+ a_CE2_CD2_CG = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CE2", "CD2", "CG");
+ a_CB_CG_CD2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CB", "CG", "CD2");
+ a_CD2_CG_CD1 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CD2", "CG", "CD1");
+
+ a_CE2_CD2_CE3 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CE2", "CD2", "CE3");
+ a_CD2_CE2_CZ2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CD2", "CE2", "CZ2");
+ a_CD2_CE3_CZ3 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CD2", "CE3", "CZ3");
+ a_CE3_CZ3_HZ3 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CE3", "CZ3", "HZ3");
+ a_CZ2_CH2_HH2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CZ2", "CH2", "HH2");
+ a_CE2_CZ2_HZ2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CE2", "CZ2", "HZ2");
+ a_CE2_CZ2_CH2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CE2", "CZ2", "CH2");
+ a_CG_CD1_HD1 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CG", "CD1", "HD1");
+ a_HE1_NE1_CE2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "HE1", "NE1", "CE2");
+ a_CD2_CE3_HE3 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TRP", "CD2", "CE3", "HE3");
/* Calculate local coordinates.
* y-axis (x=0) is the bond CD2-CE2.
bond_ch = get_ddb_bond(vsitetop, "TYR", "CD1", "HD1");
bond_co = get_ddb_bond(vsitetop, "TYR", "CZ", "OH");
bond_oh = get_ddb_bond(vsitetop, "TYR", "OH", "HH");
- angle_coh = DEG2RAD * get_ddb_angle(vsitetop, "TYR", "CZ", "OH", "HH");
+ angle_coh = gmx::c_deg2Rad * get_ddb_angle(vsitetop, "TYR", "CZ", "OH", "HH");
xi[atCG] = -bond_cc + bond_cc * std::cos(ANGLE_6RING);
xi[atCD1] = -bond_cc;
b_CE1_NE2 = get_ddb_bond(vsitetop, resname, "CE1", "NE2");
b_CG_CD2 = get_ddb_bond(vsitetop, resname, "CG", "CD2");
b_CD2_NE2 = get_ddb_bond(vsitetop, resname, "CD2", "NE2");
- a_CG_ND1_CE1 = DEG2RAD * get_ddb_angle(vsitetop, resname, "CG", "ND1", "CE1");
- a_CG_CD2_NE2 = DEG2RAD * get_ddb_angle(vsitetop, resname, "CG", "CD2", "NE2");
- a_ND1_CE1_NE2 = DEG2RAD * get_ddb_angle(vsitetop, resname, "ND1", "CE1", "NE2");
- a_CE1_NE2_CD2 = DEG2RAD * get_ddb_angle(vsitetop, resname, "CE1", "NE2", "CD2");
+ a_CG_ND1_CE1 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "CG", "ND1", "CE1");
+ a_CG_CD2_NE2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "CG", "CD2", "NE2");
+ a_ND1_CE1_NE2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "ND1", "CE1", "NE2");
+ a_CE1_NE2_CD2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "CE1", "NE2", "CD2");
if (ats[atHD1] != NOTSET)
{
b_ND1_HD1 = get_ddb_bond(vsitetop, resname, "ND1", "HD1");
- a_CE1_ND1_HD1 = DEG2RAD * get_ddb_angle(vsitetop, resname, "CE1", "ND1", "HD1");
+ a_CE1_ND1_HD1 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "CE1", "ND1", "HD1");
}
if (ats[atHE2] != NOTSET)
{
b_NE2_HE2 = get_ddb_bond(vsitetop, resname, "NE2", "HE2");
- a_CE1_NE2_HE2 = DEG2RAD * get_ddb_angle(vsitetop, resname, "CE1", "NE2", "HE2");
+ a_CE1_NE2_HE2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "CE1", "NE2", "HE2");
}
if (ats[atHD2] != NOTSET)
{
b_CD2_HD2 = get_ddb_bond(vsitetop, resname, "CD2", "HD2");
- a_NE2_CD2_HD2 = DEG2RAD * get_ddb_angle(vsitetop, resname, "NE2", "CD2", "HD2");
+ a_NE2_CD2_HD2 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "NE2", "CD2", "HD2");
}
if (ats[atHE1] != NOTSET)
{
b_CE1_HE1 = get_ddb_bond(vsitetop, resname, "CE1", "HE1");
- a_NE2_CE1_HE1 = DEG2RAD * get_ddb_angle(vsitetop, resname, "NE2", "CE1", "HE1");
+ a_NE2_CE1_HE1 = gmx::c_deg2Rad * get_ddb_angle(vsitetop, resname, "NE2", "CE1", "HE1");
}
/* constraints between CG, CE1 and NE1 */
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
{
/* Compute number of ions to be added */
vol = det(box);
- nsalt = gmx::roundToInt(conc * vol * AVOGADRO / 1e24);
+ nsalt = gmx::roundToInt(conc * vol * gmx::c_avogadro / 1e24);
p_num = abs(nsalt * n_q);
n_num = abs(nsalt * p_q);
}
nrdf += ir->opts.nrdf[g];
}
- return sum_mv2 / (nrdf * BOLTZ);
+ return sum_mv2 / (nrdf * gmx::c_boltz);
}
static real get_max_reference_temp(const t_inputrec* ir, warninp* wi)
* to be on the safe side with constraints.
*/
const real totalEnergyDriftPerAtomPerPicosecond =
- 2 * BOLTZ * buffer_temp / (ir->nsteps * ir->delta_t);
+ 2 * gmx::c_boltz * buffer_temp / (ir->nsteps * ir->delta_t);
if (ir->verletbuf_tol > 1.1 * driftTolerance * totalEnergyDriftPerAtomPerPicosecond)
{
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
d2 = distance2(x[i], x[j]);
rvec_sub(x[i], xh, nh);
rvec_sub(x[aj], xh, oh);
- a = RAD2DEG * acos(cos_angle(nh, oh));
+ a = gmx::c_rad2Deg * acos(cos_angle(nh, oh));
if ((d2 < dist2) && (a > angle))
{
hbond[i] = TRUE;
#include "gromacs/gmxpreprocess/toputil.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/calc_verletbuf.h"
#include "gromacs/mdrun/mdmodules.h"
* of errors. The factor 0.5 is because energy distributes
* equally over Ekin and Epot.
*/
- max_T_error = 0.5 * tau * ir->verletbuf_tol / (nrdf_at * BOLTZ * T);
+ max_T_error = 0.5 * tau * ir->verletbuf_tol / (nrdf_at * gmx::c_boltz * T);
if (max_T_error > T_error_warn)
{
sprintf(warn_buf,
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
vol = det(box);
fprintf(stderr, "Volume : %10g (nm^3)\n", vol);
- fprintf(stderr, "Density : %10g (g/l)\n", (mtot * 1e24) / (AVOGADRO * vol));
+ fprintf(stderr, "Density : %10g (g/l)\n", (mtot * 1e24) / (gmx::c_avogadro * vol));
fprintf(stderr, "Number of solvent molecules: %5d \n\n", nsol);
/* open topology file and append sol molecules */
*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
}
if (bFound)
{
- real param = std::sqrt(b_ij * b_ij + b_jk * b_jk
- - 2.0 * b_ij * b_jk
- * cos(DEG2RAD * ang->c0()));
+ real param = std::sqrt(
+ b_ij * b_ij + b_jk * b_jk
+ - 2.0 * b_ij * b_jk * cos(gmx::c_deg2Rad * ang->c0()));
std::vector<real> forceParm = { param, param };
if (ftype == F_CONNBONDS || ftype_a == F_CONNBONDS)
{
if (((ai == ang.ai()) && (aj == ang.aj()) && (ak == ang.ak()))
|| ((ai == ang.ak()) && (aj == ang.aj()) && (ak == ang.ai())))
{
- angle = DEG2RAD * ang.parameterValue();
+ angle = gmx::c_deg2Rad * ang.parameterValue();
break;
}
}
bError = (bij == NOTSET) || (aijk == NOTSET);
vsite->setForceParameter(1, bij);
- vsite->setForceParameter(0, RAD2DEG * aijk);
+ vsite->setForceParameter(0, gmx::c_rad2Deg * aijk);
if (bSwapParity)
{
dH = bCN - bNH * std::cos(aCNH);
rH = bNH * std::sin(aCNH);
/* we assume the H's are symmetrically distributed */
- rHx = rH * std::cos(DEG2RAD * 30);
- rHy = rH * std::sin(DEG2RAD * 30);
+ rHx = rH * std::cos(gmx::c_deg2Rad * 30);
+ rHy = rH * std::sin(gmx::c_deg2Rad * 30);
rM = 0.5 * bMM;
dM = std::sqrt(gmx::square(bCM) - gmx::square(rM));
a = 0.5 * ((dH / dM) - (rHy / rM));
.appendTextFormatted(
"virtual site %d: angle ijk = %f, angle ijl = %f, angle ijm = %f",
vsite->ai() + 1,
- RAD2DEG * aijk,
- RAD2DEG * aijl,
- RAD2DEG * aijm);
+ gmx::c_rad2Deg * aijk,
+ gmx::c_rad2Deg * aijl,
+ gmx::c_rad2Deg * aijm);
gmx_fatal(FARGS,
"invalid construction in calc_vsite4fd for atom %d: "
"cosakl=%f, cosakm=%f\n",
.appendTextFormatted(
"virtual site %d: angle ijk = %f, angle ijl = %f, angle ijm = %f",
vsite->ai() + 1,
- RAD2DEG * aijk,
- RAD2DEG * aijl,
- RAD2DEG * aijm);
+ gmx::c_rad2Deg * aijk,
+ gmx::c_rad2Deg * aijl,
+ gmx::c_rad2Deg * aijm);
gmx_fatal(FARGS,
"invalid construction in calc_vsite4fdn for atom %d: "
"pl=%f, pm=%f\n",
int ai = angle.ai();
int aj = angle.aj();
int ak = angle.ak();
- real th = RAD2DEG
+ real th = gmx::c_rad2Deg
* bond_angle(x[ai], x[aj], x[ak], bPBC ? &pbc : nullptr, r_ij, r_kj, &costh, &t1, &t2);
angle.setForceParameter(0, th);
}
int ak = dihedral.ak();
int al = dihedral.al();
real ph =
- RAD2DEG
+ gmx::c_rad2Deg
* dih_angle(x[ai], x[aj], x[ak], x[al], bPBC ? &pbc : nullptr, r_ij, r_kj, r_kl, m, n, &t1, &t2, &t3);
dihedral.setForceParameter(0, ph);
}
fill_header(reinterpret_cast<IMDHeader*>(buffer), IMD_FCOORDS, static_cast<int32_t>(nat));
for (i = 0; i < nat; i++)
{
- sendx[0] = static_cast<float>(x[i][0]) * NM2A;
- sendx[1] = static_cast<float>(x[i][1]) * NM2A;
- sendx[2] = static_cast<float>(x[i][2]) * NM2A;
+ sendx[0] = static_cast<float>(x[i][0]) * gmx::c_nm2A;
+ sendx[1] = static_cast<float>(x[i][1]) * gmx::c_nm2A;
+ sendx[2] = static_cast<float>(x[i][2]) * gmx::c_nm2A;
memcpy(buffer + c_headerSize + i * tuplesize, sendx, tuplesize);
}
void ImdSession::Impl::copyToMDForces()
{
int i;
- real conversion = CAL2JOULE * NM2A;
+ real conversion = gmx::c_cal2Joule * gmx::c_nm2A;
for (i = 0; i < nforces; i++)
type = forceatoms[i++];
ai = forceatoms[i++];
aj = forceatoms[i++];
- ksh = gmx::square(md->chargeA[aj]) * ONE_4PI_EPS0 / forceparams[type].polarize.alpha;
+ ksh = gmx::square(md->chargeA[aj]) * gmx::c_one4PiEps0 / forceparams[type].polarize.alpha;
ki = pbc_rvec_sub(pbc, x[ai], x[aj], dx); /* 3 */
dr2 = iprod(dx, dx); /* 5 */
type = forceatoms[i++];
ai = forceatoms[i++];
aj = forceatoms[i++];
- ksh = gmx::square(md->chargeA[aj]) * ONE_4PI_EPS0 / forceparams[type].anharm_polarize.alpha; /* 7*/
+ ksh = gmx::square(md->chargeA[aj]) * gmx::c_one4PiEps0 / forceparams[type].anharm_polarize.alpha; /* 7*/
khyp = forceparams[type].anharm_polarize.khyp;
drcut = forceparams[type].anharm_polarize.drcut;
type0 = forceatoms[0];
aS = forceatoms[5];
qS = md->chargeA[aS];
- kk[XX] = gmx::square(qS) * ONE_4PI_EPS0 / forceparams[type0].wpol.al_x;
- kk[YY] = gmx::square(qS) * ONE_4PI_EPS0 / forceparams[type0].wpol.al_y;
- kk[ZZ] = gmx::square(qS) * ONE_4PI_EPS0 / forceparams[type0].wpol.al_z;
+ kk[XX] = gmx::square(qS) * gmx::c_one4PiEps0 / forceparams[type0].wpol.al_x;
+ kk[YY] = gmx::square(qS) * gmx::c_one4PiEps0 / forceparams[type0].wpol.al_y;
+ kk[ZZ] = gmx::square(qS) * gmx::c_one4PiEps0 / forceparams[type0].wpol.al_z;
r_HH = 1.0 / forceparams[type0].wpol.rHH;
for (i = 0; (i < nbonds); i += 6)
{
r12sq = iprod(r12, r12); /* 5 */
r12_1 = gmx::invsqrt(r12sq); /* 5 */
r12bar = afac / r12_1; /* 5 */
- v0 = qq * ONE_4PI_EPS0 * r12_1; /* 2 */
+ v0 = qq * gmx::c_one4PiEps0 * r12_1; /* 2 */
ebar = std::exp(-r12bar); /* 5 */
v1 = (1 - (1 + 0.5 * r12bar) * ebar); /* 4 */
fscal = ((v0 * r12_1) * v1 - v0 * 0.5 * afac * ebar * (r12bar + 1)) * r12_1; /* 9 */
*dvdlambda += harmonic(forceparams[type].harmonic.krA,
forceparams[type].harmonic.krB,
- forceparams[type].harmonic.rA * DEG2RAD,
- forceparams[type].harmonic.rB * DEG2RAD,
+ forceparams[type].harmonic.rA * gmx::c_deg2Rad,
+ forceparams[type].harmonic.rB * gmx::c_deg2Rad,
theta,
lambda,
&va,
alignas(GMX_SIMD_ALIGNMENT) std::int32_t aj[GMX_SIMD_REAL_WIDTH];
alignas(GMX_SIMD_ALIGNMENT) std::int32_t ak[GMX_SIMD_REAL_WIDTH];
alignas(GMX_SIMD_ALIGNMENT) real coeff[2 * GMX_SIMD_REAL_WIDTH];
- SimdReal deg2rad_S(DEG2RAD);
+ SimdReal deg2rad_S(gmx::c_deg2Rad);
SimdReal xi_S, yi_S, zi_S;
SimdReal xj_S, yj_S, zj_S;
SimdReal xk_S, yk_S, zk_S;
ai = forceatoms[i++];
aj = forceatoms[i++];
ak = forceatoms[i++];
- th0A = forceparams[type].u_b.thetaA * DEG2RAD;
+ th0A = forceparams[type].u_b.thetaA * gmx::c_deg2Rad;
kthA = forceparams[type].u_b.kthetaA;
r13A = forceparams[type].u_b.r13A;
kUBA = forceparams[type].u_b.kUBA;
- th0B = forceparams[type].u_b.thetaB * DEG2RAD;
+ th0B = forceparams[type].u_b.thetaB * gmx::c_deg2Rad;
kthB = forceparams[type].u_b.kthetaB;
r13B = forceparams[type].u_b.r13B;
kUBB = forceparams[type].u_b.kUBB;
SimdReal rikz_S = zi_S - zk_S;
const SimdReal ktheta_S = load<SimdReal>(coeff);
- const SimdReal theta0_S = load<SimdReal>(coeff + GMX_SIMD_REAL_WIDTH) * DEG2RAD;
+ const SimdReal theta0_S = load<SimdReal>(coeff + GMX_SIMD_REAL_WIDTH) * gmx::c_deg2Rad;
const SimdReal kUB_S = load<SimdReal>(coeff + 2 * GMX_SIMD_REAL_WIDTH);
const SimdReal r13_S = load<SimdReal>(coeff + 3 * GMX_SIMD_REAL_WIDTH);
theta = bond_angle(x[ai], x[aj], x[ak], pbc, r_ij, r_kj, &cos_theta, &t1, &t2); /* 41 */
- dt = theta - forceparams[type].qangle.theta * DEG2RAD; /* 2 */
+ dt = theta - forceparams[type].qangle.theta * gmx::c_deg2Rad; /* 2 */
dVdt = 0;
va = forceparams[type].qangle.c[0];
real dopdihs(real cpA, real cpB, real phiA, real phiB, int mult, real phi, real lambda, real* V, real* dvdlambda)
{
const real L1 = 1.0 - lambda;
- const real ph0 = (L1 * phiA + lambda * phiB) * DEG2RAD;
- const real dph0 = (phiB - phiA) * DEG2RAD;
+ const real ph0 = (L1 * phiA + lambda * phiB) * gmx::c_deg2Rad;
+ const real dph0 = (phiB - phiA) * gmx::c_deg2Rad;
const real cp = L1 * cpA + lambda * cpB;
const real mdphi = mult * phi - ph0;
{
real v, dvdlambda, mdphi, v1, sdphi, ddphi;
real L1 = 1.0 - lambda;
- real ph0 = (L1 * phiA + lambda * phiB) * DEG2RAD;
- real dph0 = (phiB - phiA) * DEG2RAD;
+ real ph0 = (L1 * phiA + lambda * phiB) * gmx::c_deg2Rad;
+ real dph0 = (phiB - phiA) * gmx::c_deg2Rad;
real cp = L1 * cpA + lambda * cpB;
mdphi = mult * (phi - ph0);
alignas(GMX_SIMD_ALIGNMENT) std::int32_t al[GMX_SIMD_REAL_WIDTH];
alignas(GMX_SIMD_ALIGNMENT) real buf[3 * GMX_SIMD_REAL_WIDTH];
real * cp, *phi0, *mult;
- SimdReal deg2rad_S(DEG2RAD);
+ SimdReal deg2rad_S(gmx::c_deg2Rad);
SimdReal p_S, q_S;
SimdReal phi0_S, phi_S;
SimdReal mx_S, my_S, mz_S;
pB = forceparams[type].harmonic.rB;
kk = L1 * kA + lambda * kB;
- phi0 = (L1 * pA + lambda * pB) * DEG2RAD;
- dphi0 = (pB - pA) * DEG2RAD;
+ phi0 = (L1 * pA + lambda * pB) * gmx::c_deg2Rad;
+ dphi0 = (pB - pA) * gmx::c_deg2Rad;
dp = phi - phi0;
L1 = 1.0 - lambda;
- d2r = DEG2RAD;
+ d2r = gmx::c_deg2Rad;
for (i = 0; (i < nbonds);)
{
/* Switch to degrees */
dx = 360.0 / cmap_grid->grid_spacing;
- xphi1 = xphi1 * RAD2DEG;
- xphi2 = xphi2 * RAD2DEG;
+ xphi1 = xphi1 * gmx::c_rad2Deg;
+ xphi2 = xphi2 * gmx::c_rad2Deg;
for (i = 0; i < 4; i++) /* 16 */
{
df2 = tt * df2 + (3.0 * tc[i * 4 + 3] * tu + 2.0 * tc[i * 4 + 2]) * tu + tc[i * 4 + 1];
}
- fac = RAD2DEG / dx;
+ fac = gmx::c_rad2Deg / dx;
df1 = df1 * fac;
df2 = df2 * fac;
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2014,2015,2017,2019, by the GROMACS development team, led by
+ * Copyright (c) 2014,2015,2017,2019,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
double term_theta_theta_equil;
k_bending = forceparams[type].harmonic.krA;
- theta_equil = forceparams[type].harmonic.rA * DEG2RAD;
+ theta_equil = forceparams[type].harmonic.rA * gmx::c_deg2Rad;
theta_equil = M_PI - theta_equil;
cosine_theta_equil = cos(theta_equil);
real norm_phi;
/* Read parameters phi0 and k_torsion */
- phi0 = forceparams[type].pdihs.phiA * DEG2RAD;
+ phi0 = forceparams[type].pdihs.phiA * gmx::c_deg2Rad;
cosine_phi0 = cos(phi0);
k_torsion = forceparams[type].pdihs.cpA;
* should be negligible (unless nstlist is extremely large, which
* you wouldn't do anyhow).
*/
- kT_fac = 2 * BOLTZ * temperature * timePeriod;
+ kT_fac = 2 * gmx::c_boltz * temperature * timePeriod;
if (ir.bd_fric > 0)
{
/* This is directly sigma^2 of the displacement */
}
else
{
- kT_fac = BOLTZ * temperature * gmx::square(timePeriod);
+ kT_fac = gmx::c_boltz * temperature * gmx::square(timePeriod);
}
return kT_fac;
"interactions");
}
- elfac = ONE_4PI_EPS0 / ir.epsilon_r;
+ elfac = gmx::c_one4PiEps0 / ir.epsilon_r;
// Determine the 1st and 2nd derivative for the electostatics
pot_derivatives_t elec = { 0, 0, 0 };
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2012,2014,2015,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
for (m = 0; (m < DIM); m++)
{
- mu[m] *= ENM2DEBYE;
+ mu[m] *= gmx::c_enm2Debye;
}
if (nChargePerturbed)
mu_y += qB[i] * x[i][YY];
mu_z += qB[i] * x[i][ZZ];
}
- mu_B[XX] = mu_x * ENM2DEBYE;
- mu_B[YY] = mu_y * ENM2DEBYE;
- mu_B[ZZ] = mu_z * ENM2DEBYE;
+ mu_B[XX] = mu_x * gmx::c_enm2Debye;
+ mu_B[YY] = mu_y * gmx::c_enm2Debye;
+ mu_B[ZZ] = mu_z * gmx::c_enm2Debye;
}
else
{
Ekin = 2 * trace(tcstat->ekinh) * tcstat->ekinscaleh_nhc;
}
}
- kT = BOLTZ * reft;
+ kT = gmx::c_boltz * reft;
for (mi = 0; mi < mstepsi; mi++)
{
pscal = calc_pres(ir->pbcType, nwall, box, ekinmod, vir, localpres) + pcorr;
vol = det(box);
- GW = (vol * (MassQ->Winv / PRESFAC)) * (DIM * pscal - trace(ir->ref_p)); /* W is in ps^2 * bar * nm^3 */
+ GW = (vol * (MassQ->Winv / gmx::c_presfac))
+ * (DIM * pscal - trace(ir->ref_p)); /* W is in ps^2 * bar * nm^3 */
*veta += 0.5 * dt * GW;
}
* het systeem...
*/
- fac = PRESFAC * 2.0 / det(box);
+ fac = gmx::c_presfac * 2.0 / det(box);
for (n = 0; (n < DIM); n++)
{
for (m = 0; (m < DIM); m++)
{
if (nrdf > 0)
{
- return (2.0 * ekin) / (nrdf * BOLTZ);
+ return (2.0 * ekin) / (nrdf * gmx::c_boltz);
}
else
{
}
}
-/*! \brief Sets 1/mass for Parrinello-Rahman in wInv; NOTE: PRESFAC is not included, so not in GROMACS units! */
+/*! \brief Sets 1/mass for Parrinello-Rahman in wInv; NOTE: c_presfac is not included, so not in GROMACS units! */
static void calcParrinelloRahmanInvMass(const t_inputrec* ir, const matrix box, tensor wInv)
{
real maxBoxLength;
if (!bFirstStep)
{
- /* Note that PRESFAC does not occur here.
+ /* Note that c_presfac does not occur here.
* The pressure and compressibility always occur as a product,
* therefore the pressure unit drops out.
*/
}
real gauss = 0;
real gauss2 = 0;
- real kt = ir->opts.ref_t[0] * BOLTZ;
+ real kt = ir->opts.ref_t[0] * gmx::c_boltz;
if (kt < 0.0)
{
kt = 0.0;
{
const real factor = compressibilityFactor(d, d, ir, dt);
mu[d][d] = std::exp(-factor * (ir->ref_p[d][d] - scalar_pressure) / DIM
- + std::sqrt(2.0 * kt * factor * PRESFAC / vol) * gauss / DIM);
+ + std::sqrt(2.0 * kt * factor * gmx::c_presfac / vol) * gauss / DIM);
}
break;
case PressureCouplingType::SemiIsotropic:
{
const real factor = compressibilityFactor(d, d, ir, dt);
mu[d][d] = std::exp(-factor * (ir->ref_p[d][d] - xy_pressure) / DIM
- + std::sqrt((DIM - 1) * 2.0 * kt * factor * PRESFAC / vol / DIM)
+ + std::sqrt((DIM - 1) * 2.0 * kt * factor * gmx::c_presfac / vol / DIM)
/ (DIM - 1) * gauss);
}
{
const real factor = compressibilityFactor(ZZ, ZZ, ir, dt);
mu[ZZ][ZZ] = std::exp(-factor * (ir->ref_p[ZZ][ZZ] - pres[ZZ][ZZ]) / DIM
- + std::sqrt(2.0 * kt * factor * PRESFAC / vol / DIM) * gauss2);
+ + std::sqrt(2.0 * kt * factor * gmx::c_presfac / vol / DIM) * gauss2);
}
break;
case PressureCouplingType::SurfaceTension:
/* Notice: we here use ref_p[ZZ][ZZ] as isotropic pressure and ir->ref_p[d][d] as surface tension */
mu[d][d] = std::exp(
-factor * (ir->ref_p[ZZ][ZZ] - ir->ref_p[d][d] / box[ZZ][ZZ] - xy_pressure) / DIM
- + std::sqrt(4.0 / 3.0 * kt * factor * PRESFAC / vol) / (DIM - 1) * gauss);
+ + std::sqrt(4.0 / 3.0 * kt * factor * gmx::c_presfac / vol) / (DIM - 1) * gauss);
}
{
const real factor = compressibilityFactor(ZZ, ZZ, ir, dt);
mu[ZZ][ZZ] = std::exp(-factor * (ir->ref_p[ZZ][ZZ] - pres[ZZ][ZZ]) / DIM
- + std::sqrt(2.0 / 3.0 * kt * factor * PRESFAC / vol) * gauss2);
+ + std::sqrt(2.0 / 3.0 * kt * factor * gmx::c_presfac / vol) * gauss2);
}
break;
default:
/* units are nm^3 * ns^2 / (nm^3 * bar / kJ/mol) = kJ/mol */
/* Consider evaluating eventually if this the right mass to use. All are correct, some might be more stable */
- MassQ->Winv = (PRESFAC * trace(ir->compress) * BOLTZ * opts->ref_t[0])
+ MassQ->Winv = (gmx::c_presfac * trace(ir->compress) * gmx::c_boltz * opts->ref_t[0])
/ (DIM * state->vol0 * gmx::square(ir->tau_p / M_2PI));
/* An alternate mass definition, from Tuckerman et al. */
- /* MassQ->Winv = 1.0/(gmx::square(ir->tau_p/M_2PI)*(opts->nrdf[0]+DIM)*BOLTZ*opts->ref_t[0]); */
+ /* MassQ->Winv = 1.0/(gmx::square(ir->tau_p/M_2PI)*(opts->nrdf[0]+DIM)*c_boltz*opts->ref_t[0]); */
for (d = 0; d < DIM; d++)
{
for (n = 0; n < DIM; n++)
{
- MassQ->Winvm[d][n] =
- PRESFAC * ir->compress[d][n] / (state->vol0 * gmx::square(ir->tau_p / M_2PI));
+ MassQ->Winvm[d][n] = gmx::c_presfac * ir->compress[d][n]
+ / (state->vol0 * gmx::square(ir->tau_p / M_2PI));
/* not clear this is correct yet for the anisotropic case. Will need to reevaluate
before using MTTK for anisotropic states.*/
}
{
reft = std::max<real>(0, opts->ref_t[i]);
nd = opts->nrdf[i];
- kT = BOLTZ * reft;
+ kT = gmx::c_boltz * reft;
for (j = 0; j < nh; j++)
{
if (j == 0)
if ((ir->tau_p > 0) && (opts->ref_t[0] > 0))
{
reft = std::max<real>(0, opts->ref_t[0]);
- kT = BOLTZ * reft;
+ kT = gmx::c_boltz * reft;
for (i = 0; i < nnhpres; i++)
{
for (j = 0; j < nh; j++)
real nd = ir->opts.nrdf[i];
real reft = std::max<real>(ir->opts.ref_t[i], 0);
- real kT = BOLTZ * reft;
+ real kT = gmx::c_boltz * reft;
if (nd > 0.0)
{
}
else /* Other non Trotter temperature NH control -- no chains yet. */
{
- energy += 0.5 * BOLTZ * nd * gmx::square(ivxi[0]) / iQinv[0];
+ energy += 0.5 * gmx::c_boltz * nd * gmx::square(ivxi[0]) / iQinv[0];
energy += nd * ixi[0] * kT;
}
}
{
/* note -- assumes only one degree of freedom that is thermostatted in barostat */
real reft = std::max<real>(ir->opts.ref_t[0], 0.0); /* using 'System' temperature */
- real kT = BOLTZ * reft;
+ real kT = gmx::c_boltz * reft;
for (int j = 0; j < nh; j++)
{
{
if (invMass[d][n] > 0)
{
- energyNPT += 0.5 * gmx::square(state->boxv[d][n]) / (invMass[d][n] * PRESFAC);
+ energyNPT += 0.5 * gmx::square(state->boxv[d][n])
+ / (invMass[d][n] * gmx::c_presfac);
}
}
}
* track of unwrapped box diagonal elements. This case is
* excluded in integratorHasConservedEnergyQuantity().
*/
- energyNPT += vol * trace(ir->ref_p) / (DIM * PRESFAC);
+ energyNPT += vol * trace(ir->ref_p) / (DIM * gmx::c_presfac);
break;
}
case PressureCoupling::Mttk:
energyNPT += 0.5 * gmx::square(state->veta) / MassQ->Winv;
/* contribution from the PV term */
- energyNPT += vol * trace(ir->ref_p) / (DIM * PRESFAC);
+ energyNPT += vol * trace(ir->ref_p) / (DIM * gmx::c_presfac);
if (ir->epc == PressureCoupling::Mttk)
{
if (opts->tau_t[i] >= 0 && opts->nrdf[i] > 0 && Ek > 0)
{
- Ek_ref1 = 0.5 * opts->ref_t[i] * BOLTZ;
+ Ek_ref1 = 0.5 * opts->ref_t[i] * gmx::c_boltz;
Ek_ref = Ek_ref1 * opts->nrdf[i];
Ek_new = vrescale_resamplekin(Ek, Ek_ref, opts->nrdf[i], opts->tau_t[i] / dt, step, ir->ld_seed);
corr.virial += numCorr * density * avctwelve * iParams_.virdifftwelve_ / 3.0;
}
/* The factor 2 is because of the Gromacs virial definition */
- corr.pressure = -2.0 * invvol * corr.virial * PRESFAC;
+ corr.pressure = -2.0 * invvol * corr.virial * gmx::c_presfac;
}
if (eFep_ != FreeEnergyPerturbationType::No)
nboxs = boxs_nm.size();
}
vol = box[XX][XX] * box[YY][YY] * box[ZZ][ZZ];
- dens = (tmass * AMU) / (vol * NANO * NANO * NANO);
+ dens = (tmass * gmx::c_amu) / (vol * gmx::c_nano * gmx::c_nano * gmx::c_nano);
add_ebin(ebin_, ib_, nboxs, bs, bSum);
add_ebin(ebin_, ivol_, 1, &vol, bSum);
add_ebin(ebin_, idens_, 1, &dens, bSum);
{
/* This is pV (in kJ/mol). The pressure is the reference pressure,
not the instantaneous pressure */
- pv = vol * ref_p_ / PRESFAC;
+ pv = vol * ref_p_ / gmx::c_presfac;
add_ebin(ebin_, ipv_, 1, &pv, bSum);
enthalpy = pv + enerd->term[F_ETOT];
if (ekind && ekind->cosacc.cos_accel != 0)
{
vol = box[XX][XX] * box[YY][YY] * box[ZZ][ZZ];
- dens = (tmass * AMU) / (vol * NANO * NANO * NANO);
+ dens = (tmass * gmx::c_amu) / (vol * gmx::c_nano * gmx::c_nano * gmx::c_nano);
add_ebin(ebin_, ivcos_, 1, &(ekind->cosacc.vcos), bSum);
/* 1/viscosity, unit 1/(kg m^-1 s^-1) */
tmp = 1
- / (ekind->cosacc.cos_accel / (ekind->cosacc.vcos * PICO) * dens
- * gmx::square(box[ZZ][ZZ] * NANO / (2 * M_PI)));
+ / (ekind->cosacc.cos_accel / (ekind->cosacc.vcos * gmx::c_pico) * dens
+ * gmx::square(box[ZZ][ZZ] * gmx::c_nano / (2 * M_PI)));
add_ebin(ebin_, ivisc_, 1, &tmp, bSum);
}
if (nE_ > 1)
#include "gromacs/listed_forces/orires.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/calcmu.h"
#include "gromacs/mdlib/constr.h"
if (ir->bSimTemp)
{
/* Note -- this assumes no mass changes, since kinetic energy is not added . . . */
- scaled_lamee[i] = enerd->foreignLambdaTerms.deltaH(i) / (simtemp->temperatures[i] * BOLTZ)
- + enerd->term[F_EPOT]
- * (1.0 / (simtemp->temperatures[i])
- - 1.0 / (simtemp->temperatures[fep_state]))
- / BOLTZ;
+ scaled_lamee[i] =
+ enerd->foreignLambdaTerms.deltaH(i) / (simtemp->temperatures[i] * gmx::c_boltz)
+ + enerd->term[F_EPOT]
+ * (1.0 / (simtemp->temperatures[i])
+ - 1.0 / (simtemp->temperatures[fep_state]))
+ / gmx::c_boltz;
}
else
{
- scaled_lamee[i] = enerd->foreignLambdaTerms.deltaH(i) / (expand->mc_temp * BOLTZ);
+ scaled_lamee[i] = enerd->foreignLambdaTerms.deltaH(i) / (expand->mc_temp * gmx::c_boltz);
/* mc_temp is currently set to the system reft unless otherwise defined */
}
{
scaled_lamee[i] =
enerd->term[F_EPOT]
- * (1.0 / simtemp->temperatures[i] - 1.0 / simtemp->temperatures[fep_state]) / BOLTZ;
+ * (1.0 / simtemp->temperatures[i] - 1.0 / simtemp->temperatures[fep_state])
+ / gmx::c_boltz;
}
}
}
/* Set the Coulomb energy conversion factor */
if (interactionConst.epsilon_r != 0)
{
- interactionConst.epsfac = ONE_4PI_EPS0 / interactionConst.epsilon_r;
+ interactionConst.epsfac = gmx::c_one4PiEps0 / interactionConst.epsilon_r;
}
else
{
real wfac;
- wfac = std::cos(DEG2RAD * wangle);
+ wfac = std::cos(gmx::c_deg2Rad * wangle);
wfac = wfac * wfac;
for (int iter = 0; iter < lincsd->nIter; iter++)
int maxwarn,
int* warncount)
{
- real wfac = std::cos(DEG2RAD * wangle);
+ real wfac = std::cos(gmx::c_deg2Rad * wangle);
fprintf(stderr,
"bonds that rotated more than %g degrees:\n"
" %6d %6d %5.1f %8.4f %8.4f %8.4f\n",
ddglatnr(dd, i),
ddglatnr(dd, j),
- RAD2DEG * std::acos(cosine),
+ gmx::c_rad2Deg * std::acos(cosine),
d0,
d1,
bllen[b]);
Rc,
*krf,
*crf,
- ONE_4PI_EPS0 / eps_r);
+ gmx::c_one4PiEps0 / eps_r);
// Make sure we don't lose resolution in pow() by casting real arg to double
real rmin = gmx::invcbrt(static_cast<double>(*krf * 2.0));
fprintf(fplog, "The electrostatics potential has its minimum at r = %g\n", rmin);
if (groupOptions.tau_t[g] >= 0)
{
nrdfCoupled += groupOptions.nrdf[g];
- kineticEnergy += groupOptions.nrdf[g] * 0.5 * groupOptions.ref_t[g] * BOLTZ;
+ kineticEnergy += groupOptions.nrdf[g] * 0.5 * groupOptions.ref_t[g] * gmx::c_boltz;
}
else
{
&& (reft > 0)) /* tau_t or ref_t = 0 means that no randomization is done */
{
randomize_group[gt] = true;
- boltzfac[gt] = BOLTZ * opts->ref_t[gt];
+ boltzfac[gt] = gmx::c_boltz * opts->ref_t[gt];
}
else
{
{
for (int gt = 0; gt < inputRecord.opts.ngtc; gt++)
{
- sd_.bd_rf[gt] = std::sqrt(2.0 * BOLTZ * inputRecord.opts.ref_t[gt]
+ sd_.bd_rf[gt] = std::sqrt(2.0 * gmx::c_boltz * inputRecord.opts.ref_t[gt]
/ (inputRecord.bd_fric * inputRecord.delta_t));
}
}
{
for (int gt = 0; gt < inputRecord.opts.ngtc; gt++)
{
- sd_.bd_rf[gt] = std::sqrt(2.0 * BOLTZ * inputRecord.opts.ref_t[gt]);
+ sd_.bd_rf[gt] = std::sqrt(2.0 * gmx::c_boltz * inputRecord.opts.ref_t[gt]);
}
}
}
{
for (int gt = 0; gt < inputRecord.opts.ngtc; gt++)
{
- real kT = BOLTZ * inputRecord.opts.ref_t[gt];
+ real kT = gmx::c_boltz * inputRecord.opts.ref_t[gt];
/* The mass is accounted for later, since this differs per atom */
sd_.sdsig[gt].V = std::sqrt(kT * (1 - sd_.sdc[gt].em * sd_.sdc[gt].em));
}
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
/* Set number of stddevs such that change of exceeding < 10^-9 */
constexpr real c_numSigma = 6.0;
/* Compute the maximally stretched angle */
- const real eqAngle = angleParams.harmonic.rA * DEG2RAD;
+ const real eqAngle = angleParams.harmonic.rA * gmx::c_deg2Rad;
const real fc = angleParams.harmonic.krA;
- const real maxAngle = eqAngle + c_numSigma * BOLTZ * temperature / ((numPartnerAtoms - 1) * fc);
+ const real maxAngle =
+ eqAngle + c_numSigma * gmx::c_boltz * temperature / ((numPartnerAtoms - 1) * fc);
if (maxAngle >= M_PI)
{
return -1;
if (re->bNPT)
{
/* revist the calculation for 5.0. Might be some improvements. */
- dpV = (beta[ap] * re->pres[ap] - beta[bp] * re->pres[bp]) * (Vol[b] - Vol[a]) / PRESFAC;
+ dpV = (beta[ap] * re->pres[ap] - beta[bp] * re->pres[bp]) * (Vol[b] - Vol[a]) / gmx::c_presfac;
if (bPrint)
{
fprintf(fplog, " dpV = %10.3e d = %10.3e\n", dpV, delta + dpV);
/* temperatures of different states*/
for (i = 0; i < re->nrepl; i++)
{
- re->beta[i] = 1.0 / (re->q[ereTEMP][i] * BOLTZ);
+ re->beta[i] = 1.0 / (re->q[ereTEMP][i] * gmx::c_boltz);
}
}
else
{
for (i = 0; i < re->nrepl; i++)
{
- re->beta[i] = 1.0 / (re->temp * BOLTZ); /* we have a single temperature */
+ re->beta[i] = 1.0 / (re->temp * gmx::c_boltz); /* we have a single temperature */
}
}
if (re->type == ereLAMBDA || re->type == ereTL)
#include "gromacs/gmxlib/network.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/math/vecdump.h"
#include "gromacs/mdlib/constr.h"
aS + 1,
mb + 1);
}
- shell[nsi].k += gmx::square(qS) * ONE_4PI_EPS0
+ shell[nsi].k += gmx::square(qS) * gmx::c_one4PiEps0
/ ffparams->iparams[type].polarize.alpha;
break;
case F_WATER_POL:
+ ffparams->iparams[type].wpol.al_y
+ ffparams->iparams[type].wpol.al_z)
/ 3.0;
- shell[nsi].k += gmx::square(qS) * ONE_4PI_EPS0 / alpha;
+ shell[nsi].k += gmx::square(qS) * gmx::c_one4PiEps0 / alpha;
break;
default: gmx_fatal(FARGS, "Death Horror: %s, %d", __FILE__, __LINE__);
}
}
fprintf(fplog, "\n The temperature for test particle insertion is %.3f K\n\n", temp);
}
- beta = 1.0 / (BOLTZ * temp);
+ beta = 1.0 / (gmx::c_boltz * temp);
/* Number of insertions per frame */
nsteps = inputrec->nsteps;
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
bonds.push_back(offset + at1 + 1);
bonds.push_back(offset + at2 + 1);
bondLengths.push_back(static_cast<double>(mtop->ffparams.iparams[contype].constr.dA)
- / BOHR2NM);
+ / gmx::c_bohr2Nm);
}
for (int ncon = 0; ncon < nsettle; ++ncon)
bonds.push_back(offset + h1 + 1);
bonds.push_back(offset + h2 + 1);
bondLengths.push_back(static_cast<double>(mtop->ffparams.iparams[contype].constr.dA)
- / BOHR2NM);
+ / gmx::c_bohr2Nm);
bondLengths.push_back(static_cast<double>(mtop->ffparams.iparams[contype].constr.dA)
- / BOHR2NM);
+ / gmx::c_bohr2Nm);
bondLengths.push_back(static_cast<double>(mtop->ffparams.iparams[contype].constr.dB)
- / BOHR2NM);
+ / gmx::c_bohr2Nm);
}
nAtomsMol.push_back(type->atoms.nr);
std::vector<double> convertedCoords;
for (auto& coord : coords)
{
- convertedCoords.push_back(static_cast<double>(coord[0]) / BOHR2NM);
- convertedCoords.push_back(static_cast<double>(coord[1]) / BOHR2NM);
- convertedCoords.push_back(static_cast<double>(coord[2]) / BOHR2NM);
+ convertedCoords.push_back(static_cast<double>(coord[0]) / gmx::c_bohr2Nm);
+ convertedCoords.push_back(static_cast<double>(coord[1]) / gmx::c_bohr2Nm);
+ convertedCoords.push_back(static_cast<double>(coord[2]) / gmx::c_bohr2Nm);
}
// sending array of coordinates
MCL_receive(&*coords.begin(), 3 * natoms, TYPE_DOUBLE, 0);
for (int j = 0; j < natoms; ++j)
{
- (*x)[j][0] = static_cast<real>(coords[j * 3] * BOHR2NM);
- (*x)[j][1] = static_cast<real>(coords[j * 3 + 1] * BOHR2NM);
- (*x)[j][2] = static_cast<real>(coords[j * 3 + 2] * BOHR2NM);
+ (*x)[j][0] = static_cast<real>(coords[j * 3] * gmx::c_bohr2Nm);
+ (*x)[j][1] = static_cast<real>(coords[j * 3 + 1] * gmx::c_bohr2Nm);
+ (*x)[j][2] = static_cast<real>(coords[j * 3 + 2] * gmx::c_bohr2Nm);
}
}
void gmx::MimicCommunicator::sendEnergies(real energy)
{
- double convertedEnergy = energy / (HARTREE2KJ * AVOGADRO);
+ double convertedEnergy = energy / (gmx::c_hartree2Kj * gmx::c_avogadro);
MCL_send(&convertedEnergy, 1, TYPE_DOUBLE, 0);
}
std::vector<double> convertedForce;
for (int j = 0; j < natoms; ++j)
{
- convertedForce.push_back(static_cast<real>(forces[j][0]) / HARTREE_BOHR2MD);
- convertedForce.push_back(static_cast<real>(forces[j][1]) / HARTREE_BOHR2MD);
- convertedForce.push_back(static_cast<real>(forces[j][2]) / HARTREE_BOHR2MD);
+ convertedForce.push_back(static_cast<real>(forces[j][0]) / gmx::c_hartreeBohr2Md);
+ convertedForce.push_back(static_cast<real>(forces[j][1]) / gmx::c_hartreeBohr2Md);
+ convertedForce.push_back(static_cast<real>(forces[j][2]) / gmx::c_hartreeBohr2Md);
}
MCL_send(&*convertedForce.begin(), convertedForce.size(), TYPE_DOUBLE, 0);
}
{
for (int j = 0; j <= i; j++)
{
- real invMass = PRESFAC * (4 * M_PI * M_PI * inputrec_->compress[i][j])
+ real invMass = c_presfac * (4 * M_PI * M_PI * inputrec_->compress[i][j])
/ (3 * inputrec_->tau_p * inputrec_->tau_p * maxBoxLength);
if (invMass > 0)
{
* track of unwrapped box diagonal elements. This case is
* excluded in integratorHasConservedEnergyQuantity().
*/
- energy += volume * trace(inputrec_->ref_p) / (DIM * PRESFAC);
+ energy += volume * trace(inputrec_->ref_p) / (DIM * c_presfac);
return energy;
}
}
const real referenceKineticEnergy =
- 0.5 * temperatureCouplingData.referenceTemperature[temperatureGroup] * BOLTZ
+ 0.5 * temperatureCouplingData.referenceTemperature[temperatureGroup] * gmx::c_boltz
* temperatureCouplingData.numDegreesOfFreedom[temperatureGroup];
const real newKineticEnergy =
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
void matrix_convert(matrix box, const rvec vec, const rvec angleInDegrees)
{
rvec angle;
- svmul(DEG2RAD, angleInDegrees, angle);
+ svmul(gmx::c_deg2Rad, angleInDegrees, angle);
box[XX][XX] = vec[XX];
box[YY][XX] = vec[YY] * cos(angle[ZZ]);
box[YY][YY] = vec[YY] * sin(angle[ZZ]);
#include "gromacs/gmxlib/network.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/mdlib/gmx_omp_nthreads.h"
{
if (pull_coordinate_is_angletype(pcrd))
{
- return DEG2RAD;
+ return gmx::c_deg2Rad;
}
else
{
{
if (pull_coordinate_is_angletype(pcrd))
{
- return RAD2DEG;
+ return gmx::c_rad2Deg;
}
else
{
*
* Copyright (c) 2009,2010,2011,2012,2013 by the GROMACS development team.
* Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#include <string>
#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/utility/cstringutil.h"
#include "gromacs/utility/exceptions.h"
*
* Copyright (c) 2009,2010,2011,2012,2013 by the GROMACS development team.
* Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
GMX_THROW(gmx::InvalidInputError("Angle cutoff should be > 0"));
}
- surf->angcut *= DEG2RAD;
+ surf->angcut *= gmx::c_deg2Rad;
surf->distccut = -std::cos(surf->angcut);
surf->targetbinsize = surf->angcut / 2;
surf->ntbins = static_cast<int>(M_PI / surf->targetbinsize);
- surf->tbinsize = (180.0 / surf->ntbins) * DEG2RAD;
+ surf->tbinsize = (180.0 / surf->ntbins) * gmx::c_deg2Rad;
snew(surf->tbin, static_cast<int>(M_PI / surf->tbinsize) + 1);
surf->maxbins = 0;
/* Convert the table from degrees to radians */
for (i = 0; i < td.nx; i++)
{
- td.x[i] *= DEG2RAD;
- td.f[i] *= RAD2DEG;
+ td.x[i] *= gmx::c_deg2Rad;
+ td.f[i] *= gmx::c_rad2Deg;
}
- td.tabscale *= RAD2DEG;
+ td.tabscale *= gmx::c_rad2Deg;
}
tab.n = td.nx;
tab.scale = td.tabscale;
ekin += 0.5 * atoms->atom[i].m * v[i][j] * v[i][j];
}
}
- temp1 = (2.0 * ekin) / (natom * DIM * BOLTZ);
- temp2 = (2.0 * ekin) / (natom * (DIM - 1) * BOLTZ);
+ temp1 = (2.0 * ekin) / (natom * DIM * gmx::c_boltz);
+ temp2 = (2.0 * ekin) / (natom * (DIM - 1) * gmx::c_boltz);
fprintf(stderr, "Kinetic energy: %g (kJ/mol)\n", ekin);
fprintf(stderr,
"Assuming the number of degrees of freedom to be "
*
* Copyright (c) 2010,2011,2012,2013,2014 by the GROMACS development team.
* Copyright (c) 2015,2016,2017,2018,2019 by the GROMACS development team.
- * Copyright (c) 2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
e_dir = 2.0 * info->q2all * gmx::invsqrt(info->q2allnr * r_coulomb * info->volume);
e_dir *= std::exp(-beta * beta * r_coulomb * r_coulomb);
- return ONE_4PI_EPS0 * e_dir;
+ return gmx::c_one4PiEps0 * e_dir;
}
#define SUMORDER 6
e_rec = std::sqrt(e_rec1 + e_rec2 + e_rec3);
- return ONE_4PI_EPS0 * e_rec;
+ return gmx::c_one4PiEps0 * e_rec;
}
break;
default: GMX_THROW(InternalError("invalid -g1 value"));
}
- dh.setPoint(n, angle * RAD2DEG, bPresent);
+ dh.setPoint(n, angle * gmx::c_rad2Deg, bPresent);
}
}
dh.finishFrame();
printf("Total volume %.2f +/- %.2f nm^3\n", Vaver, Verror);
printf("Number of molecules %d total mass %.2f Dalton\n", nmol_, mtot_);
- double RhoAver = mtot_ / (Vaver * 1e-24 * AVOGADRO);
+ double RhoAver = mtot_ / (Vaver * 1e-24 * gmx::c_avogadro);
double RhoError = gmx::square(RhoAver / Vaver) * Verror;
printf("Average molar mass: %.2f Dalton\n", mtot_ / nmol_);
{
totmass += surfaceSel.position(i).mass();
}
- const real density = totmass * AMU / (totvolume * NANO * NANO * NANO);
+ const real density = totmass * gmx::c_amu / (totvolume * gmx::c_nano * gmx::c_nano * gmx::c_nano);
vh.startFrame(frnr, fr.time);
vh.setPoint(0, totvolume);
vh.setPoint(1, density);
biod.pnpi.spb.ru / alexxy / gromacs.git / commitdiff