We report a measurement of $e^+e^-$ pairs from semileptonic heavy-flavor decays in $p$+$p$ collisions at $sqrt{s_{NN}}=200$~GeV. The $e^+e^-$ pair yield from $bbar{b}$ and $cbar{c}$ is separated by exploiting a double differential fit done simultaneously in dielectron invariant mass and $p_T$. We used three different event generators, {sc pythia}, {sc mc@nlo}, and {sc powheg}, to simulate the $e^+e^-$ spectra from $cbar{c}$ and $bbar{b}$ production. The data can be well described by all three generators within the detector acceptance. However, when using the generators to extrapolate to $4pi$, significant differences are observed for the total cross section. These difference are less pronounced for $bbar{b}$ than for $cbar{c}$. The same model dependence was observed in already published $d$+$A$ data. The $p$+$p$ data are also directly compared with $d$+$A$ data in mass and $p_T$, and within the statistical accuracy no nuclear modification is seen.