We investigated the swelling behavior and density distributions of polyion segments and multivalent counterions in polyelectrolyte brushes within a primitive model that represents both polyion segments and salt ions as charged spheres and the solvent as a continuous dielectric medium. The numerical performance of a previously proposed nonlocal density functional theory (NLDFT) was tested with Monte Carlo simulations for multivalent systems and compared with the Poisson-Boltzmann equation that ignores the ionic size and electrostatic correlations. The NLDFT predicts that in the limit of low salt concentration, trivalent counterions can self-organize into a one-dimensional wavelike structure within the oppositely charged polyelectrolyte brush. Multivalent counterions introduce electrostatic correlations responsible for the nonmonotonic swelling of the brush which is qualitatively different from that in a monovalent solution. While in a monovalent solution the brush thickness increases monotonically with the grafting density, in a trivalent solution the swelling behavior resembles that for a neutral brush in a poor solvent. The NLDFT was also used to investigate the effects of ionic size on the brush swelling and counterion organization.