Dislocation slip-based crystal plasticity models have been a great success in connecting the fundamental physics with the macroscopic deformation of crystalline materials. Pioneered by Taylor in his work on “plastic strain in metals” (Taylor, 1938), and further advanced by Bishop and Hill (1951a, 1951b), the Taylor–Bishop–Hill theory laid the foundation of today’s constitutive models on crystal plasticity. An intriguing part of those modeling is to determine the active slip systems—which system to be involved in and how much it contributes to the deformation. In this paper, we developed a machine learning-based algorithm to determine accurately and efficiently the active slip systems in crystal plasticity constitutive models. Applications to the common three polycrystalline metals, face-centered cubic (FCC) copper, body-centered cubic (BCC) α-iron, and hexagonal close-packed (HCP) AZ31B, demonstrate that even a simple neural network could give rise to accurate and efficient results in comparing with traditional routines. There seems to be plenty of space for further reducing the computation time and hence scaling up the simulating samples.