<p>Variable-section leaf springs are widely used in modern vehicle suspensions for their lightweight design, low friction, and noise reduction. Their dynamic behavior directly affects ride comfort, braking, and handling stability. This paper develops a rigid-flexible coupled dynamic modeling approach for suspension systems with variable-section leaf springs using the natural coordinate formulation (NCF) and the absolute nodal coordinate formulation (ANCF). A higher-order ANCF shell element with segmented thickness functions is adopted to capture thickness variation and accurately represent distributed inertia and elasticity. The prestressed configuration of the assembled leaf springs is determined through the quasi-static simulation, and an improved normal contact force model is derived to characterize the rubber spacer’s hyper-elastic behavior. By applying the NCF-ANCF framework to vehicle suspension systems, rigid and flexible components can be modeled within a unified generalized coordinate system. The resulting system dynamic equations feature a constant mass matrix and a concise, naturally expressed set of constraint equations, enabling seamless integration with existing multibody dynamics algorithms. The developed model is validated against experimental data, showing good agreement in predicting the static stiffness and stress distribution. Further dynamic simulations demonstrate that the proposed rigid-flexible coupled suspension model effectively captures the system’s complex dynamic responses and spatial axle attitudes under asymmetric excitation. This study provides an effective modeling and integration method for high-fidelity simulation of vehicle dynamics involving variable-section leaf springs.</p>