'Covariance matrix of the multinomial distribution.
Parameters
%(_doc_default_callparams)s
Returns
cov : ndarray
The covariance matrix of the distribution'

def cov(self, n, p):

    (n, p, npcond) = self._process_parameters(n, p)
    nn = n[..., np.newaxis, np.newaxis]
    result = (nn * np.einsum('...j,...k->...jk', (- p), p))
    for i in range(p.shape[(-1)]):
        result[..., i, i] += (n * p[..., i])
    return self._checkresult(result, npcond, np.nan)
