We introduce the large average subtensor problem: given an order-$p$ tensor over $\mathbb{R}^{N\times \cdots \times N}$ with i.i.d. standard normal entries and a $k\in\mathbb{N}$, algorithmically find a $k\times \cdots \times k$ subtensor with a large average entry. This generalizes the large average submatrix problem, a key model closely related to biclustering and high-dimensional data analysis, to tensors. For the submatrix case, Bhamidi, Dey, and Nobel~\cite{bhamidi2017energy} explicitly highlight the regime $k=Θ(N)$ as an intriguing open question. Addressing the regime $k=Θ(N)$ for tensors, we establish that the largest average entry concentrates around an explicit value $E_{\mathrm{max}}$, provided that the tensor order $p$ is sufficiently large. Furthermore, we prove that for any $γ>0$ and large $p$, this model exhibits multi Overlap Gap Property ($m$-OGP) above the threshold $γE_{\mathrm{max}}$. The $m$-OGP serves as a rigorous barrier for a broad class of algorithms exhibiting input stability. These results hold for both $k=Θ(N)$ and $k=o(N)$. Moreover, for small $k$, specifically $k=o(\log^{1.5}N)$, we show that a certain polynomial-time algorithm identifies a subtensor with average entry $\frac{2\sqrt{p}}{p+1}E_{\mathrm{max}}$. In particular, the $m$-OGP is asymptotically sharp: onset of the $m$-OGP and the algorithmic threshold match as $p$ grows. Our results show that while the case $k=Θ(N)$ remains open for submatrices, it can be rigorously analyzed for tensors in the large $p$ regime. This is achieved by interpreting the model as a Boolean spin glass and drawing on insights from recent advances in the Ising $p$-spin glass model.