#ifndef MTS_MD_H
#define MTS_MD_H

#ifndef MTS_H
#include "MTS.h"
#endif

#ifndef MD_H
#include "MD.h"
#endif

namespace MM
{

class MTS_MD : public MTS
{
    double      gamma_;
    double      R2_;     // A

public:
    explicit MTS_MD (Prototype const &) {adopt_prototype (this);}
    explicit MTS_MD () : gamma_(0.0048), R2_(8.) {inner_steps_ = 20;}

//    explicit MTS_MD (MD_subject & subject) : MD_method (subject) { }

    virtual void adopt_inner_step (){}

protected:
    /*virtual Text const& method () const
    {
        static const Text name = MD::method() 
            + ": Multiple time step";
        return name;
    }//*/

    char const*     class_name () const {return "MTS_MD";}
    char const*     title      () const {return "Multiple time step";}
    char const*     alias      () const {return "MTS";}

    MD_method *     clone () const      {return new MTS_MD;}

    virtual void start ()
    {
        Lock<MD> lock(simulation());
        MD_method::start ();
        //F (MD_subject::Whole);
        //F (MD_subject::Long3, 1000000.);
        F (MD_subject::Long2, 1000000., S_);
        F (MD_subject::Long, R2_, S_);
        F (MD_subject::Short, 1000, S_);
    }


    virtual void step ()
    {
        //int last_step = inner_steps_-1;

        //drift (gamma_ * .25);
        //F (MD_subject::Long2, 1000000.);
        //F (MD_subject::Long, R2_);
        //kick  (.5, MD_subject::Long2);
        //drift (gamma_ * .25);

        kick  (.5, MD_subject::Long2);
        drift (gamma_ * .5);

        simulation_->algorithm().description ("{");

        for (int i=0;  i<inner_steps_;  ++i)
        {
            drift ((1-gamma_) * .5 / inner_steps_);
            F (MD_subject::Long, R2_, S_);
            F (MD_subject::Short, 1000, S_);
            kick  (            1. / inner_steps_, MD_subject::Long);
            kick  (            1. / inner_steps_, MD_subject::Short);
            drift ((1-gamma_) * .5 / inner_steps_);

            //kick  (            .5 / inner_steps_, MD_subject::Long);
            //kick  (            .5 / inner_steps_, MD_subject::Short);
            //drift ((1-gamma_) * 1. / inner_steps_);
            //F (MD_subject::Long, R2_);
            //F (MD_subject::Short, 1000);
            //kick  (            .5 / inner_steps_, MD_subject::Long);
            //kick  (            .5 / inner_steps_, MD_subject::Short);

            if (i<inner_steps_-1)
                simulation_->algorithm().description ("");
        }

        simulation_->algorithm().description ("}");

        drift (gamma_ * .5); 
        //F (MD_subject::Whole);
        F (MD_subject::Long2, 1000000., S_);
        F (MD_subject::Long, R2_, S_);
        kick  (.5, MD_subject::Long2);

        //drift (gamma_ * .25);
        ////F (MD_subject::Whole);
        //F (MD_subject::Long2, 1000000.);
        //F (MD_subject::Long, R2_);
        //kick  (.5, MD_subject::Long2);
        //drift (gamma_ * .25);
    }


    //virtual void step ()
    //{
    //    drift (gamma_ * .25);
    //    F (MD_subject::Long);
    //    kick  (.5, MD_subject::Long);
    //    drift (gamma_ * .25);

    //    simulation_->algorithm().description ("{");

    //    for (int i=0;  i<inner_steps_;  ++i)
    //    {
    //        drift ((1-gamma_) * .5 / inner_steps_);
    //        F (MD_subject::Short);
    //        kick  (            1. / inner_steps_, MD_subject::Short);
    //        drift ((1-gamma_) * .5 / inner_steps_);

    //        //kick  (            .5 / inner_steps_, MD_subject::Short);
    //        //drift ((1-gamma_) * 1. / inner_steps_);
    //        //F (MD_subject::Short);
    //        //kick  (            .5 / inner_steps_, MD_subject::Short);

    //        if (i<inner_steps_-1)
    //            simulation_->algorithm().description ("");
    //    }

    //    simulation_->algorithm().description ("}");

    //    drift (gamma_ * .25); 
    //    //F (MD_subject::Whole);
    //    F (MD_subject::Long);
    //    kick  (.5, MD_subject::Long);
    //    drift (gamma_ * .25);
    //}

};

}//MM

#endif //MTS_MD_H
