  Spatio-temporal evolution of Buneman instability has been followed
numerically till its quasilinear quenching and beyond, using an in-house
developed electrostatic 1D particle-in-cell simulation code. For different
initial drift velocities $k_{L}v_{0}/\omega_{pe} \approx 0.1 \, - \, 1$ and for
a wide range of electron to ion mass ratios (m/M), growth rate obtained from
simulation agrees well with the numerical solution of the fourth order
dispersion relation. Quasi-linear saturation of Buneman instability occurs when
ratio of electrostatic field energy density ($\sum\limits_{k}
|E_{k}|^{2}/8{\pi}$) to initial electron drift kinetic energy density ($W_{0} =
\frac{1}{2} n_{0}m v^{2}_{0}$) reaches up to a constant value, which as
predicted by Hirose [Plasma Physics 20, 481(1978)], is independent of initial
electron drift velocity but depends on electron to ion mass ratio m/M as
$\sum\limits_{k} |E_{k}|^{2}/16{\pi}W_{0} \approx (m/M)^{1/3}$. This result
stands verified in our simulations. Growth of the instability beyond the first
saturation (quasilinear saturation ) till its final saturation [Ishihara et.
al., PRL 44, 1404(1980)] follows an algebraic scaling with time. In contrast to
the quasilinear saturation, the ratio of final saturated electrostatic field
energy density to initial kinetic energy density, is relatively independent of
electron to ion mass ratio and is found to depend only on the initial drift
velocity. Beyond the final saturation, electron phase space holes coupled to
large amplitude ion solitary waves, a state known as coupled hole-soliton , are
seen in our simulations. The propagation characteristics ( amplitude - speed
relation ) of these coherent modes is found to be consistent with the theory of
Saeki et. al. [PRL 80, 1224(1998)].
