The dynamic model of flexible multi-body system with closed loop is described with natural coordinates. The method is formulaic and especially appropriate to the modeling of system with repetitive substructures, but Lagrange multipliers as unknown variables are contained in the model. The dynamic model is reduced to a system with minimum dimensions by means of null-space orthogonal basis of constraint Jacobin. The result is correspondence with the dynamic model of rigid multi- body system while the rigid and flexible coupling is considered. The numerical method and the simulation steps are given in detail. Numerical example dealing with a flexible parallel four-bar linkage widely used in engineering illustrates the performance of the proposed method. The dynamic response of the angle between flexible linkage and horizontal axis are presented. The results indicate that there are clearly fluctuation in the angular velocity and acceleration of crank. The motion of connecting linkage contains not only translation but also narrow range rotation especially in the start instant.