For an $n\times n$ matrix $A_n$, the $r\to p$ operator norm is defined as $$\|A_n\|_{r\to p}:= \sup_{\mathbf{x}\in\mathbb{R}^n:\|\mathbf{x} \|_r\leq 1 } \|A_n\mathbf{x} \|_p\quad\text{for}\quad r,p\geq 1.$$ For different choices of $r$ and $p$, this norm corresponds to key quantities that arise in diverse applications including matrix condition number estimation, clustering of data, and construction of oblivious routing schemes in transportation networks. This article considers $r\to p$ norms of symmetric random matrices with nonnegative entries, including adjacency matrices of Erdős-Rényi random graphs, matrices with positive sub-Gaussian entries, and certain sparse matrices. For $1<p\leq r<\infty$, the asymptotic normality, as $n\to\infty$, of the appropriately centered and scaled norm $\|A_n\|_{r\to p}$ is established. When $p \geq 2$, this is shown to imply asymptotic normality of the solution to the $\ell_p$ quadratic maximization problem, also known as the $\ell_p$ Grothendieck problem. Furthermore, a sharp $\ell_\infty$-approximation bound for the unique maximizing vector in the definition of $\|A_n\|_{r\to p}$ is obtained, and may be viewed as an $\ell_\infty$-stability result of the maximizer under random perturbations of the matrix with mean entries. This result is in fact shown to hold for a broad class of deterministic sequences of matrices having certain asymptotic expansion properties. The results obtained can be viewed as a generalization of the seminal results of Füredi and Komlós (1981) on asymptotic normality of the largest singular value of a class of symmetric random matrices. In the general case with $1<p\leq r< \infty$, spectral methods are no longer applicable, and so a new approach is developed involving a refined convergence analysis of a nonlinear power method and a perturbation bound on the maximizing vector, which may be of independent interest.