We study the near-threshold incoherent phi photoproduction on the deuteron based on a model of gamma N -> phi N, consisting of Pomeron, (pi,eta) exchanges, and a J(P) = 3/2(-) resonance, which describes the low-energy gamma p -> phi p LEPS data well, including the peak in the forward differential cross section. The calculation is done up to double rescatterings, with the spin dependence of the elementary gamma N -> phi N amplitude retained throughout the calculation. The Fermi motion and final-state interactions (FSIs) are all properly treated as prescribed by realistic nucleon-nucleon interaction. The couplings of the resonance to gamma n and phi n channels are estimated with the help of a constituent quark model. The main features of the LEPS and CLAS data are described reasonably well except for some quantitative discrepancies at very low energies and low-momentum-transfer regions. It is found that contributions of Fermi motion, pn FSI, and resonance are all indispensable in bridging the differences between the single-scattering results and the data. The off-shell rescattering is found to be important because it cancels out a large portion of the on-shell contribution. The discrepancies at low-momentum-transfer regions might be related to the binning size of the data. No peak is found to be associated with the weak resonance because it gets smeared out by the Fermi motion and FSI with the deuterium target. The problem at very-low-energy regions hints at the possible contributions from other mechanisms and should be investigated in depth with the use of recent high-statistics gamma p -> phi p data from CLAS.