Now every row along the X axis has the exact same evolution and there are no cross fluxes. The issue was that originally for the term on the advected internal energy energy the value of was computed on th same kernel were the conserved variables were updated, this means that there was a synchronization error because it was using values of the neighbors velocities to compute and at the same time those velocities were being updated.
The solution was to add the term to the internal energy in a separate kernel before the Update_Conserved kernel.
Now every row has the same values ( all overlap in the same line ) but the density still drops in the boundaries of the shock.
After fixing the bug, this are the results for the cosmological simulations using different values for :