Using a data sample of proton-proton collisions at $\sqrt{s} =$ 13 TeV, corresponding to an integrated luminosity of 140 fb$^{-1}$ collected by the CMS experiment in 2016-2018, the B$^0_\mathrm{s}$ $\to$ X(3872)$\phi$ decay is observed. Decays into J/$\psi\,\pi^+\pi^-$ and K$^+$K$^-$ are used to reconstruct, respectively, the X(3872) and $\phi$. The ratio of the product of branching fractions $\mathcal{B}($B$^0_\mathrm{s}$ $\to$ X(3872)$\phi)\,\mathcal{B}($X(3872) $\to$ J$/\psi\,\pi^+\pi^-)$ to the product $\mathcal{B}($B$^0_\mathrm{s}$ $\to$ $\psi$(2S)$\phi)\,\mathcal{B}(\psi$(2S) $\to$ J/$\psi\,\pi^+\pi^-)$ is measured to be (2.21 $\pm$ 0.29 (stat) $\pm$ 0.17 (syst))%. The ratio $\mathcal{B}($B$^0_\mathrm{s}$ $\to$ X(3872)$\phi) / \mathcal{B}($B$^{0}$ $\to$ X(3872)K$^0)$ is found to be consistent with one, while the ratio $\mathcal{B}($B$^0_\mathrm{s}$ $\to$ X(3872)$\phi) / \mathcal{B}($B$^+$ $\to$ X(3872)K$^+)$ is two times smaller. This suggests a difference in the production dynamics of the X(3872) in B$^0$ and B$^0_\mathrm{s}$ meson decays compared to B$^+$. The reported observation may shed new light on the nature of the X(3872) particle.