In this paper, we consider approximation algorithms for optimizing a generic multi-variate homogeneous polynomial function, subject to homogeneous quadratic constraints. Such optimization models have wide applications, e.g., in signal processing, magnetic resonance imaging (MRI), data training, approximation theory, and portfolio selection. Since polynomial functions are nonconvex in general, the problems under consideration are all NP-hard. In this paper we shall focus on polynomial-time approximation algorithms. In particular, we first study optimization of a multi-linear tensor function over the Cartesian product of spheres. We shall propose approximation algorithms for such problem and derive worst-case performance ratios, which are shown to depend only on the dimensions of the models. The methods are then extended to optimize a generic multi-variate homogeneous polynomial function with spherical constraints. Likewise, approximation algorithms are proposed with provable relative ...