In case of multivariate normal distribution of dimension n, the squared Mahalanobis distance is subject to a chi-squared distribution with n degrees of freedom.
So one approach here is:
(1) calculate the mean vector X and covariance matrix C
(2) calculate (u0-X)'*C^-1*(u0-X) and compare it against the critical value X^2(alpha, 3).