#include "Test.h"

#ifdef MM_WIN32
#pragma warning (disable : 4786)
#endif

#include "MolMeccano_file.h"

#include "Log.h"
#include "User.h"

#include "Timer.h"
#include "Atom.h"
#include "Model.h"
#include "Model_kit.h"
#include "Path.h"
#include "Create.h"
#include "Force_field.h"
#include "Common_interactions.h"
#include "Conventional_interactions.h"
#include "List_interactions.h"

namespace MM
{

class AMBER94_interactions_test : public Test
{
    Common_interactions *link;

    Force_field     *force_field_;

public:
    AMBER94_interactions_test (const Text & suite_name) 
        : Test (suite_name) 
    {
    }
    
    void run()
    {
        //force_field_ = new Force_field ("AMBER94");//fix
        
        test_20_acids         ();       count_tests();

        //delete force_field_;
    }

    void test_20_acids ()
    {
        MolMeccano_file::apply_torsions (false);

        Own <Model>     model (prototype <Model> (). clone ());
        model().kit().set_force_field ("AMBER94");

        //Own <Common_interactions>   janitor (Common_interactions::
        //    create(model(), model(), *force_field_));
        Own <Common_interactions>   janitor 
            (new Conventional_interactions(model()/*, model(), *force_field_*/));
            //(new List_interactions(model()));
        Common_interactions & interactions = janitor();
        
        double E;
                                                        // Radzicka & Wolfenden
        model(). load (Path::aminoacids() + "Ala.mlm" );//  1.81 Alanine 
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  1.7973, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  1.9549, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4916, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  2.3116, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-31.0889, 0.1);

        model().load (Path::aminoacids() + "Arg+.mlm");
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.7450, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.2103, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4916, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  5.4251, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-154.0290, 0.1);

        model().load (Path::aminoacids() + "Asn.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  1.8521, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.9589, 0.0001);
        //DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.7953, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (), 10.1467, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (), 72.2839, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-99.9042, 0.1);

        model().load (Path::aminoacids() + "Asp-.mlm");
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  1.7069, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.4794, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.5024, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (), 39.7277, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-40.5903, 0.1);

        model().load (Path::aminoacids() + "Cys.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  1.6893, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.3669, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.5024, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  2.8744, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-29.0021, 0.1);

        model().load (Path::aminoacids() + "Gln.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.1849, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.9601, 0.0001);
        //DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.7953, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (), 10.1468, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (), 13.3893, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-94.9695, 0.1);

        model().load (Path::aminoacids() + "Glu-.mlm");
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.0397, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.4806, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.5025, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  7.2342, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-27.6643, 0.1);

        model().load (Path::aminoacids() + "Gly.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  1.4645, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  1.7753, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  3.6303, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  2.3467, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-36.8732, 0.1);

        model().load (Path::aminoacids() + "His+.mlm");
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  3.2492, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (), 21.7742, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4672, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  0.9239, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (), 18.7911, 0.1);

        model().load (Path::aminoacids() + "Ile.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.7955, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  1.8302, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4675, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (), 14.1843, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-38.8969, 0.1);

        model().load (Path::aminoacids() + "Leu.mlm" );
        E = interactions.potential ();
    
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.7956, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  1.9586, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4917, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),1729.1700,0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-58.0901, 0.1);

        model().load (Path::aminoacids() + "Lys+.mlm");
        E = interactions.potential ();
 
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  3.0984, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.0641, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4929, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  4.2341, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (), 10.2367, 0.1);

        model().load (Path::aminoacids() + "Met.mlm" );
        E = interactions.potential ();
 
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.4474, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  3.9504, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.5129, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  3.3083, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-34.2456, 0.1);

        model().load (Path::aminoacids() + "Phe.mlm" );
        E = interactions.potential ();
  
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.7084, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.3429, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4916, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  6.0444, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-35.4605, 0.1);

        model().load (Path::aminoacids() + "Pro.mlm" );
        E = interactions.potential ();
  
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.5361, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (), 17.8519, 0.0001);
        //DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  1.6151, 0.0001);

        // x-CT-CT-x tinker 0.156     AMBER 0.1555555
        DOUBLES_EQUAL (interactions.last_torsion_potential          (), 21.9335, 0.03);

        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (), 86.6211, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-22.5618, 0.1);

        model().load (Path::aminoacids() + "Ser.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  1.8825, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.8516, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.5208, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  2.6805, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-38.2704, 0.1);
 
        model().load (Path::aminoacids() + "Thr.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.2152, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.6655, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);

        // x-CT-CT-x tinker 0.156     AMBER 0.1555555
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.5096, 0.001);
        
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  8.0991, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-61.4639, 0.1);

        model().load (Path::aminoacids() + "Trp.mlm" );
        E = interactions.potential ();

        DOUBLES_EQUAL (interactions.last_bond_potential             (), 16.6909, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  5.0560, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4672, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  6.3378, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-27.2378, 0.1);

        model().load (Path::aminoacids() + "Tyr.mlm" );
        E = interactions.potential ();
  
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.8689, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  2.6094, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4916, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  5.9395, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-53.2062, 0.1);

        model().load (Path::aminoacids() + "Val.mlm" );
        E = interactions.potential ();
 
        DOUBLES_EQUAL (interactions.last_bond_potential             (),  2.4628, 0.0001);
        DOUBLES_EQUAL (interactions.last_angle_potential            (),  1.9570, 0.0001);
        DOUBLES_EQUAL (interactions.last_improper_torsion_potential (),  0.0000, 0.0001);
        DOUBLES_EQUAL (interactions.last_torsion_potential          (),  4.4916, 0.0001);
        DOUBLES_EQUAL (interactions.last_Van_der_Waals_potential    (),  6.0575, 0.0001);
        DOUBLES_EQUAL (interactions.last_electrostatic_potential    (),-55.2887, 0.1);

        MolMeccano_file::apply_torsions (true);
    }
};

AMBER94_interactions_test  test_AMBER94_interactions ("correctness");

}//MM
