A new algorithm for constrained molecular dynamics is proposed. In contrast to the standard approach, the constrained bond-length/bond-angle value is adjusted at each time step so that total energy is minimized with respect to the constrained distances. This can be viewed as modifying the equilibrium bond-length/bond-angle according to external and centripetal forces. Two approaches are constructed to implement the algorithm. Method I includes all energy terms, but it is neither holonomic nor symplectic. Method II neglects a rotational kinetic energy term, resulting in a more expensive symplectic integrator. Both integrators are reversible and well conserve total energy. Due to the iterative nature of the optimization, which requires an energy evaluation at each iteration, these methods are intended for use in conjunction with other integration techniques that require an optimization at each time step iteration, such as the self-consistent electronic polarization methods or iterative combined quantum mechanics/molecular mechanics (QM/MM) hybrid simulations. Simulation results for collisions between two diatomic molecules, two water molecules, and for a periodic water box are compared and contrasted with fully constrained and free dynamics. These results demonstrate the utility and efficacy of these new methods.