  A search for rare $B \rightarrow D \mu^+ \mu^-$ decays is performed using
proton-proton collision data collected by the LHCb experiment, corresponding to
an integrated luminosity of 9 fb$^{-1}$. No significant signals are observed in
the non-resonant $\mu^+\mu^-$ modes, and upper limits of $\mathcal{B}(B^0
\rightarrow \overline{D}^0 \mu^+ \mu^-) < 5.1 \times 10^{-8}$, $\mathcal{B}(B^+
\rightarrow D_s^+ \mu^+ \mu^-) < 3.2 \times 10^{-8}$, $\mathcal{B}(B_s^0
\rightarrow \overline{D}^0 \mu^+ \mu^-) < 1.6 \times 10^{-7}$ and $f_c/f_u
\cdot \mathcal{B}(B_c^+ \rightarrow D_s^+ \mu^+ \mu^-) < 9.6 \times 10^{-8}$
are set at the 95\% confidence level, where $f_c$ and $f_u$ are the
fragmentation fractions of a $B$ meson with a $c$ and $u$ quark respectively in
proton-proton collisions. Each result is either the first such measurement or
an improvement by three orders of magnitude on an existing limit. Separate
upper limits are calculated when the muon pair originates from a $J/\psi
\rightarrow \mu^+ \mu^-$ decay. The branching fraction of $B_c^+ \rightarrow
D_s^+ J/\psi$ multiplied by the fragmentation-fraction ratio is measured to be
$f_c/f_u \cdot \mathcal{B}(B_c^+ \rightarrow D_s^+ J/\psi ) = (1.63 \pm 0.15
\pm 0.13) \times 10^{-5}$, where the first uncertainty is statistical and the
second systematic.
