Different from previous discrete element methods (DEM), where irregular 3D particle shapes are approximated by subspheres, vertices, or voxels, this study aims to develop an innovative and computationally effective DEM method directly employing spherical harmonic functions for simulations of 3D irregular‐shaped particles. First, the discrete surface points of a 3D irregular‐shaped particle are represented by spherical harmonic functions with only a limited number of harmonic coefficients to restore the particle morphology. Then, the intrinsic physical quantities are computed directly using spherical harmonic functions. Next, specific algorithms for interparticle overlapping detection and contact resolution involving the spherical harmonic functions are developed. Subsequently, the interparticle contact forces, moments, and particle movements are computed. The feasibility and capability of the proposed 3D method are verified by simulating random deposition of superellipsoids, repose angle tests, and triaxial tests on particles with various shapes. The proposed method could pave a viable pathway for realistic modeling of granular media pertaining to various engineering and industrial processes.