This article proposes an adaptively bounded gradient descent (ABGD) algorithm for group elastic net penalized regression. Unlike previously proposed algorithms, the proposed algorithm adaptively bounds the Fisher information matrix, which results in a flexible and stable computational framework. In particular, the proposed algorithm (i) does not require orthogonalization of the predictors, and (ii) can be easily applied to any combination of exponential family response distribution and link function. The proposed algorithm is implemented in the grpnet R package (available from CRAN), which implements the approach for common response distributions (Gaussian, binomial, and Poisson), as well as several response distributions not previously considered in the group penalization literature (i.e., multinomial, negative binomial, gamma, and inverse Gaussian). Simulated and real data examples demonstrate that the proposed algorithm is as or more efficient than existing methods for Gaussian, binomial, and Poisson distributions. Furthermore, using two genomic examples, I demonstrate how the proposed algorithm can be applied to high-dimensional multinomial regression problems with grouped predictors. R code to reproduce the results is included as supplementary materials.