docs: def getUs(classname, eta_s, f, alpha_s, alpha_e, m_u, m_d, m_s, m_c, m_b, m_e, m_mu, m_tau):
      """Get the QCD evolution matrix."""
      w, v = admeig(classname, f, m_u, m_d, m_s, m_c, m_b, m_e, m_mu, m_tau)
      b0s = 11 - 2 * f / 3
      a = w / (2 * b0s)
      return v @ np.diag(eta_s**a) @ np.linalg.inv(v)