We use density functional theory (B3LYP) and the COSMO continuum solvent model to characterize the structure and stability of the hydrated Cu(II) complexes [Cu(MeNH(2))(H(2)O)(n-1)](2+) and [Cu(OH)(x)(H(2)O)(n-x)](2-x) (x = 1-3) as a function of metal coordination number (4-6) and cluster size (n = 4-8, 18). The small clusters with n = 4-8 are found to be the most stable in the nearly square-planar four-coordinate configuration, except for [Cu(OH)(3)(H(2)O)](-), which is three-coordinate. In the presence of the two full hydration shells (n = 18), however, the five-coordinate square-pyramidal geometry is the most favorable for Cu(MeNH(2))(2+) (5, 6) and Cu(OH)(+) (5, 4, 6), and the four-coordinate geometry is the most stable for Cu(OH)(2) (4, 5) and Cu(OH)(3)(-) (4). (Other possible coordination numbers for these complexes in the aqueous phase are shown in parentheses.) A small energetic difference between these structures (0.23-2.65 kcal/mol) suggests that complexes with different coordination numbers may coexist in solution. Using two full hydration shells around the Cu(2+) ion (18 ligands) gives Gibbs free energies of aqueous reactions that are in excellent agreement with experiment. The mean unsigned error is 0.7 kcal/mol for the three consecutive hydrolysis steps of Cu(2+) and the complexation of Cu(2+) with methylamine. Conversely, calculations for the complexes with only one coordination shell (four equatorial ligands) lead to a mean unsigned error that is >6.0 kcal/mol. Thus, the explicit treatment of the first and the second shells is critical for the accurate prediction of structural and thermodynamic properties of Cu(II) species in aqueous solution.