The paper proposes a numerical technique within the Lagrangian description
for propagating the quantum fluid dynamical (QFD) equations in terms of the
Madelung field variables R and S, which are connected to the wave function
via the transformation Psi= exp{(R + iS)/(h) over bar}. The technique rest
s on the QFD equations depending only on the form, not the magnitude, of th
e probability density rho = \psi\(2) and on the structure of R = (h) over b
ar/2 In rho generally being simpler and smoother than rho. The spatially sm
ooth functions R and S are especially suitable for multivariate radial basi
s function interpolation to enable the implementation of a robust numerical
scheme. Examples of two-dimensional model systems show that the method riv
als, in both efficiency and accuracy, the split-operator and Chebychev expa
nsion methods. The results on a three-dimensional model system indicates th
at the present method is superior to the existing ones, especially, for its
low storage requirement and its uniform accuracy. The advantage of the new
algorithm is expected to increase fur higher dimensional systems to provid
e a practical computational tool.