Skip to content

Commit

Permalink
Merge #3355
Browse files Browse the repository at this point in the history
3355: Remove Langevin noise correlation r=RudolfWeeber a=jngrad

Fixes #3343

Description of changes:
- decoupled translational and rotational Langevin RNG seeds
- add a test for Langevin noise correlation
- cleanup docs

Co-authored-by: Jean-Noël Grad <[email protected]>
Co-authored-by: Florian Weik <[email protected]>
  • Loading branch information
3 people authored Dec 7, 2019
2 parents 4d8b68b + f81e7b0 commit 846cc74
Show file tree
Hide file tree
Showing 10 changed files with 205 additions and 161 deletions.
20 changes: 8 additions & 12 deletions src/core/random.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -23,9 +23,8 @@
#define RANDOM_H

/** \file
A random generator
*/
* Random number generation using Philox.
*/

#include "errorhandling.hpp"

Expand All @@ -50,6 +49,7 @@ enum class RNGSalt : uint64_t {
FLUID = 0,
PARTICLES,
LANGEVIN,
LANGEVIN_ROT,
SALT_DPD,
THERMALIZED_BOND
};
Expand Down Expand Up @@ -96,8 +96,7 @@ inline void unseeded_error() {

/**
* @brief checks the seeded state and throws error if unseeded
*
**/
*/
inline void check_user_has_seeded() {
static bool unseeded_error_thrown = false;
if (!user_has_seeded && !unseeded_error_thrown) {
Expand All @@ -111,26 +110,24 @@ inline void check_user_has_seeded() {
*
* @param cnt Unused.
* @param seeds A vector of seeds, must be at least n_nodes long.
**/
*/
void mpi_random_seed(int cnt, std::vector<int> &seeds);

/**
* @brief Gets a string representation of the state of all
* the nodes.
* @brief Gets a string representation of the state of all the nodes.
*/
std::string mpi_random_get_stat();

/**
* @brief Set the seeds on all the node to the state represented
* by the string.
* The string representation must be one that was returned by
* mpi_random_get_stat.
* @ref mpi_random_get_stat.
*/
void mpi_random_set_stat(const std::vector<std::string> &stat);

/**
* @brief
* Get the state size of the random number generator
* @brief Get the state size of the random number generator
*/
int get_state_size_of_generator();

Expand All @@ -152,7 +149,6 @@ void init_random_seed(int seed);
* @brief Draws a random real number from the uniform distribution in the range
* [0,1)
*/

inline double d_random() {
using namespace Random;
check_user_has_seeded();
Expand Down
10 changes: 1 addition & 9 deletions src/core/thermostat.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -35,18 +35,14 @@
#include <iostream>
#include <unistd.h>

/** thermostat switch */
int thermo_switch = THERMO_OFF;
/** Temperature */
double temperature = 0.0;

/** True if the thermostat should act on virtual particles. */
bool thermo_virtual = true;

using Thermostat::GammaType;

namespace {
/** @name Langevin parameters sentinel
/** @name Langevin parameters sentinels.
* These functions return the sentinel value for the Langevin
* parameters, indicating that they have not been set yet.
*/
Expand All @@ -56,11 +52,7 @@ Utils::Vector3d sentinel(Utils::Vector3d) { return {-1.0, -1.0, -1.0}; }
/*@}*/
} // namespace

/* LANGEVIN THERMOSTAT */

/** Langevin friction coefficient gamma for translation. */
GammaType langevin_gamma = sentinel(GammaType{});
/** Friction coefficient gamma for rotation. */
GammaType langevin_gamma_rotation = sentinel(GammaType{});
GammaType langevin_pref1;
GammaType langevin_pref2;
Expand Down
79 changes: 43 additions & 36 deletions src/core/thermostat.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
#ifndef CORE_THERMOSTAT_HPP
#define CORE_THERMOSTAT_HPP
/** \file
* Implementation in \ref thermostat.cpp.
*/

#include "config.hpp"
Expand All @@ -38,7 +39,6 @@
#include <tuple>

/** \name Thermostat switches */
/************************************************************/
/*@{*/
#define THERMO_OFF 0
#define THERMO_LANGEVIN 1
Expand All @@ -62,55 +62,60 @@ using GammaType = double;
* exported variables
************************************************/

/** Switch determining which thermostat to use. This is a or'd value
of the different possible thermostats (defines: \ref THERMO_OFF,
\ref THERMO_LANGEVIN, \ref THERMO_DPD \ref THERMO_NPT_ISO). If it
is zero all thermostats are switched off and the temperature is
set to zero. */
/** Switch determining which thermostat(s) to use. This is a or'd value
* of the different possible thermostats (defines: \ref THERMO_OFF,
* \ref THERMO_LANGEVIN, \ref THERMO_DPD \ref THERMO_NPT_ISO). If it
* is zero all thermostats are switched off and the temperature is
* set to zero.
*/
extern int thermo_switch;

/** temperature. */
/** Temperature of the thermostat. */
extern double temperature;

/** True if the thermostat should act on virtual particles. */
extern bool thermo_virtual;

/** Langevin friction coefficient gamma. */
/** Langevin friction coefficient gamma for translation. */
extern Thermostat::GammaType langevin_gamma;
/** Langevin friction coefficient gamma. */
/** Langevin friction coefficient gamma for rotation. */
extern Thermostat::GammaType langevin_gamma_rotation;

/** Friction coefficient for nptiso-thermostat's inline-function
* friction_therm0_nptiso */
/** Friction coefficient for nptiso-thermostat's function
* @ref friction_therm0_nptiso
*/
extern double nptiso_gamma0;
/** Friction coefficient for nptiso-thermostat's inline-function
* friction_thermV_nptiso */
/** Friction coefficient for nptiso-thermostat's function
* @ref friction_thermV_nptiso
*/
extern double nptiso_gammav;

/** Langevin RNG counter, used for both translation and rotation. */
extern std::unique_ptr<Utils::Counter<uint64_t>> langevin_rng_counter;

/************************************************
* functions
************************************************/

/** only require seed if rng is not initialized */
/** Only require seed if rng is not initialized. */
bool langevin_is_seed_required();

/** philox functionality: increment, get/set */
/** @name philox functionality: increment, get/set */
/*@{*/
void langevin_rng_counter_increment();
void langevin_set_rng_state(uint64_t counter);
uint64_t langevin_get_rng_state();
/*@}*/

/** initialize constants of the thermostat on
start of integration */
/** Initialize constants of the thermostat at the start of integration */
void thermo_init();

#ifdef NPT
/** add velocity-dependent noise and friction for NpT-sims to the particle's
velocity
@param vj j-component of the velocity
@return j-component of the noise added to the velocity, also scaled by
dt (contained in prefactors)
/** Add velocity-dependent noise and friction for NpT-sims to the particle's
* velocity
* @param vj j-component of the velocity
* @return j-component of the noise added to the velocity, also scaled by
* dt (contained in prefactors)
*/
inline double friction_therm0_nptiso(double vj) {
extern double nptiso_pref1, nptiso_pref2;
Expand All @@ -123,8 +128,9 @@ inline double friction_therm0_nptiso(double vj) {
return 0.0;
}

/** add p_diff-dependent noise and friction for NpT-sims to \ref
* nptiso_struct::p_diff */
/** Add p_diff-dependent noise and friction for NpT-sims to \ref
* nptiso_struct::p_diff
*/
inline double friction_thermV_nptiso(double p_diff) {
extern double nptiso_pref3, nptiso_pref4;
if (thermo_switch & THERMO_NPT_ISO) {
Expand All @@ -137,11 +143,11 @@ inline double friction_thermV_nptiso(double p_diff) {
}
#endif

/** Langevin thermostat core function.
Collects the particle velocity (different for ENGINE, PARTICLE_ANISOTROPY).
Collects the langevin parameters kt, gamma (different for
LANGEVIN_PER_PARTICLE). Applies the noise and friction term.
*/
/** Langevin thermostat for particle translational velocities.
* Collects the particle velocity (different for ENGINE, PARTICLE_ANISOTROPY).
* Collects the langevin parameters kT, gamma (different for
* LANGEVIN_PER_PARTICLE). Applies the noise and friction term.
*/
inline Utils::Vector3d friction_thermo_langevin(Particle const &p) {
// Early exit for virtual particles without thermostat
if (p.p.is_virtual && !thermo_virtual) {
Expand All @@ -167,7 +173,7 @@ inline Utils::Vector3d friction_thermo_langevin(Particle const &p) {
}
#endif /* LANGEVIN_PER_PARTICLE */

// Get velocity effective in the thermostatting
// Get effective velocity in the thermostatting
#ifdef ENGINE
auto const &velocity = (p.p.swim.v_swim != 0)
? p.m.v - p.p.swim.v_swim * p.r.calc_director()
Expand Down Expand Up @@ -201,10 +207,11 @@ inline Utils::Vector3d friction_thermo_langevin(Particle const &p) {
}

#ifdef ROTATION
/** set the particle torques to the friction term, i.e. \f$\tau_i=-\gamma w_i +
\xi_i\f$.
The same friction coefficient \f$\gamma\f$ is used as that for translation.
*/
/** Langevin thermostat for particle angular velocities.
* Collects the particle velocity (different for PARTICLE_ANISOTROPY).
* Collects the langevin parameters kT, gamma_rot (different for
* LANGEVIN_PER_PARTICLE). Applies the noise and friction term.
*/
inline Utils::Vector3d friction_thermo_langevin_rotation(const Particle &p) {
extern Thermostat::GammaType langevin_pref2_rotation;

Expand All @@ -228,8 +235,8 @@ inline Utils::Vector3d friction_thermo_langevin_rotation(const Particle &p) {
using Random::v_noise;

// Here the thermostats happens
auto const noise =
v_noise<RNGSalt::LANGEVIN>(langevin_rng_counter->value(), p.p.identity);
auto const noise = v_noise<RNGSalt::LANGEVIN_ROT>(
langevin_rng_counter->value(), p.p.identity);
#ifdef PARTICLE_ANISOTROPY
return hadamard_product(langevin_pref_friction_buf, p.m.omega) +
hadamard_product(langevin_pref_noise_buf, noise);
Expand Down
Loading

0 comments on commit 846cc74

Please sign in to comment.