In this study an advanced approach for describing diffusion fluxes in a stochastically nonhomogeneous two-phase body with randomly located spherical inclusions is presented. The research is focused on modeling the mass transfer of an impurity substance, considering the dominant phase within the body. The model accounts for arbitrary probability distributions of inclusions and applies Fick’s first law to determine diffusion flux. Based on the representation of the impurity concentration in the form of the Neumann integral series, calculation formulas for the averaged flux and the amount of impurity substance that passed through the given cross-section of the two-phase body with uniformly distributed spherical inclusions are obtained. Software modules are developed to perform numerical simulations, offering a computational tool for further analysis of diffusion process in randomly nonhomogeneous structures.