  We study thermalization of transverse field Ising chain with power law
decaying interaction $\sim 1/r^{\alpha}$ following a global quantum quench of
the transverse field to two different dynamical regimes. We quantify the
thermalization behavior by comparing the full probability distribution function
(PDF) of the evolving states with the corresponding thermal state given by the
Gibbs canonical ensemble (GCE). To this end, we use matrix product state (MPS)
based time dependent variational principle (TDVP) algorithm to simulate both
real time evolution following a global quantum quench and the finite
temperature density operator. We observe that thermalization is strongly
suppressed in the region with strong confinement for all the interaction
strength $\alpha$ considered whereas thermalization occurs in the region with
weak confinement.
