In this paper we proposed a new symbolic, non-standard recursive and fast orthonormalization procedure of linearly independent vectors but as in other approaches not orthonormal based on the Gram-Schmidt orthonormalization algorithm. Our adaptation of the Gram-Schmidt orthonormalization procedure provide simple analytic formulas for the SU(3) Bargmann-Moshinsky basis orthonormalization coefficients and do not involve any square root operation on the expressions coming from the previous iterative computation steps. This distinct features of the proposed orthonormalization algorithm may make the large scale symbolic calculations feasible. We demonstrate efficiency of our procedure by benchmark large-scale calculations of the non-canonical BM basis with the highest weight vectors of SO(3) irreducible representations.