code: def _genA(self):
          
          p, df = self._p, self.df
          A = np.zeros((p, p))
  
          for i in range(p):
              A[i, i] = sqrt(st.chi2.rvs(df - i))
  
          for j in range(p-1):
              for i in range(j+1, p):
                  A[i, j] = np.random.randn()
  
          return A
docs: Generate the matrix A in the Bartlett decomposition
  
            A is a lower triangular matrix, with
  
                A(i, j) ~ sqrt of Chisq(df - i + 1) when i == j
                        ~ Normal()                  when i > j