Summary: Bi-quadratic programming over unit spheres is a fundamental problem in quantum mechanics introduced by pioneer work of Einstein, Schrödinger, and others. It has been shown to be NP-hard; so it must be solve by efficient heuristic algorithms such as the block improvement method (BIM). This paper focuses on the maximization of bi-quadratic forms with nonnegative coefficient tensors, which leads to a rank-one approximation problem that is equivalent to computing the M-spectral radius and its corresponding eigenvectors. Specifically, we propose a tight upper bound of the M-spectral radius for nonnegative fourth-order partially symmetric (PS) tensors. This bound, serving as an improved shift parameter, significantly enhances the convergence speed of BIM while maintaining computational complexity aligned with the initial shift parameter of BIM. Moreover, we elucidate that the computation cost of such upper bound can be further simplified for certain sets and delve into the nature of these sets. Building on the insights gained from the proposed bounds, we derive the exact solutions of the M-spectral radius and its corresponding M-eigenvectors for certain classes of fourth-order PS-tensors and discuss the nature of this specific category. Lastly, as a practical application, we introduce a testable sufficient condition for the strong ellipticity in the field of solid mechanics. Numerical experiments demonstrate the utility of the proposed results.