In this paper the authors study the problem of testing the hypothesis of a doubly exchangeable covariance matrix for three-level multivariate observations, taken on m variables over u sites and over v time/space points. Through the decomposition of the main hypothesis into a set of three sub-hypotheses, the likelihood ratio test statistic is defined, its exact moments are determined, and its exact distribution is studied. Because this distribution is very much intricate, a very precise near-exact distribution is developed. Numerical studies conducted to evaluate the closeness between this near-exact distribution and the exact distribution show the very good performance of this approximation even for very small sample sizes. A simulation study is also conducted and two real-data examples are presented.