Solution (source code)

= Solution

Define the physical <velocity potential> by $\psi=\sqrt n\exp(im\phi/\hbar)$, so the <superfluid velocity> is $\mathbf v=\nabla\phi$. The <condensate number current> is
$$
\mathbf j=\frac\hbar m\operatorname{Im}(\psi^*\nabla\psi)=n\nabla\phi.
$$
Multiplying the pumped <Gross–Pitaevskii equation> by $\psi^*$ and subtracting its <complex conjugate> gives the density equation. Equivalently, substituting the amplitude and phase and equating imaginary parts yields
$$
\boxed{n_t+\nabla\cdot(n\nabla\phi)=\frac2\hbar(\gamma-\Gamma n)n}.
$$
Equating real parts gives the other <gain-loss Madelung equations> relation,
$$
\boxed{m\phi_t+\frac m2|\nabla\phi|^2+V_0n
-\frac{\hbar^2}{2m}\frac{\nabla^2\sqrt n}{\sqrt n}=0}.
$$
The last term is the <quantum potential>, retained wherever $n>0$. Taking a <gradient> gives the equivalent velocity equation
$$
 m[\mathbf v_t+(\mathbf v\cdot\nabla)\mathbf v]
=-\nabla\left(V_0n-\frac{\hbar^2}{2m}\frac{\nabla^2\sqrt n}{\sqrt n}\right)
$$
in an irrotational region. In equilibrium $n_t=0$ and $\phi(\mathbf x,t)=\Phi(\mathbf x)-\mu t/m$, so
$$
\boxed{\nabla\cdot(n\nabla\Phi)=\frac{2n}{\hbar}(\gamma-\Gamma n)},\qquad
\mu=\frac m2|\nabla\Phi|^2+V_0n-
\frac{\hbar^2}{2m}\frac{\nabla^2\sqrt n}{\sqrt n}.
$$
This convention distinguishes the physical <velocity potential> from the dimensionless <quantum phase> $m\phi/\hbar$; using the latter instead would change the divergence prefactor.