<p>Clearance joint systems exhibit various nonlinear factors, such as freeplay, variable stiffness, bilinearity, which give rise to multiple motion states. Transitions between these motion states are typically non-differentiable, resulting in abrupt changes in the system’s dynamic characteristics. Traditional parameter continuation strategies in the shooting method become invalid at these state transition points. This paper presents an improved shooting method for nonlinear normal mode analysis of clearance joint systems, incorporating an adaptive step size adjustment strategy. First, a hyperbolic tangent kernel function is introduced to smooth the local non-differentiable points in the one-dimensional clearance joint system. Next, a shooting method combined with a pseudo-arclength continuation algorithm is applied to compute and track the relationship between system frequency characteristics and energy. At state transition points, the adaptive step size adjustment strategy facilitates the continuation of periodic solutions. Finally, the stability of the solutions is evaluated using Floquet multipliers. Numerical examples demonstrate that clearances reduce structural modal frequencies and induce multistage behavior in frequency–energy curves. Additionally, clearance with bilinear characteristics generates novel internal resonance tongues.</p>