The nonlinear dynamic equations of a multibody system composed of ?exible beams are derived by using the Lagrange multiplier method. The nonlinear Euler beam theory with inclusion of axial deformation e?ect is employed and its deformation ?eld is described by exact vibration modes. A numerical procedure for solving the dynamic equations is presented based on the Newmark direct integration method combined with Newton-Raphson iterative method. The results of numerical examples prove the correctness and e?ciency of the method proposed.