As far as I know, the only significant result to speed up these calculations is that $E_2(\frac{p}{q}) = \frac{1}{2}|\lbrace d: d | q^2, d \equiv -q (mod p) \rbrace|$, where $E_2(p/q)$ represents the number of decompositions into 2 unit fractions, and each matching $d$ represents the decomposition $\frac{p}{q} = \frac{qp}{q(q+d)} + \frac{dp}{q(q+d)}$. (Take floor() or ceil() depending on whether you want to allow repeats.)

When I've coded this in the past, I called one of 4 different functions depending on a) whether $p=1$ or not, and b) whether $q/p \ge min$ or not, where $min$ is the greatest denominator I'm already using. When $p=1$ and $q \ge min$, in particular, we can just calculate $\tau(q^2)/2$ from the factorisation of $q$; in the other cases I actually walked the factors from $q/p$ to $\sqrt{q}$.

So: yes, you can count the number of matching sets without generating the 7 elements of each set, but computationally the elements are just a whisker away.

Comput. Math.69(1999), no. 229, 407--420]. – Wadim Zudilin Jun 12 '10 at 11:23