Working in a subspace with dimensionality much smaller than the dimension of the full Hilbert space, we deduce exact 4-particle ground states in 2D samples containing hexagonal repeat units and described by Hubbard type of models. The procedure identifies first a small subspace ${cal{S}}$ in which the ground state $|Psi_grangle$ is placed, than deduces $|Psi_grangle$ by exact diagonalization in ${cal{S}}$. The small subspace is obtained by the repeated application of the Hamiltonian $hat H$ on a carefully chosen starting wave vector describing the most interacting particle configuration, and the wave vectors resulting from the application of $hat H$, till the obtained system of equations closes in itself. The procedure which can be applied in principle at fixed but arbitrary system size and number of particles, is interesting by its own since provides exact information for the numerical approximation techniques which use a similar strategy, but apply non-complete basis for ${cal{S}}$. The diagonalization inside ${cal{S}}$ provides an incomplete image about the low lying part of the excitation spectrum, but provides the exact $|Psi_grangle$. Once the exact ground state is obtained, its properties can be easily analyzed. The $|Psi_grangle$ is found always as a singlet state whose energy, interestingly, saturates in the $U to infty$ limit. The unapproximated results show that the emergence probabilities of different particle configurations in the ground state present Zittern (trembling) characteristics which are absent in 2D square Hubbard systems. Consequently, the manifestation of the local Coulomb repulsion in 2D square and honeycomb types of systems presents differences, which can be a real source in the differences in the many-body behavior.