The paper examines the variational equations of motion of a mechanical system of variable mass of the form $\ddot{\xi}^\gamma=A^\gamma_\delta(t)\xi^\gamma+B^\gamma_\delta(t)\dot{\xi}^\delta$, $(\xi,\delta=1,\dots,n)$. Here it is shown how a discrete model of a linear system $x(t_{n+1})=E(t_n)x(t_n)+F(t_n)U$ be used to solve the variations $\xi^\gamma(t)$.