Numerical solution of $Q^2$ evolution equations for fragmentation functions

Semi-inclusive hadron-production processes are becoming important in high-energy hadron reactions. They are used for investigating properties of quark-hadron matters in heavy-ion collisions, for finding the origin of nucleon spin in polarized lepton-nucleon and nucleon-nucleon reactions, and possibly for finding exotic hadrons. In describing the hadron-production cross sections in high-energy reactions, fragmentation functions are essential quantities. A fragmentation function indicates the probability of producing a hadron from a parton in the leading order of the running coupling constant $\alpha_s$. Its $Q^2$ dependence is described by the standard DGLAP (Dokshitzer-Gribov-Lipatov-Altarelli-Parisi) evolution equations, which are often used in theoretical and experimental analyses of the fragmentation functions and in calculating semi-inclusive cross sections. The DGLAP equations are complicated integro-differential equations, which cannot be solved in an analytical method. In this work, a simple method is employed for solving the evolution equations by using Gauss-Legendre quadrature for evaluating integrals, and a useful code is provided for calculating the $Q^2$ evolution of the fragmentation functions in the leading order (LO) and next-to-leading order (NLO) of $\alpha_s$. The renormalization scheme is $\overline{MS}$ in the NLO evolution. Our evolution code is explained for using it in one's studies on the fragmentation functions.