Ignore:
Timestamp:
Sep 9, 2008, 8:37:07 AM (17 years ago)
Author:
Frederik Heber <heber@…>
Children:
8e52aa
Parents:
ec4e84
Message:

new function: molecule::OutputTemperatureFromTrajectories() stores velocities as temperature

In order to use VibrAlyzer we need a temperature file, which we don't have anymore as we do the M
D in molecuilder now. Hence, there is a new command line switch (-L) and a menu item in MeasureAt
oms() in order to print temperature or to store it to a given file.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • molecuilder/src/molecules.cpp

    rec4e84 r79cd79  
    46644664};
    46654665
     4666/** Stores the temperature evaluated from velocities in molecule::Trajectories.
     4667 * We simply use the formula equivaleting temperature and kinetic energy:
     4668 * \f$k_B T = \sum_i m_i v_i^2\f$
     4669 * \param *out output stream for debugging
     4670 * \param startstep first MD step in molecule::Trajectories
     4671 * \param endstep last plus one MD step in molecule::Trajectories
     4672 * \param *output output stream of temperature file
     4673 * \return file written (true), failure on writing file (false)
     4674 */
     4675bool molecule::OutputTemperatureFromTrajectories(ofstream *out, int startstep, int endstep, ofstream *output)
     4676{
     4677  double temperature;
     4678  atom *Walker = NULL;
     4679  // test stream
     4680  if (output == NULL)
     4681    return false;
     4682  else
     4683    *output << "# Step Temperature [K] Temperature [a.u.]" << endl;
     4684  for (int step=startstep;step < endstep; step++) { // loop over all time steps
     4685    temperature = 0.;
     4686    Walker = start;
     4687    while (Walker->next != end) {
     4688      Walker = Walker->next;
     4689      for (int i=NDIM;i--;)
     4690        temperature += Walker->type->mass * Trajectories[Walker].U.at(step).x[i]* Trajectories[Walker].U.at(step).x[i];
     4691    }
     4692    *output << step << "\t" << temperature*AtomicEnergyToKelvin << "\t" << temperature << endl;
     4693  }
     4694  return true;
     4695};
Note: See TracChangeset for help on using the changeset viewer.