Computing CP-violating nucleonic matrix elements on the lattice allows one to place theoretical constraints on the couplings of effective interactions related to BSM sources of CP-violation. These interactions are related to local operators that mix under renormalization. Typically, this mixing is parametrized by the only scale available, the lattice spacing, and induces local divergences in the coefficients of lower-dimensional operators, obscuring the continuum limit. The gradient flow has become an attractive method to circumvent this problem. In adopting the flow to define renormalized operators, the renormalization and mixing scales are disentangled, allowing for a clean computation of the corresponding matching (Wilson) coefficients. Perturbative calculations within the gradient flow formalism can be used to fix the high-energy behavior of the matching coefficients, so that the matrix elements are renormalized across a wide range of energy scales. We present results on the renormalization and mixing of the gluon chromoelectric dipole moment (gCEDM) operator to one-loop order in perturbation theory. These include the power-divergent mixing of the gCEDM with the topological charge density and the logarithmic mixing with various dimension-six operators. We also discuss the construction of a basis compatible with the chiral anomaly.