We propose an optimization technique designed specifically for molecular structure optimization performed on an ab initio level. This gradient-based technique is a modification of quadratically convergent quasi-Newton method, and although it requires more energy evaluations than the conventional method, each of these energy evaluations is much cheaper due to O(N-3) scaling of the two-electron integrals evaluation. Statistics obtained from numerous optimization runs with Lennard-Jones molecules shows that the number of energy and gradient evaluations for the proposed technique is only 1.5-5 times (for 3-27 atoms, respectively) larger than that for conventional method. Given the great advantage of O(N-3) scaling of the two-electron integrals in the former, a substantial speedup of the overall computation can be achieved in certain cases. We consider the factors which affect the performance of the proposed technique and we also present timings and other details of several molecular structure optimization tests of the method on the ab initio level. Additionally, a novel approach to numerical Hessian evaluation during optimization is proposed, where the quality of the Hessian so obtained can be assessed.