In this paper we are interested in the following nonlinear Choquard equation −∆u + (λV (x) − β)u = |x| −µ * |u| 2 * µ |u| 2 * µ −2 u in R N , where λ, β ∈ R + , 0 < µ < N , N ≥ 4, 2 * µ = (2N −µ)/(N −2) is the upper critical exponent due to the Hardy-Littlewood-Sobolev inequality and the nonnegative potential function V ∈ C(R N , R) such that Ω := intV −1 (0) is a nonempty bounded set with smooth boundary. If β > 0 is a constant such that the operator −∆ + λV (x) − β is non-degenerate, we prove the existence of ground state solutions which localize near the potential well int V −1 (0) for λ large enough and also characterize the asymptotic behavior of the solutions as the parameter λ goes to infinity. Furthermore, for any 0 < β < β 1 , we are able to prove the existence of multiple solutions by the Lusternik-Schnirelmann category theory, where β 1 is the first eigenvalue of −∆ on Ω with Dirichlet boundary condition.