diff --git a/doc/sphinx/advanced_methods.rst b/doc/sphinx/advanced_methods.rst index a6f699890d9..95a7f5ec087 100644 --- a/doc/sphinx/advanced_methods.rst +++ b/doc/sphinx/advanced_methods.rst @@ -1486,7 +1486,7 @@ Particle polarizability with thermalized cold Drude oscillators .. note:: - Requires features ``THOLE``, ``P3M``, ``LANGEVIN_PER_PARTICLE``. + Requires features ``THOLE``, ``P3M``, ``THERMOSTAT_PER_PARTICLE``. .. note:: @@ -1498,7 +1498,7 @@ charge) to a particle (Drude core) that mimics an electron cloud which can be elongated to create a dynamically inducible dipole. The energetic minimum of the Drude charge can be obtained self-consistently, which requires several iterations of the system's electrostatics and is usually considered -computational expensive. However, with thermalized cold Drude oscillators, the +computationally expensive. However, with thermalized cold Drude oscillators, the distance between Drude charge and core is coupled to a thermostat so that it fluctuates around the SCF solution. This thermostat is kept at a low temperature compared to the global temperature to minimize the heat flow into @@ -1515,8 +1515,8 @@ In |es|, the basic ingredients to simulate such a system are split into three bo The system-wide thermostat has to be applied to the centre of mass and not to the core particle directly. Therefore, the particles have to be excluded from -global thermostatting. With ``LANGEVIN_PER_PARTICLE`` enabled, we set the -temperature and friction coefficient of the Drude complex to zero, which allows +global thermostatting. With ``THERMOSTAT_PER_PARTICLE`` enabled, we set the +friction coefficient of the Drude complex to zero, which allows to still use a global Langevin thermostat for non-polarizable particles. As the Drude charge should not alter the *charge* or *mass* of the Drude diff --git a/doc/sphinx/installation.rst b/doc/sphinx/installation.rst index 8ca57ab431b..3d4949fe3f1 100644 --- a/doc/sphinx/installation.rst +++ b/doc/sphinx/installation.rst @@ -317,11 +317,8 @@ General features additional degrees of freedom, which for example means that the kinetic energy changes at constant temperature is twice as large. -- ``LANGEVIN_PER_PARTICLE`` Allows to choose the Langevin temperature and friction coefficient - per particle. - -- ``BROWNIAN_PER_PARTICLE`` Allows to choose the Brownian temperature and friction coefficient - per particle. +- ``THERMOSTAT_PER_PARTICLE`` Allows setting a per-particle friction + coefficient for the Langevin and Brownian thermostats. - ``ROTATIONAL_INERTIA`` diff --git a/doc/sphinx/inter_non-bonded.rst b/doc/sphinx/inter_non-bonded.rst index 51061e3bb78..441bea995aa 100644 --- a/doc/sphinx/inter_non-bonded.rst +++ b/doc/sphinx/inter_non-bonded.rst @@ -659,7 +659,7 @@ To use the script, compile espresso with the following features: #define EXTERNAL_FORCES #define MASS - #define LANGEVIN_PER_PARTICLE + #define THERMOSTAT_PER_PARTICLE #define ROTATION #define ROTATIONAL_INERTIA #define ELECTROSTATICS diff --git a/doc/sphinx/running.rst b/doc/sphinx/running.rst index 773325d982b..5d94a2d98af 100644 --- a/doc/sphinx/running.rst +++ b/doc/sphinx/running.rst @@ -499,7 +499,7 @@ the anisotropic diffusion of anisotropic colloids (rods, etc.). Using the Langevin thermostat, it is possible to set a temperature and a friction coefficient for every particle individually via the feature -``LANGEVIN_PER_PARTICLE``. Consult the reference of the ``part`` command +``THERMOSTAT_PER_PARTICLE``. Consult the reference of the ``part`` command (chapter :ref:`Setting up particles`) for information on how to achieve this. .. _Brownian thermostat: @@ -532,7 +532,7 @@ Best explained in an example:: where ``gamma`` (hereinafter :math:`\gamma`) is a viscous friction coefficient. In terms of the Python interface and setup, the Brownian thermostat is very similar to the :ref:`Langevin thermostat`. The feature -``BROWNIAN_PER_PARTICLE`` is used to control the per-particle +``THERMOSTAT_PER_PARTICLE`` is used to control the per-particle temperature and the friction coefficient setup. The major differences are its internal integrator implementation and other temporal constraints. The integrator is still a symplectic velocity Verlet-like one. diff --git a/maintainer/configs/maxset.hpp b/maintainer/configs/maxset.hpp index cc8c78ffc92..74368803f10 100644 --- a/maintainer/configs/maxset.hpp +++ b/maintainer/configs/maxset.hpp @@ -32,8 +32,7 @@ #define BOND_CONSTRAINT #define COLLISION_DETECTION -#define LANGEVIN_PER_PARTICLE -#define BROWNIAN_PER_PARTICLE +#define THERMOSTAT_PER_PARTICLE #define NPT diff --git a/maintainer/configs/no_rotation.hpp b/maintainer/configs/no_rotation.hpp index 8e2fd495f81..5f88f292a1e 100644 --- a/maintainer/configs/no_rotation.hpp +++ b/maintainer/configs/no_rotation.hpp @@ -29,8 +29,7 @@ // Geometry, equation of motion, thermostat/barostat #define MASS #define EXTERNAL_FORCES -#define LANGEVIN_PER_PARTICLE -#define BROWNIAN_PER_PARTICLE +#define THERMOSTAT_PER_PARTICLE #define BOND_CONSTRAINT #define NPT #define DPD diff --git a/samples/chamber_game.py b/samples/chamber_game.py index 55ae8f40891..544124508bc 100644 --- a/samples/chamber_game.py +++ b/samples/chamber_game.py @@ -32,7 +32,7 @@ from espressomd.visualization_opengl import openGLLive, KeyboardButtonEvent, KeyboardFireEvent required_features = ["LENNARD_JONES", "WCA", "MASS", - "EXTERNAL_FORCES", "LANGEVIN_PER_PARTICLE"] + "EXTERNAL_FORCES", "THERMOSTAT_PER_PARTICLE"] espressomd.assert_features(required_features) print("""THE CHAMBER GAME @@ -61,13 +61,11 @@ # PARAMETERS # PHYSICS -temperature_snake = 0.0 gamma_snake_head = 1.0 gamma_snake_bead = 15.0 -temperature_bubbles = 10000.0 -temp_l = temperature_bubbles -temp_r = temperature_bubbles +temp_l = 10000.0 +temp_r = temp_l temp_max = 1e5 gamma_bubbles = 0.5 @@ -180,7 +178,6 @@ type=snake_head_type, fix=[False, False, True], mass=snake_head_mass, - temp=temperature_snake, gamma=gamma_snake_head) else: system.part.add( @@ -192,7 +189,6 @@ type=snake_bead_type, fix=[False, False, True], mass=snake_bead_mass, - temp=temperature_snake, gamma=gamma_snake_bead) # NB INTER @@ -262,7 +258,6 @@ type=bubble_type, fix=[False, False, True], mass=bubble_mass, - temp=temperature_bubbles, gamma=gamma_bubbles) testid = len(system.part) - 1 n += 1 @@ -453,10 +448,8 @@ def T_to_g(temp): Nl = len(pl) Nr = len(pr) for p in pl: - p.temp = temp_l p.gamma = T_to_g(temp_l) for p in pr: - p.temp = temp_r p.gamma = T_to_g(temp_r) w = visualizer.specs['window_size'] diff --git a/samples/drude_bmimpf6.py b/samples/drude_bmimpf6.py index 2b2767e07a0..6f8fbee21d0 100644 --- a/samples/drude_bmimpf6.py +++ b/samples/drude_bmimpf6.py @@ -27,7 +27,7 @@ import espressomd required_features = ["LENNARD_JONES", "P3M", "MASS", "ROTATION", "ROTATIONAL_INERTIA", "VIRTUAL_SITES_RELATIVE", - "THOLE", "LANGEVIN_PER_PARTICLE"] + "THOLE", "THERMOSTAT_PER_PARTICLE"] espressomd.assert_features(required_features) import espressomd.observables @@ -215,7 +215,7 @@ def combination_rule_sigma(rule, sig1, sig2): system.part.add( id=rid, type=types["BMIM_COM"], pos=pos_com, mass=masses["BMIM_COM"], rinertia=[646.284, 585.158, 61.126], - temp=0, gamma=0, rotation=[1, 1, 1]) + gamma=0, rotation=[1, 1, 1]) cation_com_ids.append(rid) com_id = rid rid += 1 diff --git a/src/config/features.def b/src/config/features.def index 13e57d755f0..6a203625f25 100644 --- a/src/config/features.def +++ b/src/config/features.def @@ -18,8 +18,7 @@ EXTERNAL_FORCES MASS EXCLUSIONS BOND_CONSTRAINT -LANGEVIN_PER_PARTICLE -BROWNIAN_PER_PARTICLE +THERMOSTAT_PER_PARTICLE COLLISION_DETECTION NPT ENGINE implies ROTATION, EXTERNAL_FORCES diff --git a/src/config/myconfig-default.hpp b/src/config/myconfig-default.hpp index 8eb7d1285fb..72a77a8ad11 100644 --- a/src/config/myconfig-default.hpp +++ b/src/config/myconfig-default.hpp @@ -32,8 +32,7 @@ #define MASS #define PARTICLE_ANISOTROPY #define EXTERNAL_FORCES -#define LANGEVIN_PER_PARTICLE -#define BROWNIAN_PER_PARTICLE +#define THERMOSTAT_PER_PARTICLE #define BOND_CONSTRAINT #define NPT #define DPD diff --git a/src/core/Particle.hpp b/src/core/Particle.hpp index 9bffa329a70..82bfe54cdc3 100644 --- a/src/core/Particle.hpp +++ b/src/core/Particle.hpp @@ -146,8 +146,7 @@ struct ParticleProperties { static constexpr bool is_virtual = false; #endif /* VIRTUAL_SITES */ -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) - double T = -1.; +#ifdef THERMOSTAT_PER_PARTICLE #ifndef PARTICLE_ANISOTROPY double gamma = -1.; #else @@ -161,7 +160,7 @@ struct ParticleProperties { Utils::Vector3d gamma_rot = {-1., -1., -1.}; #endif // ROTATIONAL_INERTIA #endif // ROTATION -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef EXTERNAL_FORCES /** flag whether to fix a particle in space. @@ -220,13 +219,12 @@ struct ParticleProperties { #endif #endif /* VIRTUAL_SITES */ -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) - ar &T; +#ifdef THERMOSTAT_PER_PARTICLE ar γ #ifdef ROTATION ar &gamma_rot; #endif -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef EXTERNAL_FORCES ar &ext_flag; ar &ext_force; diff --git a/src/core/particle_data.cpp b/src/core/particle_data.cpp index 43b5910d6a3..73db3218e78 100644 --- a/src/core/particle_data.cpp +++ b/src/core/particle_data.cpp @@ -121,8 +121,7 @@ using UpdatePropertyMessage = boost::variant &Prop::vs_relative> #endif #endif -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) - , UpdateProperty +#ifdef THERMOSTAT_PER_PARTICLE #ifndef PARTICLE_ANISOTROPY , UpdateProperty #else @@ -135,7 +134,7 @@ using UpdatePropertyMessage = boost::variant , UpdateProperty #endif // PARTICLE_ANISOTROPY #endif // ROTATION -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef EXTERNAL_FORCES , UpdateProperty , UpdateProperty @@ -871,11 +870,7 @@ void set_particle_torque_lab(int part, const Utils::Vector3d &torque_lab) { } #endif -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) -void set_particle_temperature(int part, double T) { - mpi_update_particle_property(part, T); -} - +#ifdef THERMOSTAT_PER_PARTICLE #ifndef PARTICLE_ANISOTROPY void set_particle_gamma(int part, double gamma) { mpi_update_particle_property(part, gamma); @@ -900,7 +895,7 @@ void set_particle_gamma_rot(int part, Utils::Vector3d gamma_rot) { } #endif // PARTICLE_ANISOTROPY #endif // ROTATION -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef EXTERNAL_FORCES #ifdef ROTATION @@ -1249,7 +1244,7 @@ void pointer_to_fix(Particle const *p, const uint8_t *&res) { } #endif -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) +#ifdef THERMOSTAT_PER_PARTICLE void pointer_to_gamma(Particle const *p, double const *&res) { #ifndef PARTICLE_ANISOTROPY res = &(p->p.gamma); @@ -1268,10 +1263,7 @@ void pointer_to_gamma_rot(Particle const *p, double const *&res) { } #endif // ROTATION -void pointer_to_temperature(Particle const *p, double const *&res) { - res = &(p->p.T); -} -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef ENGINE void pointer_to_swimming(Particle const *p, diff --git a/src/core/particle_data.hpp b/src/core/particle_data.hpp index b2410d6c37e..8fab4251da9 100644 --- a/src/core/particle_data.hpp +++ b/src/core/particle_data.hpp @@ -248,13 +248,7 @@ void set_particle_vs_relative(int part, int vs_relative_to, double vs_distance, Utils::Quaternion const &rel_ori); #endif -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) -/** Call only on the master node: set particle temperature. - * @param part the particle. - * @param T its new temperature. - */ -void set_particle_temperature(int part, double T); - +#ifdef THERMOSTAT_PER_PARTICLE /** Call only on the master node: set particle frictional coefficient. * @param part the particle. * @param gamma its new frictional coefficient. @@ -271,7 +265,7 @@ void set_particle_gamma_rot(int part, double gamma); void set_particle_gamma_rot(int part, Utils::Vector3d gamma_rot); #endif #endif -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef EXTERNAL_FORCES #ifdef ROTATION @@ -396,13 +390,12 @@ void pointer_to_ext_torque(Particle const *p, double const *&res2); void pointer_to_fix(Particle const *p, const uint8_t *&res); #endif -#if defined(LANGEVIN_PER_PARTICLE) || defined(BROWNIAN_PER_PARTICLE) +#ifdef THERMOSTAT_PER_PARTICLE void pointer_to_gamma(Particle const *p, double const *&res); -void pointer_to_temperature(Particle const *p, double const *&res); #ifdef ROTATION void pointer_to_gamma_rot(Particle const *p, double const *&res); #endif -#endif // LANGEVIN_PER_PARTICLE || BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE #ifdef ENGINE void pointer_to_swimming(Particle const *p, diff --git a/src/core/thermostats/brownian_inline.hpp b/src/core/thermostats/brownian_inline.hpp index f3ed6af7f1d..733d5c4510c 100644 --- a/src/core/thermostats/brownian_inline.hpp +++ b/src/core/thermostats/brownian_inline.hpp @@ -44,7 +44,7 @@ inline Utils::Vector3d bd_drag(Thermostat::GammaType const &brownian_gamma, // The friction tensor Z from the Eq. (14.31) of Schlick2010: Thermostat::GammaType gamma; -#ifdef BROWNIAN_PER_PARTICLE +#ifdef THERMOSTAT_PER_PARTICLE if (p.p.gamma >= Thermostat::GammaType{}) { gamma = p.p.gamma; } else @@ -106,7 +106,7 @@ inline Utils::Vector3d bd_drag_vel(Thermostat::GammaType const &brownian_gamma, // The friction tensor Z from the eq. (14.31) of Schlick2010: Thermostat::GammaType gamma; -#ifdef BROWNIAN_PER_PARTICLE +#ifdef THERMOSTAT_PER_PARTICLE if (p.p.gamma >= Thermostat::GammaType{}) { gamma = p.p.gamma; } else @@ -170,43 +170,18 @@ inline Utils::Vector3d bd_random_walk(BrownianThermostat const &brownian, // skip the translation thermalizing for virtual sites unless enabled if (p.p.is_virtual && !thermo_virtual) return {}; - // first, set defaults - Thermostat::GammaType sigma_pos = brownian.sigma_pos; - // Override defaults if per-particle values for T and gamma are given -#ifdef BROWNIAN_PER_PARTICLE + Thermostat::GammaType sigma_pos = brownian.sigma_pos; +#ifdef THERMOSTAT_PER_PARTICLE + // override default if particle-specific gamma if (p.p.gamma >= Thermostat::GammaType{}) { - // Is a particle-specific temperature also specified? - if (p.p.T >= 0.) { - if (p.p.T > 0.0) { - sigma_pos = BrownianThermostat::sigma(p.p.T, p.p.gamma); - } else { - // just an indication of the infinity - sigma_pos = Thermostat::GammaType{}; - } - } else - // default temperature but particle-specific gamma - if (temperature > 0.0) { + if (temperature > 0.0) { sigma_pos = BrownianThermostat::sigma(temperature, p.p.gamma); } else { sigma_pos = Thermostat::GammaType{}; } - } // particle-specific gamma - else { - // No particle-specific gamma, but is there particle-specific temperature - if (p.p.T >= 0.) { - if (p.p.T > 0.0) { - sigma_pos = BrownianThermostat::sigma(p.p.T, brownian.gamma); - } else { - // just an indication of the infinity - sigma_pos = Thermostat::GammaType{}; - } - } else { - // default values for both - sigma_pos = brownian.sigma_pos; - } } -#endif // BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE // Eq. (14.37) is factored by the Gaussian noise (12.22) with its squared // magnitude defined in the second eq. (14.38), Schlick2010. @@ -270,22 +245,6 @@ inline Utils::Vector3d bd_random_walk_vel(BrownianThermostat const &brownian, // skip the translation thermalizing for virtual sites unless enabled if (p.p.is_virtual && !thermo_virtual) return {}; - // Just a square root of kT, see eq. (10.2.17) and comments in 2 paragraphs - // afterwards, Pottier2010 - double sigma_vel; - - // Override defaults if per-particle values for T and gamma are given -#ifdef BROWNIAN_PER_PARTICLE - // Is a particle-specific temperature specified? - if (p.p.T >= 0.) { - sigma_vel = BrownianThermostat::sigma(p.p.T); - } else { - sigma_vel = brownian.sigma_vel; - } -#else - // defaults - sigma_vel = brownian.sigma_vel; -#endif // BROWNIAN_PER_PARTICLE auto const noise = Random::noise_gaussian( brownian.rng_counter(), brownian.rng_seed(), p.identity()); @@ -302,7 +261,7 @@ inline Utils::Vector3d bd_random_walk_vel(BrownianThermostat const &brownian, // (14.31) of Schlick2010. A difference is the mass factor to the friction // tensor. The noise is Gaussian according to the convention at p. 237 // (last paragraph), Pottier2010. - velocity[j] += sigma_vel * noise[j] / sqrt(p.p.mass); + velocity[j] += brownian.sigma_vel * noise[j] / sqrt(p.p.mass); } } return velocity; @@ -321,7 +280,7 @@ bd_drag_rot(Thermostat::GammaType const &brownian_gamma_rotation, Particle &p, double dt) { Thermostat::GammaType gamma; -#ifdef BROWNIAN_PER_PARTICLE +#ifdef THERMOSTAT_PER_PARTICLE if (p.p.gamma_rot >= Thermostat::GammaType{}) { gamma = p.p.gamma_rot; } else @@ -363,7 +322,7 @@ bd_drag_vel_rot(Thermostat::GammaType const &brownian_gamma_rotation, Particle const &p) { Thermostat::GammaType gamma; -#ifdef BROWNIAN_PER_PARTICLE +#ifdef THERMOSTAT_PER_PARTICLE if (p.p.gamma_rot >= Thermostat::GammaType{}) { gamma = p.p.gamma_rot; } else @@ -398,43 +357,18 @@ bd_drag_vel_rot(Thermostat::GammaType const &brownian_gamma_rotation, inline Utils::Quaternion bd_random_walk_rot(BrownianThermostat const &brownian, Particle const &p, double dt) { - // first, set defaults - Thermostat::GammaType sigma_pos = brownian.sigma_pos_rotation; - // Override defaults if per-particle values for T and gamma are given -#ifdef BROWNIAN_PER_PARTICLE + Thermostat::GammaType sigma_pos = brownian.sigma_pos_rotation; +#ifdef THERMOSTAT_PER_PARTICLE + // override default if particle-specific gamma if (p.p.gamma_rot >= Thermostat::GammaType{}) { - // Is a particle-specific temperature also specified? - if (p.p.T >= 0.) { - if (p.p.T > 0.0) { - sigma_pos = BrownianThermostat::sigma(p.p.T, p.p.gamma_rot); - } else { - // just an indication of the infinity - sigma_pos = {}; - } - } else if (temperature > 0.) { - // Default temperature but particle-specific gamma + if (temperature > 0.) { sigma_pos = BrownianThermostat::sigma(temperature, p.p.gamma_rot); } else { - // just an indication of the infinity - sigma_pos = {}; - } - } // particle-specific gamma - else { - // No particle-specific gamma, but is there particle-specific temperature - if (p.p.T >= 0.) { - if (p.p.T > 0.0) { - sigma_pos = BrownianThermostat::sigma(p.p.T, brownian.gamma_rotation); - } else { - // just an indication of the infinity - sigma_pos = {}; - } - } else { - // Default values for both - sigma_pos = brownian.sigma_pos_rotation; + sigma_pos = {}; // just an indication of the infinity } } -#endif // BROWNIAN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE Utils::Vector3d dphi = {}; auto const noise = Random::noise_gaussian( @@ -472,20 +406,7 @@ bd_random_walk_rot(BrownianThermostat const &brownian, Particle const &p, */ inline Utils::Vector3d bd_random_walk_vel_rot(BrownianThermostat const &brownian, Particle const &p) { - double sigma_vel; - - // Override defaults if per-particle values for T and gamma are given -#ifdef BROWNIAN_PER_PARTICLE - // Is a particle-specific temperature specified? - if (p.p.T >= 0.) { - sigma_vel = BrownianThermostat::sigma(p.p.T); - } else { - sigma_vel = brownian.sigma_vel_rotation; - } -#else - // set defaults - sigma_vel = brownian.sigma_vel_rotation; -#endif // BROWNIAN_PER_PARTICLE + auto const sigma_vel = brownian.sigma_vel_rotation; Utils::Vector3d domega{}; auto const noise = Random::noise_gaussian( diff --git a/src/core/thermostats/langevin_inline.hpp b/src/core/thermostats/langevin_inline.hpp index 3dd7eb467af..008d83b8a0f 100644 --- a/src/core/thermostats/langevin_inline.hpp +++ b/src/core/thermostats/langevin_inline.hpp @@ -35,7 +35,7 @@ /** 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. + * THERMOSTAT_PER_PARTICLE). Applies the noise and friction term. * @param[in] langevin Parameters * @param[in] p %Particle * @param[in] time_step Time step @@ -49,19 +49,17 @@ friction_thermo_langevin(LangevinThermostat const &langevin, Particle const &p, } // Determine prefactors for the friction and the noise term - // first, set defaults Thermostat::GammaType pref_friction = langevin.pref_friction; Thermostat::GammaType pref_noise = langevin.pref_noise; - // Override defaults if per-particle values for T and gamma are given -#ifdef LANGEVIN_PER_PARTICLE - if (p.p.gamma >= Thermostat::GammaType{} or p.p.T >= 0.) { - auto const kT = p.p.T >= 0. ? p.p.T : temperature; +#ifdef THERMOSTAT_PER_PARTICLE + // override default if particle-specific gamma + if (p.p.gamma >= Thermostat::GammaType{}) { auto const gamma = p.p.gamma >= Thermostat::GammaType{} ? p.p.gamma : langevin.gamma; pref_friction = -gamma; - pref_noise = LangevinThermostat::sigma(kT, time_step, gamma); + pref_noise = LangevinThermostat::sigma(temperature, time_step, gamma); } -#endif // LANGEVIN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE // Get effective velocity in the thermostatting #ifdef ENGINE @@ -98,7 +96,7 @@ friction_thermo_langevin(LangevinThermostat const &langevin, Particle const &p, /** 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. + * THERMOSTAT_PER_PARTICLE). Applies the noise and friction term. * @param[in] langevin Parameters * @param[in] p %Particle * @param[in] time_step Time step @@ -110,17 +108,16 @@ friction_thermo_langevin_rotation(LangevinThermostat const &langevin, auto pref_friction = -langevin.gamma_rotation; auto pref_noise = langevin.pref_noise_rotation; - // Override defaults if per-particle values for T and gamma are given -#ifdef LANGEVIN_PER_PARTICLE - if (p.p.gamma_rot >= Thermostat::GammaType{} or p.p.T >= 0.) { - auto const kT = p.p.T >= 0. ? p.p.T : temperature; +#ifdef THERMOSTAT_PER_PARTICLE + // override default if particle-specific gamma + if (p.p.gamma_rot >= Thermostat::GammaType{}) { auto const gamma = p.p.gamma_rot >= Thermostat::GammaType{} ? p.p.gamma_rot : langevin.gamma_rotation; pref_friction = -gamma; - pref_noise = LangevinThermostat::sigma(kT, time_step, gamma); + pref_noise = LangevinThermostat::sigma(temperature, time_step, gamma); } -#endif // LANGEVIN_PER_PARTICLE +#endif // THERMOSTAT_PER_PARTICLE auto const noise = Random::noise_uniform( langevin.rng_counter(), langevin.rng_seed(), p.p.identity); diff --git a/src/python/espressomd/drude_helpers.py b/src/python/espressomd/drude_helpers.py index 10273b9a2bf..21ed4b08e91 100644 --- a/src/python/espressomd/drude_helpers.py +++ b/src/python/espressomd/drude_helpers.py @@ -76,7 +76,7 @@ def add_drude_particle_to_core(system, harmonic_bond, thermalized_bond, gamma_off = 0.0 system.part.add(id=id_drude, pos=p_core.pos, type=type_drude, - q=q_drude, mass=mass_drude, temp=0, gamma=gamma_off) + q=q_drude, mass=mass_drude, gamma=gamma_off) if verbose: print( @@ -88,7 +88,6 @@ def add_drude_particle_to_core(system, harmonic_bond, thermalized_bond, p_core.add_bond((harmonic_bond, id_drude)) p_core.add_bond((thermalized_bond, id_drude)) - p_core.temp = 0.0 p_core.gamma = gamma_off if type_drude in drude_dict and not ( diff --git a/src/python/espressomd/particle_data.pxd b/src/python/espressomd/particle_data.pxd index be7e83ac72c..dc50d78194b 100644 --- a/src/python/espressomd/particle_data.pxd +++ b/src/python/espressomd/particle_data.pxd @@ -130,10 +130,7 @@ cdef extern from "particle_data.hpp": void set_particle_virtual(int part, int isVirtual) void pointer_to_virtual(const particle * P, const bint * & res) - IF LANGEVIN_PER_PARTICLE or BROWNIAN_PER_PARTICLE: - void set_particle_temperature(int part, double T) - void pointer_to_temperature(const particle * p, const double * & res) - + IF THERMOSTAT_PER_PARTICLE: IF PARTICLE_ANISOTROPY: void set_particle_gamma(int part, Vector3d gamma) ELSE: diff --git a/src/python/espressomd/particle_data.pyx b/src/python/espressomd/particle_data.pyx index 08cbe8c207d..3fa32c849a5 100644 --- a/src/python/espressomd/particle_data.pyx +++ b/src/python/espressomd/particle_data.pyx @@ -892,7 +892,7 @@ cdef class ParticleHandle: self.particle_data, ext_t) return array_locked([ext_t[0], ext_t[1], ext_t[2]]) - IF LANGEVIN_PER_PARTICLE or BROWNIAN_PER_PARTICLE: + IF THERMOSTAT_PER_PARTICLE: IF PARTICLE_ANISOTROPY: property gamma: """ @@ -902,8 +902,8 @@ cdef class ParticleHandle: gamma : :obj:`float` or (3,) array_like of :obj:`float` .. note:: - This needs features ``PARTICLE_ANISOTROPY`` and either - ``LANGEVIN_PER_PARTICLE`` or ``BROWNIAN_PER_PARTICLE``. + This needs features ``PARTICLE_ANISOTROPY`` and + ``THERMOSTAT_PER_PARTICLE``. See Also ---------- @@ -940,7 +940,7 @@ cdef class ParticleHandle: gamma : :obj:`float` .. note:: - This needs the feature ``LANGEVIN_PER_PARTICLE`` or ``BROWNIAN_PER_PARTICLE``. + This needs the feature ``THERMOSTAT_PER_PARTICLE``. See Also ---------- @@ -969,7 +969,7 @@ cdef class ParticleHandle: .. note:: This needs features ``ROTATION``, ``PARTICLE_ANISOTROPY`` - and either ``LANGEVIN_PER_PARTICLE`` or ``BROWNIAN_PER_PARTICLE``. + and ``THERMOSTAT_PER_PARTICLE``. """ @@ -1002,8 +1002,8 @@ cdef class ParticleHandle: gamma_rot : :obj:`float` .. note:: - This needs features ``ROTATION`` and either - ``LANGEVIN_PER_PARTICLE`` or ``BROWNIAN_PER_PARTICLE``. + This needs features ``ROTATION`` and + ``THERMOSTAT_PER_PARTICLE``. """ @@ -1019,29 +1019,6 @@ cdef class ParticleHandle: self.particle_data, gamma_rot) return gamma_rot[0] - property temp: - """ - Particle's temperature in the Langevin and Brownian thermostats. - - temp: :obj:`float` - - .. note:: - This needs the feature ``LANGEVIN_PER_PARTICLE`` or - ``BROWNIAN_PER_PARTICLE``. - - """ - - def __set__(self, _temp): - check_type_or_throw_except( - _temp, 1, float, "temp has to be a float.") - set_particle_temperature(self._id, _temp) - - def __get__(self): - self.update_particle_data() - cdef const double * temp = NULL - pointer_to_temperature(self.particle_data, temp) - return temp[0] - IF ROTATION: property rotation: """ diff --git a/testsuite/python/array_properties.py b/testsuite/python/array_properties.py index 2f28546f828..6e450f1cffa 100644 --- a/testsuite/python/array_properties.py +++ b/testsuite/python/array_properties.py @@ -175,7 +175,7 @@ def test_external_forces(self): self.assert_copy_is_writable(self.system.part[0].ext_force) self.assert_copy_is_writable(self.system.part[0].fix) - @utx.skipIfMissingFeatures(["ROTATION", "LANGEVIN_PER_PARTICLE", + @utx.skipIfMissingFeatures(["ROTATION", "THERMOSTAT_PER_PARTICLE", "PARTICLE_ANISOTROPY"]) def test_rot_aniso(self): self.assert_operator_usage_raises(self.system.part[0].gamma_rot) @@ -193,9 +193,9 @@ def test_lb(self): self.assert_operator_usage_raises(lbf[0, 0, 0].pressure_tensor_neq) self.assert_operator_usage_raises(lbf[0, 0, 0].population) - @utx.skipIfMissingFeatures(["LANGEVIN_PER_PARTICLE", + @utx.skipIfMissingFeatures(["THERMOSTAT_PER_PARTICLE", "PARTICLE_ANISOTROPY"]) - def test_langevinpp_aniso(self): + def test_thermostat_per_particle_aniso(self): self.assert_operator_usage_raises(self.system.part[0].gamma) check_array_writable(self.system.part[0].gamma) diff --git a/testsuite/python/brownian_dynamics.py b/testsuite/python/brownian_dynamics.py index c14dd5fb9ab..b8e25c149eb 100644 --- a/testsuite/python/brownian_dynamics.py +++ b/testsuite/python/brownian_dynamics.py @@ -40,7 +40,7 @@ def tearDown(self): self.system.thermostat.turn_off() self.system.integrator.set_vv() - def check_rng(self, per_particle_gamma=False, per_particle_temp=False): + def check_rng(self, per_particle_gamma=False): """Test for RNG consistency.""" kT = 1.1 @@ -51,17 +51,14 @@ def reset_particle(): p = system.part.add(pos=[0, 0, 0]) if espressomd.has_features("ROTATION"): p.rotation = [1, 1, 1] - if per_particle_gamma or per_particle_temp: - assert espressomd.has_features("BROWNIAN_PER_PARTICLE") if per_particle_gamma: + assert espressomd.has_features("THERMOSTAT_PER_PARTICLE") if espressomd.has_features("PARTICLE_ANISOTROPY"): p.gamma = 3 * [gamma / 2] else: p.gamma = gamma / 2 if espressomd.has_features("ROTATION"): p.gamma_rot = p.gamma * 1.5 - if per_particle_temp: - p.temp = 2 * kT return p system = self.system @@ -125,12 +122,10 @@ def test_01__rng(self): self.system.thermostat.set_brownian(kT=1, gamma=2) self.check_rng() - @utx.skipIfMissingFeatures("BROWNIAN_PER_PARTICLE") + @utx.skipIfMissingFeatures("THERMOSTAT_PER_PARTICLE") def test_01__rng_per_particle(self): """Test for RNG consistency.""" - self.check_rng(False, True) - self.check_rng(True, False) - self.check_rng(True, True) + self.check_rng(True) @utx.skipIfMissingFeatures("VIRTUAL_SITES") def test_07__virtual(self): diff --git a/testsuite/python/brownian_dynamics_stats.py b/testsuite/python/brownian_dynamics_stats.py index d5433a24d50..b73073a522f 100644 --- a/testsuite/python/brownian_dynamics_stats.py +++ b/testsuite/python/brownian_dynamics_stats.py @@ -77,7 +77,7 @@ def test_vel_dist_global_temp_initial_forces(self): """ self.check_vel_dist_global_temp(True, loops=170) - @utx.skipIfMissingFeatures("BROWNIAN_PER_PARTICLE") + @utx.skipIfMissingFeatures("THERMOSTAT_PER_PARTICLE") def test_vel_dist_per_particle(self): """Test Brownian dynamics with particle-specific kT and gamma. Covers all combinations of particle-specific gamma and temp set or not set. @@ -88,23 +88,13 @@ def test_vel_dist_per_particle(self): kT = 0.9 gamma = 3.2 gamma2 = 4.3 - kT2 = 1.5 system.thermostat.set_brownian(kT=kT, gamma=gamma, seed=41) loops = 200 v_minmax = 5 bins = 4 error_tol = 0.012 self.check_per_particle( - N, kT, kT2, gamma2, loops, v_minmax, bins, error_tol) - - def setup_diff_mass_rinertia(self, p): - if espressomd.has_features("MASS"): - p.mass = 0.5 - if espressomd.has_features("ROTATION"): - p.rotation = [1, 1, 1] - # Make sure rinertia does not change diff coeff - if espressomd.has_features("ROTATIONAL_INERTIA"): - p.rinertia = [0.4, 0.4, 0.4] + N, kT, gamma2, loops, v_minmax, bins, error_tol) def test_msd_global_temp(self): """Tests diffusion via MSD for global gamma and temperature""" diff --git a/testsuite/python/drude.py b/testsuite/python/drude.py index bcd41de687b..627cb3ecb57 100644 --- a/testsuite/python/drude.py +++ b/testsuite/python/drude.py @@ -25,7 +25,7 @@ @utx.skipIfMissingFeatures(["P3M", "ELECTROSTATICS", "THOLE", - "LANGEVIN_PER_PARTICLE", "MASS"]) + "THERMOSTAT_PER_PARTICLE", "MASS"]) class Drude(ut.TestCase): def test(self): diff --git a/testsuite/python/langevin_thermostat.py b/testsuite/python/langevin_thermostat.py index e1e033053ae..2129031ffe8 100644 --- a/testsuite/python/langevin_thermostat.py +++ b/testsuite/python/langevin_thermostat.py @@ -40,7 +40,7 @@ def tearDown(self): self.system.thermostat.turn_off() self.system.integrator.set_vv() - def check_rng(self, per_particle_gamma=False, per_particle_temp=False): + def check_rng(self, per_particle_gamma=False): """Test for RNG consistency.""" kT = 1.1 @@ -51,17 +51,14 @@ def reset_particle(): p = system.part.add(pos=[0, 0, 0]) if espressomd.has_features("ROTATION"): p.rotation = [1, 1, 1] - if per_particle_gamma or per_particle_temp: - assert espressomd.has_features("LANGEVIN_PER_PARTICLE") if per_particle_gamma: + assert espressomd.has_features("THERMOSTAT_PER_PARTICLE") if espressomd.has_features("PARTICLE_ANISOTROPY"): p.gamma = 3 * [gamma / 2] else: p.gamma = gamma / 2 if espressomd.has_features("ROTATION"): p.gamma_rot = p.gamma * 1.5 - if per_particle_temp: - p.temp = 2 * kT return p system = self.system @@ -140,12 +137,10 @@ def test_01__rng(self): self.system.thermostat.set_langevin(kT=1, gamma=2) self.check_rng() - @utx.skipIfMissingFeatures("LANGEVIN_PER_PARTICLE") + @utx.skipIfMissingFeatures("THERMOSTAT_PER_PARTICLE") def test_01__rng_per_particle(self): """Test for RNG consistency.""" - self.check_rng(False, True) - self.check_rng(True, False) - self.check_rng(True, True) + self.check_rng(True) def test_02__friction_trans(self): """Tests the translational friction-only part of the thermostat.""" diff --git a/testsuite/python/langevin_thermostat_stats.py b/testsuite/python/langevin_thermostat_stats.py index bfe4d865797..35ee6534346 100644 --- a/testsuite/python/langevin_thermostat_stats.py +++ b/testsuite/python/langevin_thermostat_stats.py @@ -76,7 +76,7 @@ def test_vel_dist_global_temp_initial_forces(self): """ self.check_vel_dist_global_temp(True, loops=170) - @utx.skipIfMissingFeatures("LANGEVIN_PER_PARTICLE") + @utx.skipIfMissingFeatures("THERMOSTAT_PER_PARTICLE") def test_vel_dist_per_particle(self): """Test Langevin dynamics with particle-specific kT and gamma. Covers all combinations of particle-specific gamma and temp set or not set. @@ -87,14 +87,13 @@ def test_vel_dist_per_particle(self): kT = 0.9 gamma = 3.2 gamma2 = 4.3 - kT2 = 1.5 system.thermostat.set_langevin(kT=kT, gamma=gamma, seed=41) loops = 300 v_minmax = 5 bins = 4 error_tol = 0.016 self.check_per_particle( - N, kT, kT2, gamma2, loops, v_minmax, bins, error_tol) + N, kT, gamma2, loops, v_minmax, bins, error_tol) def setup_diff_mass_rinertia(self, p): if espressomd.has_features("MASS"): @@ -140,8 +139,6 @@ def test_06__diffusion(self): gamma_rot_a = [4.2, 1, 1.2] # If we have langevin per particle: - # per particle kT - per_part_kT = 1.6 # Translation per_part_gamma = 1.63 # Rotational @@ -154,7 +151,7 @@ def test_06__diffusion(self): self.setup_diff_mass_rinertia(p_global) # particle specific gamma, kT, and both - if espressomd.has_features("LANGEVIN_PER_PARTICLE"): + if espressomd.has_features("THERMOSTAT_PER_PARTICLE"): p_gamma = system.part.add(pos=(0, 0, 0)) self.setup_diff_mass_rinertia(p_gamma) if espressomd.has_features("PARTICLE_ANISOTROPY"): @@ -166,22 +163,6 @@ def test_06__diffusion(self): if espressomd.has_features("ROTATION"): p_gamma.gamma_rot = per_part_gamma_rot_i - p_kT = system.part.add(pos=(0, 0, 0)) - self.setup_diff_mass_rinertia(p_kT) - p_kT.temp = per_part_kT - - p_both = system.part.add(pos=(0, 0, 0)) - self.setup_diff_mass_rinertia(p_both) - p_both.temp = per_part_kT - if espressomd.has_features("PARTICLE_ANISOTROPY"): - p_both.gamma = per_part_gamma, per_part_gamma, per_part_gamma - if espressomd.has_features("ROTATION"): - p_both.gamma_rot = per_part_gamma_rot_a - else: - p_both.gamma = per_part_gamma - if espressomd.has_features("ROTATION"): - p_both.gamma_rot = per_part_gamma_rot_i - # Thermostat setup if espressomd.has_features("ROTATION"): if espressomd.has_features("PARTICLE_ANISOTROPY"): @@ -205,10 +186,8 @@ def test_06__diffusion(self): corr_vel = {} corr_omega = {} all_particles = [p_global] - if espressomd.has_features("LANGEVIN_PER_PARTICLE"): + if espressomd.has_features("THERMOSTAT_PER_PARTICLE"): all_particles.append(p_gamma) - all_particles.append(p_kT) - all_particles.append(p_both) # linear vel vel_obs = ParticleVelocities(ids=system.part[:].id) @@ -239,11 +218,8 @@ def test_06__diffusion(self): gamma = np.ones(3) * gamma per_part_gamma = np.ones(3) * per_part_gamma self.verify_diffusion(p_global, corr_vel, kT, gamma) - if espressomd.has_features("LANGEVIN_PER_PARTICLE"): + if espressomd.has_features("THERMOSTAT_PER_PARTICLE"): self.verify_diffusion(p_gamma, corr_vel, kT, per_part_gamma) - self.verify_diffusion(p_kT, corr_vel, per_part_kT, gamma) - self.verify_diffusion( - p_both, corr_vel, per_part_kT, per_part_gamma) # Rotation if espressomd.has_features("ROTATION"): @@ -258,13 +234,9 @@ def test_06__diffusion(self): eff_per_part_gamma_rot = per_part_gamma_rot_i * np.ones(3) self.verify_diffusion(p_global, corr_omega, kT, eff_gamma_rot) - if espressomd.has_features("LANGEVIN_PER_PARTICLE"): + if espressomd.has_features("THERMOSTAT_PER_PARTICLE"): self.verify_diffusion( p_gamma, corr_omega, kT, eff_per_part_gamma_rot) - self.verify_diffusion( - p_kT, corr_omega, per_part_kT, eff_gamma_rot) - self.verify_diffusion(p_both, corr_omega, - per_part_kT, eff_per_part_gamma_rot) def test_08__noise_correlation(self): """Checks that the Langevin noise is uncorrelated""" diff --git a/testsuite/python/mass-and-rinertia_per_particle.py b/testsuite/python/mass-and-rinertia_per_particle.py index 38bd60a07f2..42c4e986a8f 100644 --- a/testsuite/python/mass-and-rinertia_per_particle.py +++ b/testsuite/python/mass-and-rinertia_per_particle.py @@ -23,7 +23,7 @@ @utx.skipIfMissingFeatures(["MASS", "PARTICLE_ANISOTROPY", - "ROTATIONAL_INERTIA", "LANGEVIN_PER_PARTICLE"]) + "ROTATIONAL_INERTIA", "THERMOSTAT_PER_PARTICLE"]) class ThermoTest(ut.TestCase): longMessage = True # Handle for espresso system @@ -235,7 +235,7 @@ def check_dissipation(self, n): """ - tol = 1.2E-3 + tol = 1.3E-3 for _ in range(100): self.system.integrator.run(2) for i in range(n): @@ -384,23 +384,6 @@ def set_particle_specific_gamma(self, n): if espressomd.has_features("ROTATION"): self.system.part[ind].gamma_rot = self.gamma_rot_p[k, :] - def set_particle_specific_temperature(self, n): - """ - Set the particle-specific temperature. - - Parameters - ---------- - n : :obj:`int` - Number of particles of the each type. There are 2 types. - - """ - - for k in range(2): - self.halfkT_p_validate[k] = self.kT_p[k] / 2.0 - for i in range(n): - ind = i + k * n - self.system.part[ind].temp = self.kT_p[k] - def set_diffusivity_tran(self): """ Set the translational diffusivity to validate further. @@ -412,7 +395,7 @@ def set_diffusivity_tran(self): self.D_tran_p_validate[k, :] = 2.0 * \ self.halfkT_p_validate[k] / self.gamma_tran_p_validate[k, :] - # Test case 0.0: no particle specific values / dissipation only + # Test case 0.0: no particle specific gamma / dissipation only def test_case_00(self): # Each of 2 kind of particles will be represented by n instances: n = 1 @@ -424,7 +407,7 @@ def test_case_00(self): # Actual integration and validation run self.check_dissipation(n) - # Test case 0.1: no particle specific values / fluctuation & dissipation + # Test case 0.1: no particle specific gamma / fluctuation & dissipation def test_case_01(self): # Each of 2 kind of particles will be represented by n instances: n = 500 @@ -437,7 +420,7 @@ def test_case_01(self): # Actual integration and validation run self.check_fluctuation_dissipation(n) - # Test case 1.0: particle specific gamma but not temperature / dissipation + # Test case 1.0: particle specific gamma / dissipation # only def test_case_10(self): # Each of 2 kind of particles will be represented by n instances: @@ -451,7 +434,7 @@ def test_case_10(self): # Actual integration and validation run self.check_dissipation(n) - # Test case 1.1: particle specific gamma but not temperature / fluctuation + # Test case 1.1: particle specific gamma / fluctuation # & dissipation def test_case_11(self): # Each of 2 kind of particles will be represented by n instances: @@ -466,67 +449,7 @@ def test_case_11(self): # Actual integration and validation run self.check_fluctuation_dissipation(n) - # Test case 2.0: particle specific temperature but not gamma / dissipation - # only - def test_case_20(self): - # Each of 2 kind of particles will be represented by n instances: - n = 1 - self.dissipation_param_setup(n) - self.set_langevin_global_defaults() - # The test case-specific thermostat and per-particle parameters - self.system.thermostat.set_langevin( - kT=self.kT, gamma=self.gamma_global, seed=42) - self.set_particle_specific_temperature(n) - # Actual integration and validation run - self.check_dissipation(n) - - # Test case 2.1: particle specific temperature but not gamma / fluctuation - # & dissipation - def test_case_21(self): - # Each of 2 kind of particles will be represented by n instances: - n = 500 - self.fluctuation_dissipation_param_setup(n) - self.set_langevin_global_defaults() - # The test case-specific thermostat and per-particle parameters - self.system.thermostat.set_langevin( - kT=self.kT, gamma=self.gamma_global, seed=42) - self.set_particle_specific_temperature(n) - self.set_diffusivity_tran() - # Actual integration and validation run - self.check_fluctuation_dissipation(n) - - # Test case 3.0: both particle specific gamma and temperature / - # dissipation only - def test_case_30(self): - # Each of 2 kind of particles will be represented by n instances: - n = 1 - self.dissipation_param_setup(n) - self.set_langevin_global_defaults() - # The test case-specific thermostat and per-particle parameters - self.system.thermostat.set_langevin( - kT=self.kT, gamma=self.gamma_global, seed=42) - self.set_particle_specific_gamma(n) - self.set_particle_specific_temperature(n) - # Actual integration and validation run - self.check_dissipation(n) - - # Test case 3.1: both particle specific gamma and temperature / - # fluctuation & dissipation - def test_case_31(self): - # Each of 2 kind of particles will be represented by n instances: - n = 500 - self.fluctuation_dissipation_param_setup(n) - self.set_langevin_global_defaults() - # The test case-specific thermostat and per-particle parameters - self.system.thermostat.set_langevin( - kT=self.kT, gamma=self.gamma_global, seed=42) - self.set_particle_specific_gamma(n) - self.set_particle_specific_temperature(n) - self.set_diffusivity_tran() - # Actual integration and validation run - self.check_fluctuation_dissipation(n) - - # Test case 4.0: no particle specific values / rotational specific global + # Test case 4.0: no particle specific gamma / rotational specific global # thermostat / dissipation only def test_case_40(self): # Each of 2 kind of particles will be represented by n instances: @@ -542,7 +465,7 @@ def test_case_40(self): # Actual integration and validation run self.check_dissipation(n) - # Test case 4.1: no particle specific values / rotational specific global + # Test case 4.1: no particle specific gamma / rotational specific global # thermostat / fluctuation & dissipation def test_case_41(self): # Each of 2 kind of particles will be represented by n instances: diff --git a/testsuite/python/particle.py b/testsuite/python/particle.py index 243e6d76d10..320c823557e 100644 --- a/testsuite/python/particle.py +++ b/testsuite/python/particle.py @@ -135,7 +135,7 @@ def test_director(self): sample_vector_normalized ) - if espressomd.has_features(["LANGEVIN_PER_PARTICLE"]): + if espressomd.has_features(["THERMOSTAT_PER_PARTICLE"]): if espressomd.has_features(["PARTICLE_ANISOTROPY"]): test_gamma = generateTestForVectorProperty( "gamma", np.array([2., 9., 0.23])) diff --git a/testsuite/python/thermostats_common.py b/testsuite/python/thermostats_common.py index ea8c35d1c14..6df58369616 100644 --- a/testsuite/python/thermostats_common.py +++ b/testsuite/python/thermostats_common.py @@ -88,6 +88,15 @@ def check_global(self, N, kT, steps, v_minmax, # Warmup system.integrator.run(20) + vel_obs = ParticleVelocities(ids=system.part[:].id) + vel_acc = TimeSeries(obs=vel_obs) + system.auto_update_accumulators.add(vel_acc) + + if espressomd.has_features("ROTATION"): + omega_obs = ParticleBodyAngularVelocities(ids=system.part[:].id) + omega_acc = TimeSeries(obs=omega_obs) + system.auto_update_accumulators.add(omega_acc) + # Sampling v_stored = np.zeros((steps, N, 3)) omega_stored = np.zeros((steps, N, 3)) @@ -97,14 +106,15 @@ def check_global(self, N, kT, steps, v_minmax, if espressomd.has_features("ROTATION"): omega_stored[i] = system.part[:].omega_body - self.check_velocity_distribution( - v_stored.reshape((-1, 3)), v_minmax, n_bins, error_tol, kT) + vel = vel_acc.time_series().reshape((-1, 3)) + self.check_velocity_distribution(vel, v_minmax, n_bins, error_tol, kT) if espressomd.has_features("ROTATION"): + omega = omega_acc.time_series().reshape((-1, 3)) self.check_velocity_distribution( - omega_stored.reshape((-1, 3)), v_minmax, n_bins, error_tol, kT) + omega, v_minmax, n_bins, error_tol, kT) - def check_per_particle(self, N, kT_global, kT_local, gamma_local, - steps, v_minmax, n_bins, error_tol): + def check_per_particle(self, N, kT, gamma_local, steps, v_minmax, n_bins, + error_tol): """ Test Langevin/Brownian dynamics velocity distribution with global and particle-specific kT/gamma. Covers all combinations of particle @@ -114,10 +124,8 @@ def check_per_particle(self, N, kT_global, kT_local, gamma_local, ---------- N : :obj:`int` Number of particles. - kT_global : :obj:`float` + kT : :obj:`float` Global temperature. - kT_local : :obj:`float` - Per-particle temperature. gamma_local : :obj:`float` Per-particle gamma. steps : :obj:`int` @@ -137,45 +145,28 @@ def check_per_particle(self, N, kT_global, kT_local, gamma_local, if espressomd.has_features("PARTICLE_ANISOTROPY"): gamma_local = 3 * [gamma_local] - # Set different kT on 2nd half of particles - system.part[N // 2:].temp = kT_local - # Set different gamma on half of the particles (overlap over both kTs) - system.part[N // 4:3 * N // 4].gamma = gamma_local + # Set different gamma on 2nd half of particles + system.part[N // 2:].gamma = gamma_local system.integrator.run(50) - vel_obs_global = ParticleVelocities(ids=system.part[:N // 2].id) - vel_obs_local = ParticleVelocities(ids=system.part[N // 2:].id) - vel_acc_global = TimeSeries(obs=vel_obs_global) - vel_acc_local = TimeSeries(obs=vel_obs_local) - system.auto_update_accumulators.add(vel_acc_global) - system.auto_update_accumulators.add(vel_acc_local) + vel_obs = ParticleVelocities(ids=system.part[:].id) + vel_acc = TimeSeries(obs=vel_obs) + system.auto_update_accumulators.add(vel_acc) if espressomd.has_features("ROTATION"): - omega_obs_global = ParticleBodyAngularVelocities( - ids=system.part[:N // 2].id) - omega_obs_local = ParticleBodyAngularVelocities( - ids=system.part[N // 2:].id) - omega_acc_global = TimeSeries(obs=omega_obs_global) - omega_acc_local = TimeSeries(obs=omega_obs_local) - system.auto_update_accumulators.add(omega_acc_global) - system.auto_update_accumulators.add(omega_acc_local) + omega_obs = ParticleBodyAngularVelocities(ids=system.part[:].id) + omega_acc = TimeSeries(obs=omega_obs) + system.auto_update_accumulators.add(omega_acc) system.integrator.run(steps) - vel_global = vel_acc_global.time_series().reshape((-1, 3)) - vel_local = vel_acc_local.time_series().reshape((-1, 3)) - self.check_velocity_distribution( - vel_global, v_minmax, n_bins, error_tol, kT_global) - self.check_velocity_distribution( - vel_local, v_minmax, n_bins, error_tol, kT_local) + vel = vel_acc.time_series().reshape((-1, 3)) + self.check_velocity_distribution(vel, v_minmax, n_bins, error_tol, kT) if espressomd.has_features("ROTATION"): - omega_global = omega_acc_global.time_series().reshape((-1, 3)) - omega_local = omega_acc_local.time_series().reshape((-1, 3)) - self.check_velocity_distribution( - omega_global, v_minmax, n_bins, error_tol, kT_global) + omega = omega_acc.time_series().reshape((-1, 3)) self.check_velocity_distribution( - omega_local, v_minmax, n_bins, error_tol, kT_local) + omega, v_minmax, n_bins, error_tol, kT) def check_noise_correlation(self, kT, steps, delta): """Test the Langevin/Brownian noise is uncorrelated. diff --git a/testsuite/python/virtual_sites_tracers.py b/testsuite/python/virtual_sites_tracers.py index c18bd798cf8..5d3d3786a6a 100644 --- a/testsuite/python/virtual_sites_tracers.py +++ b/testsuite/python/virtual_sites_tracers.py @@ -23,7 +23,8 @@ from virtual_sites_tracers_common import VirtualSitesTracersCommon -@utx.skipIfMissingFeatures(['VIRTUAL_SITES_INERTIALESS_TRACERS']) +@utx.skipIfMissingFeatures( + ['VIRTUAL_SITES_INERTIALESS_TRACERS', 'LB_BOUNDARIES']) class VirtualSitesTracers(ut.TestCase, VirtualSitesTracersCommon): def setUp(self): diff --git a/testsuite/python/virtual_sites_tracers_gpu.py b/testsuite/python/virtual_sites_tracers_gpu.py index 1e982c7904f..7b00383bca2 100644 --- a/testsuite/python/virtual_sites_tracers_gpu.py +++ b/testsuite/python/virtual_sites_tracers_gpu.py @@ -24,7 +24,8 @@ @utx.skipIfMissingGPU() -@utx.skipIfMissingFeatures(['VIRTUAL_SITES_INERTIALESS_TRACERS']) +@utx.skipIfMissingFeatures( + ['VIRTUAL_SITES_INERTIALESS_TRACERS', 'LB_BOUNDARIES']) class VirtualSitesTracers(ut.TestCase, VirtualSitesTracersCommon): def setUp(self):