#include "MD.h"

#include "Defs.h"
#include "Coordinates.h"
#include "Velocities.h"
#include "Force.h"
#include "Mass.h"
#include "MD_method.h"
#include "Thermostat.h"
#include "Temperature_unit.h"
#include "Barostat.h"
//#include "MD_subject.h"
#include "Lock.h"

#include "Log.h"
#include "User.h"
#include "Profiler.h"
#include "Timer.h"

#include <math.h>

namespace MM
{

const int max_descriptions = 2;//fix

MD::
MD ()
:   
    step_length_(1e-015 * Constant::sec_to_time_unit),
    drift_length_(0.),
    kick_length_(0.),
    subject_ (0),
    method_(0),
    thermostat_(&Thermostat::none()),
    barostat_  (&Barostat::none()),
    stop_flag_(false)
{
}

MD::
MD (MD_subject & subject, MD_method & method)
:   
    step_length_(1e-015 * Constant::sec_to_time_unit),
    drift_length_(0.),
    kick_length_(0.),
    subject_ (&subject),
    method_(&method),
    thermostat_(&Thermostat::none()),
    barostat_  (&Barostat::none()),
    stop_flag_(false)
{
    method_->set_simulation (*this);
}


Text MD::
title () const
{
    return Text() + "MD: " + method_->title();
}

bool MD::
locked () const
{
    return subject_->locked ();
}

void MD::
set_lock (bool on)
{
    subject_->set_lock (on);
}

double MD::
step_length () const
{
    return step_length_;
}

void MD::
set_step_length (double new_length)
{
    step_length_ = new_length;
}

void MD::
set_method (MD_method & new_method)
{
    method_ = &new_method;
    method_->set_simulation (*this);
}

void MD::
set_subject (MD_subject & new_subject)
{
    subject_ = &new_subject;
}

void MD::
setup_velocities ()
{
    subject_->setup_velocities (thermostat_->desired_temperature ());
}

void MD::
setup_velocities (Temperature_unit const& t)
{
    subject_->setup_velocities (t);
}

void MD::
start ()
{
    if (algorithm().is_on())
        algorithm().description (Text("start  ") + title());

    stop_flag_ = false;
    subject_->start();
    method_->start();
}

void MD::
cycle ()
{
    Text message("MD::cycle  ");
    if (algorithm().is_on())
        algorithm().description (Text("begin ") + iterations() + " step cycle");

    //Lock<MD> lock(*this);

    Timer timer;    Text  time;
    int steps = iterations();

    set_lock (true);

    for (int i=0;  i<steps && !stop_flag_;  ++i)
    {
        algorithm().set_on (i < max_descriptions);
        if (algorithm().is_on())
        {
            algorithm().description ("");
            algorithm().description (Text(i) + " step");
            algorithm().description ("");
        }

        profile().interactions___add_force___start ();
        try
        {
            method_->step ();
        }
        catch(...) {to_user().fatal_error (message + " in step"); throw;}
        //profile().interactions___add_force___finish ();

        try
        {
            thermostat_->execute(); 
        }
        catch(...) {to_user().fatal_error (message + " in thermostat"); throw;}

        try
        {
            barostat_->execute(); 
        }
        catch(...) {to_user().fatal_error (message + " in barostat"); throw;}
        
        try
        {
            step_finished_.publish ();  
        }
        catch(...) {to_user().fatal_error (message + " in step_finished"); throw;}

        //if (timer.elapsed() > .05)
        if (timer.elapsed() > .08)
        {
            try
            {
                set_lock (false);
                update ();                      //fix
                set_lock (true);
            }
            catch(...) {to_user().fatal_error (message + " in update"); throw;}

            timer.restart();
        }

        //if (i >= max_descriptions)
        //    algorithm().off();

        //time += "iteration ";
        //time += timer.elapsed();    
        //time += " sec\n";
        //log () << time;
        profile().interactions___add_force___finish ();
    }
    set_lock (false);

    try
    {
        update (true);                      
    }
    catch(...) {to_user().fatal_error (message + " in update"); throw;}

    if (algorithm().is_on())
        algorithm().description ("end cycle");
}

int MD::
iterations () const
{
    return method_->iterations ();
}

void MD::
finish ()
{
    method_->finish ();
    subject_->finish ();
    stop_.publish ();
    algorithm().set_on (true);
    algorithm().description (Text("finish ") + title());
    ENSURE ("", fabs(drift_length_ - kick_length_) <= 0.000001 * drift_length_);
}

void MD::
zero_virial ()
{
    subject_->zero_virial ();
    if (algorithm().is_on())
        algorithm().description ("zero virial");
}

void MD::
list (Hint hint, double skin)
{
    //Timer timer;    Text  time;

    subject_->list (hint, skin);
    if (algorithm().is_on())
        algorithm().description (Text("list       ") + hint + " skin="+ skin);

    //time += "list         ";
    //time += timer.elapsed();    timer.restart();
    //time += " sec\n";
    //log () << time;
}
void MD::
list (Hint hint, double R, double skin)
{
    //Timer timer;    Text  time;

    subject_->list (hint, R, skin);
    if (algorithm().is_on())
        algorithm().description
        (Text("list       ") + hint + " R="+ R + " skin="+ skin);

    //time += "list         ";
    //time += timer.elapsed();    timer.restart();
    //time += " sec\n";
    //log () << time;
}

    
double MD::
degrees_of_freedom ()
{
    if (algorithm().is_on())
        algorithm().description (Text("degrees of freedom "));
    return subject_->degrees_of_freedom ();
}

double MD::
H ()
{
    if (algorithm().is_on())
        algorithm().description (Text("H       "));
    return subject_->H ();
}

double MD::
K ()
{
    if (algorithm().is_on())
        algorithm().description (Text("K       "));
    return subject_->K ();
}

double MD::
U () 
{
    if (algorithm().is_on())
        algorithm().description (Text("U       "));
    return subject_->U ();
}

void MD::
F (Hint hint)
{
    //Timer timer;    Text  time;

    subject_->F (hint);
    if (algorithm().is_on())
        algorithm().description (Text("F       ") + hint);

    //time += "F         ";
    //time += timer.elapsed();    timer.restart();
    //time += " sec\n";
    //log () << time;
}


void MD::
F (Hint hint, double R, double Switching)
{
    //Timer timer;    Text  time;

    subject_->F (hint, R, Switching);
    if (algorithm().is_on())
        algorithm().description
        (Text("F       ") + hint + " R="+ R + " S=" + Switching);

    //time += "F         ";
    //time += timer.elapsed();    timer.restart();
    //time += " sec\n";
    //log () << time;
}


void MD::
clear_moment (Hint hint)
{
    subject_->clear_moment (hint);
    if (algorithm().is_on())
        algorithm().description (Text("clear moment "));
}

void MD::
kick (double h, Hint hint)
{
    //Timer timer;    Text  time;

    kick_length_ += h * step_length_;
    subject_->kick (h*step_length_, hint);
    if (algorithm().is_on())
        algorithm().description
        (Text("kick    ") + h*step_length_ + "  \tsum = " + kick_length_);

    //time += "kick      ";
    //time += timer.elapsed();    timer.restart();
    //time += " sec\n";
    //log () << time;
}

void MD::
kick  (int substeps, double dT, Hint hint)
{
    kick_length_ += dT * step_length_;
    subject_->kick (substeps, dT*step_length_, hint);
    if (algorithm().is_on())
        algorithm().description
        (Text("kick    ") + dT*step_length_ +
        "  \tsum = " + kick_length_ +
        "  \tsubsteps = " + substeps );
}

void MD::
drift (double h)
{
    //Timer timer;    Text  time;

    drift_length_ += h * step_length_;
    subject_->drift (h*step_length_);

    if (algorithm().is_on())
        algorithm().description
        (Text("drift   ") + h*step_length_ + "  \tsum = " +drift_length_);

    //time += "drift     ";
    //time += timer.elapsed();    timer.restart();
    //time += " sec\n";
    //log () << time;
}

void MD::
keep_configuration ()
{
    subject_->keep_configuration ();
    if (algorithm().is_on())
        algorithm().description (Text("keep configuration"));
}

void MD::
restore_configuration ()
{
    subject_->restore_configuration ();
    if (algorithm().is_on())
        algorithm().description (Text("restore configuration"));
}

void MD::
update (bool force)
{
    subject_->update (force);
}

/*
bool MD::
invariant () const
{
    bool result = coordinats_ != 0 &&
                  velocities_ != 0 &&
                  forces_     != 0 &&
                  masses_     != 0;
    if (!result)    
        return false;

    int coordinats = coordinats_-> size ();
    int velocities = velocities_-> size ();
    int forces     = forces_->     size ();
    int masses     = masses_->     size ();

    result = coordinats == velocities &&
             coordinats == forces &&
             coordinats == masses;

    return result;
}//*/

}//MM








