We have studied nuclear fusion reactions from the Coulomb explosion of deuterium clusters induced by high-intensity laser beams within a transport approach. By incorporating the D+D → n + He 3 channel as inelastic collisions based on the stochastic method, we have calibrated the neutron yield from the simulation in a box system with that from the reaction rate equation, and then have investigated the dynamics from Coulomb explosions of systems of different sizes and with different deuteron numbers in clusters. We have found that the kinetic energy spectrum of deuterons at final stages, which depends on the system size and deuteron number in clusters, is different from that when neutrons are abundantly produced, corresponding to significantly different reaction rates. We have also extrapolated the neutron yield result from small systems to large systems based on an intuitive parameterized form. The present framework can be extended by incorporating more channels, and useful for further studies of nuclear fusion reactions in plasma systems at higher energies reached in more recent experiments.