code: def orth_chol(order, dist, normed=True, sort="GR", cross_truncation=1., **kws):
      dim = len(dist)
      basis = chaospy.poly.basis(
          start=1, stop=order, dim=dim, sort=sort,
          cross_truncation=cross_truncation,
      )
      length = len(basis)
      cholmat = chaospy.chol.gill_king(chaospy.descriptives.Cov(basis, dist))
      cholmat_inv = numpy.linalg.inv(cholmat.T).T
      if not normed:
          diag_mesh = numpy.repeat(numpy.diag(cholmat_inv), len(cholmat_inv))
          cholmat_inv /= diag_mesh.reshape(cholmat_inv.shape)
      coefs = numpy.empty((length+1, length+1))
      coefs[1:, 1:] = cholmat_inv
      coefs[0, 0] = 1
      coefs[0, 1:] = 0
      expected = -numpy.sum(
          cholmat_inv*chaospy.descriptives.E(basis, dist, **kws), -1)
      coefs[1:, 0] = expected
      coefs = coefs.T
      out = {}
      out[(0,)*dim] = coefs[0]
      for idx in range(length):
          index = basis[idx].keys[0]
          out[index] = coefs[idx+1]
      polynomials = chaospy.poly.Poly(out, dim, coefs.shape[1:], float)
      return polynomials