Numerical implementation of evolution equations for twist-3 collinear PDFs

Twist-3 collinear parton distribution functions (PDFs) are matrix elements of quark-gluon-quark or three-gluons light-cone operators. They depend on three momentum fraction variables, which are restricted to a hexagon region, and the evolution kernels are defined via two-dimensional convolution in these variables. We present the numerical realisation of the twist-3 evolution equations at leading order in the strong coupling for all kinds of twist-3 PDF (quark, gluon, chiral-even/odd, etc). We provide two independent codes (in C and Fortran) that have been extensively cross-checked, and are ready-to-use. We supplement the paper with a review of known properties of twist-3 PDFs.