Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/atom_trajectoryparticle.cpp

    ra3fded r112b09  
    55 *      Author: heber
    66 */
     7
     8#include "Helpers/MemDebug.hpp"
    79
    810#include "atom.hpp"
     
    1315#include "log.hpp"
    1416#include "parser.hpp"
    15 #include "ThermoStatContainer.hpp"
    1617#include "verbose.hpp"
    1718
     
    196197void TrajectoryParticle::Thermostat_Langevin(int Step, gsl_rng * r, double *ekin, config *configuration)
    197198{
    198   double sigma  = sqrt(configuration->Thermostats->TargetTemp/type->mass); // sigma = (k_b T)/m (Hartree/atomicmass = atomiclength/atomictime)
     199  double sigma  = sqrt(configuration->TargetTemp/type->mass); // sigma = (k_b T)/m (Hartree/atomicmass = atomiclength/atomictime)
    199200  Vector &U = Trajectory.U.at(Step);
    200201  if (FixedIon == 0) { // even FixedIon moves, only not by other's forces
    201202    // throw a dice to determine whether it gets hit by a heat bath particle
    202     if (((((rand()/(double)RAND_MAX))*configuration->Thermostats->TempFrequency) < 1.)) {
     203    if (((((rand()/(double)RAND_MAX))*configuration->TempFrequency) < 1.)) {
    203204      DoLog(3) && (Log() << Verbose(3) << "Particle " << *this << " was hit (sigma " << sigma << "): " << sqrt(U[0]*U[0]+U[1]*U[1]+U[2]*U[2]) << " -> ");
    204205      // pick three random numbers from a Boltzmann distribution around the desired temperature T for each momenta axis
     
    224225  if (FixedIon == 0) { // even FixedIon moves, only not by other's forces
    225226    for (int d=0; d<NDIM; d++) {
    226       U[d] *= sqrt(1+(configuration->Deltat/configuration->Thermostats->TempFrequency)*(ScaleTempFactor-1));
     227      U[d] *= sqrt(1+(configuration->Deltat/configuration->TempFrequency)*(ScaleTempFactor-1));
    227228      *ekin += 0.5*type->mass * U[d]*U[d];
    228229    }
     
    254255  if (FixedIon == 0) { // even FixedIon moves, only not by other's forces
    255256    for (int d=0; d<NDIM; d++) {
    256         U[d] += configuration->Deltat/type->mass * (configuration->Thermostats->alpha * (U[d] * type->mass));
     257        U[d] += configuration->Deltat/type->mass * (configuration->alpha * (U[d] * type->mass));
    257258        *ekin += (0.5*type->mass) * U[d]*U[d];
    258259      }
Note: See TracChangeset for help on using the changeset viewer.