  Large scale molecular dynamics simulations for bidisperse nanoparticle
suspensions with an explicit solvent are used to investigate the effects of
evaporation rates and volume fractions on the nanoparticle distribution during
drying. Our results show that ``small-on-top'' stratification can occur when
${\rm Pe}_s\phi_s \gtrsim c$ with $c\sim 1$, where ${\rm Pe}_s$ is the
P\'{e}clet number and $\phi_s$ is the volume fraction of the smaller particles.
This threshold of ${\rm Pe}_s\phi_s$ for ``small-on-top'' is larger by a factor
of $\sim\alpha^2$ than the prediction of the model treating solvent as an
implicit viscous background, where $\alpha$ is the size ratio between the large
and small particles. Our simulations further show that when the evaporation
rate of the solvent is reduced, the ``small-on-top'' stratification can be
enhanced, which is not predicted by existing theories. This unexpected behavior
is explained with thermophoresis associated with a positive gradient of solvent
density caused by evaporative cooling at the liquid-vapor interface. For
ultrafast evaporation the gradient is large and drives the nanoparticles toward
the liquid/vapor interface. This phoretic effect is stronger for larger
nanoparticles and consequently the ``small-on-top'' stratification becomes more
distinct when the evaporation rate is slower (but not too slow such that a
uniform distribution of nanoparticles in the drying film is produced), as
thermophoresis that favors larger particles on the top is mitigated. A similar
effect can lead to ``large-on-top'' stratification for ${\rm Pe}_s\phi_s$ above
the threshold when ${\rm Pe}_s$ is large but $\phi_s$ is small. Our results
reveal the importance of including the solvent explicitly when modeling
evaporation-induced particle separation and organization and point to the
important role of density gradients brought about by ultrafast evaporation.
