This paper proposes a framework for fully automatic gradient-based constrained aerodynamic shape optimization in a multirow turbomachinery environment. The concept of adjoint-based gradient calculation is discussed and the development of the discrete adjoint equations for a turbomachinery Reynolds-averaged Navier–Stokes (RANS) solver, particularly the derivation of flow-consistent adjoint boundary conditions as well as the implementation of a discrete adjoint mixing-plane formulation, are described in detail. A parallelized, automatic grid perturbation scheme utilizing radial basis functions (RBFs), which is accurate and robust as well as able to handle highly resolved complex multiblock turbomachinery grid configurations, is developed and employed to calculate the gradient from the adjoint solution. The adjoint solver is validated by comparing its sensitivities with finite-difference gradients obtained from the flow solver. A sequential quadratic programming (SQP) algorithm is then utilized to determine an improved blade shape based on the gradient information from the objective functional and the constraints. The developed optimization method is used to redesign a single-stage transonic flow compressor in both inviscid and viscous flow. The design objective is to maximize the isentropic efficiency while constraining the mass flow rate and the total pressure ratio.