code: def pdinv(A, *args):
      
      L = jitchol(A, *args)
      logdet = 2.*np.sum(np.log(np.diag(L)))
      Li = dtrtri(L)
      Ai, _ = dpotri(L, lower=1)
      # Ai = np.tril(Ai) + np.tril(Ai,-1).T
      symmetrify(Ai)
  
      return Ai, L, Li, logdet
docs: :param A: A DxD pd numpy array
  
      :rval Ai: the inverse of A
      :rtype Ai: np.ndarray
      :rval L: the Cholesky decomposition of A
      :rtype L: np.ndarray
      :rval Li: the Cholesky decomposition of Ai
      :rtype Li: np.ndarray
      :rval logdet: the log of the determinant of A
      :rtype logdet: float64