It is well-known that the Levenberg–Marquardt method is a good choice for solving nonlinear equations, especially in the cases of singular/nonisolated solutions. We first exhibit some numerical experiments with local convergence, showing that this method for “generic” equations actually also works very well when applied to the specific case of the Lagrange optimality system, i.e., to the equation given by the first-order optimality conditions for equality-constrained optimization. In particular, it appears to outperform not only the basic Newton method applied to such systems, but also its modifications supplied with dual stabilization mechanisms, intended specially for tackling problems with nonunique Lagrange multipliers. The usual globalizations of the Levenberg–Marquardt method are based on linesearch for the squared Euclidean residual of the equation being solved. In the case of the Lagrange optimality system, this residual does not involve the objective function of the underlying optimization problem (only its derivative), and in particular, the resulting globalization scheme has no preference for converging to minima versus maxima, or to any other stationary point. We thus develop a special globalization of the Levenberg–Marquardt method when it is applied to the Lagrange optimality system, based on linesearch for a smooth exact penalty function of the optimization problem, which in particular involves the objective function of the problem. The algorithm is shown to have appropriate global convergence properties, preserving also fast local convergence rate under weak assumptions. © 2018, Springer Science+Business Media, LLC, part of Springer Nature.