#include <iostream>
#include <cmath>
#include <vector>
#include <cstdlib>
#include <tuple>
// #include <ctime>

#define KELVIN 273.15
#define R 8.314
#define ctemp 30
#define temp_delta 30

using std::cout;
using std::endl;

const int hmax = 1000;
const int wmax = 0.007;
const double mmassarg = 40*pow(10,-3);
const double massarg = 40*1.66*pow(10,-27);
const double T2 = KELVIN + ctemp;
const double T1 = T2 + temp_delta;
const double V1 = sqrt((3*R*T1)/mmassarg);
const double V2 = sqrt((3*R*T2)/mmassarg);
const double time_delta = 1*pow(10,-6);

std::tuple<double,double> getrs(bool fl){
  float x = static_cast <float> (rand()) / static_cast <float> (RAND_MAX);
  float y = static_cast <float> (rand()) / static_cast <float> (RAND_MAX);
  if (fl) {
    return std::make_tuple(x*V1,y*V1);
  }else{
    return std::make_tuple(x*V2,y*V2);
  }
}

struct molecule{
  double posx,posy;
  double speedx,speedy;
  molecule(double posx,double posy,double speedx,double speedy):posx(posx),posy(posy),speedx(speedx),speedy(speedy){}
  molecule(double posx,double posy,std::tuple<double,double> tp):posx(posx),posy(posy){
    speedx = std::get<0>(tp);
    speedy = std::get<1>(tp);
  }
};

int main(int argc, char const *argv[]) {
  cout << "App speed1: " << V1 << endl;
  cout << "App speed2: " << V2 << endl;
  molecule am(0.0045,0.0045,getrs(false));
  double time = 0;
  int count = 0;
  double energy = 0;
  double mposx,mposy;
  while (time < 1) {
    mposx = am.posx + am.speedx*time_delta;
    mposy = am.posy + am.speedy*time_delta;
    if (mposx > wmax) {
      double ltime = (wmax-am.posx)/am.speedx;
      double stime = time_delta - ltime;
      am.posx = wmax;
      am.posy = ltime*am.speedy;
      auto [sx,sy] = getrs(false);
      am.speedx = -sx;
      am.speedy = sy;
      am.posx = stime*am.speedx;
      am.posy = stime*am.speedy;
    }
    if (mposx < 0) {
      double E1 = massarg*(am.speedx*am.speedx + am.speedy*am.speedy)/2;
      double ltime = (wmax-am.posx)/am.speedx;
      double stime = time_delta - ltime;
      am.posx = 0;
      am.posy = ltime*am.speedy;
      auto [sx,sy] = getrs(true);
      am.speedx = sx;
      am.speedy = sy;
      am.posx = stime*am.speedx;
      am.posy = stime*am.speedy;
      double E2 = massarg*(am.speedx*am.speedx + am.speedy*am.speedy)/2;
      energy +=abs(E1-E2);
      cout << "E1: " << E1 << " E2: " << E2 << " Delta: " << abs(E1-E2) << endl;
    }
    if (mposy > hmax) {
      am.posx = mposx;
      double ltime = (hmax-am.posy)/am.speedy;
      double stime = time_delta - ltime;
      am.posy = hmax;
      am.speedy = -am.speedy;
      am.posy = stime*am.speedy;
    }
    if (mposy < 0) {
      am.posx = mposx;
      double ltime = (hmax-am.posy)/am.speedy;
      double stime = time_delta - ltime;
      am.posy = 0;
      am.speedy = -am.speedy;
      am.posy = stime*am.speedy;
    }
    time+=time_delta;
    ++count;
  }

  cout << count << ' ' << time << endl;
  cout << energy << endl;

  return 0;
}