We study the existence and stability of localized states in the two-dimensional (2D) nonlinear Schrödinger (NLS)/Gross-Pitaevskii equation with a symmetric four-well potential. Using a fourmode approximation, we are able to trace the parametric evolution of the trapped stationary modes, starting from the corresponding linear limits, and thus derive the complete bifurcation diagram for the families of these stationary modes. The predictions based on the four-mode decomposition are found to be in good agreement with the numerical results obtained from the NLS equation. Actually, the stability properties coincide with those suggested by the corresponding discrete model in the large-amplitude limit. The dynamics of the unstable modes is explored by means of direct simulations. Finally, while we present the full analysis for the case of the focusing nonlinearity, the bifurcation diagram for the defocusing case is briefly considered too.