In this paper, an Optimal Homotopy Analysis Method (Optimal HAM) is applied to solve the linear optimal control problems (OCPs), which have a quadratic performance index. This approach contains at most two convergence-control parameters which depend on the control system and is computationally rather efficient. A squared residual error for the system is defined, which can be used to find the unknown optimal convergence-control parameters by using Mathematica package BVPh (version 2.0). The results of comparisons among the proposed method, the homotopy perturbation method (HPM), the Adomian decomposition method (ADM), the differential transform method (DTM) and the homotopy analysis method (HAM) provide verification for the validity of the proposed approach. Moreover, numerical results are presented by several examples involving scalar and 2nd-order systems to clarify the efficiency and high accuracy of the proposed approach.