We present a novel and flexible method to optimize the phase response of reflective metasurfaces towards a desired scattering profile. The scattering power is expressed as a spin-chain Hamiltonian using the radar cross section formalism. For metasurfaces reflecting an oblique plane wave, an Ising Hamiltonian is obtained. Thereby, the problem of achieving the scattering profile is recast into finding the ground-state solution of the associated Ising Hamiltonian. To rapidly explore the configuration states, we encode the Ising coefficients with quantum annealing algorithms, taking advantage of the fact that the adiabatic evolution efficiently performs energy minimization in the Ising model. Finally, the optimization problem is solved on the D-Wave 2048-qubit quantum adiabatic optimizer machine for binary phase as well as quadriphase reflective metasurfaces. Even though the work is focused on the phase modulation of metasurfaces, we believe this approach paves the way to fast optimization of reconfigurable intelligent surfaces that are mod- ulated in both amplitude and phase for multi-beam generation in and beyond 5G/6G mobile networks.