Within the framework of the density gradient theory, the interfacial tension can be calculated by finding the density profiles that minimize an integral of two terms over the system of infinite width. It is found that the two integrands exhibit a constant difference along the interface for a finite planar interface, and this difference becomes smaller as the interface width increases. These findings naturally lead to a solution procedure that consists of an inner loop and an outer loop for calculating the interfacial tension of a planar interface. The outer loop deals with the relationship between the interfacial tension and the interface width, and it permits us to obtain accurate results from finite width calculations. The inner loop minimizes the interfacial tension for a given interface width by adjusting the density profiles, in which the integrals are evaluated by a combination of Gauss-Lobatto quadrature and Lagrange interpolation based polynomial approximation. A better approximation of the interfacial tension is derived by a path integration along the density profiles. These strategies enable us to obtain accurate solutions with looser tolerance criteria and a fewer number of thermodynamic property evaluations compared to other methods. The performance of the algorithm with recommended parameters is analyzed for various systems, and the efficiency is further compared with the geometric-mean density gradient theory, which only needs to solve nonlinear algebraic equations. The results show that the algorithm is only 5-10 times less efficient than solving the geometric-mean density gradient theory. (C) 2017 Elsevier B.V. All rights reserved.