A numerical algorithm is proposed for multi-orbital slave-boson mean field approach through the integrating pattern search method, the generalized Lagrange multiplier method, and the Rosenbrock method. Since the crystal field splitting, inter-orbital hopping and realistic band structures can be considered, the proposed slave-boson mean field approach can be utilized to study realistic material. To validate our algorithm, the Mott transitions in twoorbital Hubbard models are studied with the elliptical density of states. The results are consistent with the reported ones available. Then we use this method to study the correlation effect on the three-orbital Hubbard model for NaxCoO2. It is shown that the six small Fermi surfaces constructed by the eg' orbital vanish in the intermediate Coulomb correlations. The physical reason is that the hole occupations of eg' orbital decrease with U increasing. All the calculated results verify the accuracy and the efficiency of our numerical algorithm.