void LeapfrogStep(int part)
{
	int i;
	if (part == 1) {
		for (i = 0; i < nMol; i++) {
			mol[i].rv.x += 0.5 * deltaT * mol[i].ra.x;
			mol[i].rv.y += 0.5 * deltaT * mol[i].ra.y;
			mol[i].r.x += deltaT * mol[i].rv.x;
			mol[i].r.y += deltaT * mol[i].rv.y;
		}
	}
	else {
		for (i = 0; i < nMol; i++) {
			mol[i].rv.x += 0.5 * deltaT * mol[i].ra.x;
			mol[i].rv.y += 0.5 * deltaT * mol[i].ra.y;
		}
	}
}

void ComputeForces()
{
    int i, j;
    VecR dist, forces;
    real r, force;
    for (i = 0; i < nMol - 1; i++)
        for (j = i + 1; j < nMol; j++) {
            dist.x = mol[j].r.x - mol[i].r.x;
            dist.y = mol[j].r.y - mol[i].r.y;
            r = sqrt(dist.x * dist.x + dist.y * dist.y); //distance formula
            force = 24 * (r*r*r*r*r*r - 2) / (r*r*r*r*r*r*r*r*r*r*r*r*r); //derivative of the Lennard-Jones potential
            forces.x = dist.x / r * force;
            forces.y = dist.y / r * force;
            mol[i].rv.x += forces.x; //the program works correctly when I apply the force to rv (velocity) but not ra (acceleration)
            mol[i].rv.y += forces.y;
            mol[j].rv.x -= forces.x;
            mol[j].rv.y -= forces.y;
        }
}