  We present a rapid analytic framework for predicting kilonova light curves
following neutron star (NS) mergers, where the main input parameters are
binary-based properties measurable by gravitational wave detectors (chirp mass
and mass ratio, orbital inclination) and properties dependent on the nuclear
equation of state (tidal deformability, maximum NS mass). This enables
synthesis of a kilonova sample for any NS source population, or determination
of the observing depth needed to detect a live kilonova given gravitational
wave source parameters in low latency. We validate this code, implemented in
the public MOSFiT package, by fitting it to GW170817. A Bayes factor analysis
overwhelmingly ($B>10^{10}$) favours the inclusion of an additional luminosity
source in addition to lanthanide-poor dynamical ejecta during the first day.
This is well fit by a shock-heated cocoon model, though differences in the
ejecta structure, opacity or nuclear heating rate cannot be ruled out as
alternatives. The emission thereafter is dominated by a lanthanide-rich viscous
wind. We find the mass ratio of the binary is $q=0.92\pm0.07$ (90% credible
interval). We place tight constraints on the maximum stable NS mass, $M_{\rm
TOV}=2.17^{+0.08}_{-0.11}$ M$_\odot$. For a uniform prior in tidal
deformability, the radius of a 1.4 M$_\odot$ NS is $R_{1.4}\sim 10.7$ km.
Re-weighting with a prior based on equations of state that support our credible
range in $M_{\rm TOV}$, we derive a final measurement
$R_{1.4}=11.06^{+1.01}_{-0.98}$ km. Applying our code to the second
gravitationally-detected neutron star merger, GW190425, we estimate that an
associated kilonova would have been fainter (by $\sim0.7$ mag at one day
post-merger) and declined faster than GW170817, underlining the importance of
tuning follow-up strategies individually for each GW-detected NS merger.
