#ifndef ISOKINETIC_TERMOSTAT_H
#define ISOKINETIC_TERMOSTAT_H

#ifndef THERMOSTAT_H
#include "Thermostat.h"
#endif

#ifndef ATOM_KIT_H
#include "Atom_kit.h"
#endif

#ifndef HYDROGEN_ATOMS_H
#include "Hydrogen_atoms.h"
#endif

#ifndef  MACH_EPS_H
#include "Mach_eps.h"
#endif

namespace MM
{

class Isokinetic_termostat : public Thermostat
{
public:
    Isokinetic_termostat (Model & model, Atom_group & atom_group, MD & simulation)
    :   Thermostat (model, atom_group, simulation)    { }

    char const* title () const      {return "Isokinetic";}

    virtual void    handle_execute ()
    {
        MD_simulation_.algorithm().description ("Isokinetic thermostat");
        //Hydrogen_atoms group (atom_group_);
        Atom_group & group = atom_group_;

        double T   = calc_temperature (group);   // actual temperature
        int    count = group.atom_count();

        double lambda = sqrt (T0_ / (T+d_mach_eps_2));

        for (int i=0;  i<count;  ++i)
        {
            Atom & atom = group.atom (i);
            atom.kit().scale_velocity (lambda);
        }
    }

    double              response_constant () const     {return 0.;}
    void                set_response_constant (double) {;}
};

}//MM

#endif //ISOKINETIC_TERMOSTAT_H
