The translation of experimental limits on the neutron electric dipole moment into constraints on heavy CP-violating physics beyond the Standard Model requires knowledge about non-perturbative matrix elements of effective operators, which ideally should be computed in lattice QCD. However, this necessitates a matching calculation as an interface to the effective field theory framework, which is based on dimensional regularization and renormalization by minimal subtraction. We calculate the one-loop matching between the gradient-flow and minimal-subtraction schemes for the CP-violating four-quark operators contributing to the neutron electric dipole moment. The gradient flow is a modern regularization-independent scheme amenable to lattice computations that promises, e.g., better control over power divergences than traditional momentum-subtraction schemes. Our results extend previous work on dimension-five operators and provide a necessary ingredient for future lattice-QCD computations of the contribution of four-quark operators to the neutron electric dipole moment.