Abstract: The design of shelling equipment for fresh Camellia oleifera is constrained by lack of a discrete element model that characterizes the multi-layer mechanical response during shelling. Single-layer models treat fruit as homogeneous, failing to distinguish shell from seed; double-layer models separate shell and seed but treat seed as uniform, missing the sequential failure from shell to seed coat to kernel Given that the Camellia oleifera consists of a three-layer composite structure comprising the outer shell, the seed coat, and the seed kernel, a multi-layer flexible discrete element bonding model was established based on the Hertz-Mindlin with Bonding contact model. Within each layer, the Bonding mechanism was employed to characterize material fracture behavior through bond breakage between constituent particles, while the basic Hertz-Mindlin contact model was applied between adjacent layers to simulate contact compression and extrusion interaction. This heterogeneous modeling strategy enables a physically realistic representation of intra-layer fracture and inter-layer separation, a distinction that conventional double-layer models cannot achieve. Geometric parameters—including transverse diameter, longitudinal diameter, and shell thickness—and physical properties—including density, Poisson's ratio, and shear modulus—of fresh Camellia oleifera were determined through systematic physical experiments. Contact characteristic parameters, including coefficients of restitution, static friction, and rolling friction for various material pairings, were calibrated by combining experimental measurements with EDEM simulation results. The bonding parameters for the seed coat–seed coat, kernel–kernel, and shell–shell contacts were calibrated separately using a multi-stage experimental design protocol comprising Plackett-Burman screening, steepest ascent climbing, and response surface optimization. For the tea seed model, the relative errors between simulated and experimentally measured peak cracking forces in the transverse and longitudinal directions were 3.17% and 7.13%, respectively. For the whole fruit model, the corresponding relative errors were 1.54% and 1.83%, respectively, indicating that the calibrated model accurately predicts the critical failure loads. To comprehensively evaluate model fidelity beyond single-point peak force comparisons, the experimentally measured and simulated compression load-displacement curves were compared across the entire loading history. The coefficient of determination R² was employed to quantify the overall agreement between the two curves throughout the elastic deformation, crack initiation, and crack propagation stages. For the whole fruit, the R² values for transverse and longitudinal compression were 0.92 and 0.90, respectively; for the tea seeds, the R² values were 0.89 and 0.87, respectively. These high R² values confirm that the simulation not only matches the peak failure loads but also reproduces the full force-deformation response with excellent fidelity. Furthermore, macro-cracks in the shells consistently appeared along the longitudinal valve lines in both physical experiments and simulations, corroborating the model's capacity to capture the correct failure morphology. The developed model effectively characterizes the progressive cracking behavior of the fruit shell and the layer-by-layer mechanical response within the tea seeds—from seed coat rupture to kernel deformation—thereby providing a robust theoretical basis and model support for the simulation-driven design and process optimization of dehusking equipment for Camellia oleifera.