Summary: A class of optimal control problems in which the control signals are bounded and appear linearly in the Hamiltonian function occur in many applications. In the solutions to these problems, the fact that optimal control trajectories may contain discontinuities and possibly singular arcs/sub-arcs makes the problem extremely difficult to solve, if not impossible. This paper presents a novel and efficient computational method for solving such problems. The method consists of two steps. In the first step, the original optimal control problem with possible discontinuities and singular arcs in control is converted into one with continuous and non-singular control trajectories by adding to the performance indexes a perturbed (or weighted) energy term. The resultant boundary value problem can easily be solved for an appropriately large value of the perturbation parameter. In the second step, a continuation method (imbedding method, or homotopy method) is developed to obtain the solution to the original problem by solving a set of initial value problems sequentially and/or in parallel as the perturbation parameter goes to zero. The proposed two-step method is computationally efficient since the resulting two-point boundary value problem is solved only once for a large perturbation parameter and the remaining problem becomes that of solving a set of initial value subproblems. The proposed algorithm is applicable to a large class of optimal control problems with various boundary conditions (e.g. fixed and free terminal time). The practicability of the method is demonstrated by computer simulations on several simple examples, including one with a singular arc along the optimal control trajectory.