[en] A Monte Carlo method is presented for the calculation of the QCD evolution of structure functions. Its application is discussed in detail in the framework of the LLA, but it can also be used with modified parton decay probability functions including higher-order effects. For heavy quark production, threshold constraints can be correctly taken into account, and one obtains results which at low Q2 are consistent with those of the photon-gluon fusion model. (orig.)