  We have implemented the 3+1 dimensional CLVisc hydrodynamics model with
TRENTO-3D initial conditions to investigate the spin polarization of $\Lambda$
hyperons along the beam direction in p+Pb collisions at $\sqrt{s_{NN}} = 8.16$
TeV. Following our previous theoretical framework based on quantum kinetic
theory, we consider three different scenarios: $\Lambda$ equilibrium, $s$ quark
equilibrium, and iso-thermal equilibrium scenarios. We have computed the second
Fourier sine coefficients of spin polarization along the beam direction,
denoted as $\left\langle P_{z} \sin 2(\phi_{p} - \Psi_{2}) \right\rangle$, with
$\phi_p$ being the azimuthal angle relative to the second-order event plane
$\Psi_{2}$, as functions of multiplicity and transverse momentum in the three
scenarios. Additionally, we have also computed the spin polarization along the
beam direction, $P_{z}$, as a function of the azimuthal angle. We find that the
spin polarization induced by thermal vorticity always provides an opposite
contribution compared to the shear-induced polarization in p+Pb collisions. The
total spin polarization computed by the current hydrodynamic model disagrees
with the data measured by LHC-CMS experiments. Our findings imply that other
non-flow effects may play a crucial role in p+Pb collisions.
