  We present a comprehensive, configurable open-source framework for estimating
the rate of electromagnetic detection of kilonovae (KNe) associated with
gravitational wave detections of binary neutron star (BNS) mergers. We simulate
the current LIGO-Virgo-KAGRA (LVK) observing run (O4) using up-to-date
sensitivity and up-time values as well as the next observing run (O5) using
predicted sensitivities. We find the number of discoverable kilonovae during
LVK O4 to be ${ 1}_{- 1}^{+ 4}$ or ${ 2 }_{- 2 }^{+ 3 }$, (at 90% confidence)
depending on the distribution of NS masses in coalescing binaries, with the
number increasing by an order of magnitude during O5 to ${ 19 }_{- 11 }^{+ 24
}$. Regardless of mass model, we predict at most five detectable KNe (at 95%
confidence) in O4. We also produce optical and near-infrared light curves that
correspond to the physical properties of each merging system. We have collated
important information for allocating observing resources and directing search
and follow-up observations including distributions of peak magnitudes in
several broad bands and timescales for which specific facilities can detect
each KN. The framework is easily adaptable, and new simulations can quickly be
produced as input information such as merger rates and NS mass distributions
are refined. Finally, we compare our suite of simulations to the thus-far
completed portion of O4 (as of October 14, 2023), finding a median number of
discoverable KNe of 0 and a 95-percentile upper limit of 2, consistent with no
detection so far in O4.
