We present an efficient and flexible method for solving the non-linear lasing equations of the steady-state ab initio laser theory. Our strategy is to solve the underlying system of partial differential equations directly, without the need of setting up a parametrized basis of constant flux states. We validate this approach in one-dimensional as well as in cylindrical systems, and demonstrate its scalability to full-vector three-dimensional calculations in photonic-crystal slabs. Our method paves the way for efficient and accurate simulations of microlasers which were previously inaccessible.