  Accurate determination of the full momentum-dependent spin susceptibility
$\chi(\mathbf{q}) $ is very important for the description of magnetism and
superconductivity. While in principle the formalism for calculating
$\chi(\mathbf{q})$ in the linear response density functional theory (DFT) is
well established, hardly any publicly available code includes this capability.
Here, we describe an alternative way to calculate the static
$\chi(\mathbf{q})$, which can be applied to most common DFT codes without
additional programming. The method combined standard fixed-spin-moment
calculations of $\chi(\mathbf{0}) $ with direct calculations of the energy of
spin spirals stabilized by an artificial Hubbard interaction. From these
calculations, $\chi_{DFT}(\mathbf{q} )$ can be extracted by inverting the RPA
formula. We apply this recipe to the recently discovered Ising
superconductivity in NbSe$_2$ monolayer, one of the most exciting findings in
superconductivity in recent years. It was proposed that spin fluctuations may
strongly affect the parity of the order parameter. Previous estimates suggested
proximity to ferromagnetism, $i.e.$, $\chi(\mathbf{q})$ peaked at
$\mathbf{q}=0$. We find that the structure of spin fluctuations is more
complicated, with the fluctuation spectrum sharply peaked at $\mathbf{q}\approx
(0.2,0)$. Such a spectrum would change the interband pairing interaction and
considerably affect the superconducting state.
