  This paper presents a phenomenological study for the differential cross
section of the forward Z boson production in leptonic decay channels as a
function of the angular variable $\phi^{*}_{\eta}$ in proton-proton collisions.
The $\phi^{*}_{\eta}$ distribution is predicted for the forward pseudorapidity
region $2.0 < \eta_{l} < 4.5$ of the decay leptons at center-of-mass energies
8, 13, and 14 TeV. Accurate prediction of the $\phi^{*}_{\eta}$ distribution is
achieved by means of the state-of-the-art calculations including fixed-order
perturbative QCD and large logarithmic corrections. The predicted distributions
are obtained by employing the resummation either at next-to-next-to-leading
logarithmic (NNLL) or next-to-NNLL (N$^{3}$LL) accuracy which is matched to the
perturbative QCD calculation at next-to-next-to-leading order (NNLO) accuracy,
that is at NNLO+NNLL and NNLO+N$^{3}$LL, respectively. The Z boson
$\phi^{*}_{\eta}$ variable is experimentally preferable as it has been measured
more precisely than the Z boson transverse momentum though both variables probe
the same physics, thereby an accurate description of the $\phi^{*}_{\eta}$
variable is required by using theoretical predictions at both NNLO+NNLL and
NNLO+N$^{3}$LL accuracies. The matched predictions are compared with the
available 8 TeV and 13 TeV measurements by the LHCb experiment at the LHC and
found to be in good agreement. The 14 TeV predicted distributions at both
NNLO+NNLL and NNLO+N$^{3}$LL are also reported. In all the predicted results,
the NNLO+N$^{3}$LL provides more improved accuracy for the reliable description
of the $\phi^{*}_{\eta}$ distribution throughout its entire phase space region.
