code: def rwishart_cov(n, C):
      
      # return rwishart(n, np.linalg.inv(C))
  
      p = np.shape(C)[0]
      # Need cholesky decomposition of precision matrix C^-1?
      sig = np.linalg.cholesky(C)
  
      if n <= (p-1):
          raise ValueError('Wishart parameter n must be greater '
                           'than size of matrix.')
  
      norms = np.random.normal(size=(p * (p - 1)) // 2)
      chi_sqs = np.sqrt(np.random.chisquare(df=np.arange(n, n - p, -1)))
      A = flib.expand_triangular(chi_sqs, norms)
  
      flib.dtrmm_wrap(sig, A, side='L', uplo='L', transa='N', alpha=1.)
      w = np.asmatrix(np.dot(A, A.T))
      flib.symmetrize(w)
      return w
docs: Return a Wishart random matrix.
  
      :Parameters:
        n : int
          Degrees of freedom, > 0.
        C : matrix
          Symmetric and positive definite