def METHOD_NAME(nat_mean, nat_covar):
    try:
        L_inv = psd_safe_cholesky(-2.0 * nat_covar, upper=False)
    except RuntimeError as e:
        if str(e).startswith("cholesky"):
            raise RuntimeError(
                "Non-negative-definite natural covariance. You probably "
                "updated it using an optimizer other than gpytorch.optim.NGD (such as Adam). "
                "This is not supported."
            )
        else:
            raise e
    L = _triangular_inverse(L_inv, upper=False)
    S = L.transpose(-1, -2) @ L
    mu = (S @ nat_mean.unsqueeze(-1)).squeeze(-1)
    # Two choleskys are annoying, but we don't have good support for a
    # LinearOperator of form L.T @ L
    return mu, psd_safe_cholesky(S, upper=False)