Finding Values of All Nine Constants in q-hypergeometric Equation

I suspect that an identity of the form $$\sum_{j=-n}^{j=n}(-1)^jq^{j(5j+1)/2}=\sum_{r=0}^{n}f_{ar+bn+c}(-1)^{n-r}q^{(n-r)(dr+en+f)/2}\frac{(q;q)_{gr+hn+j}}{(q;q)_{n-r}}$$ for $a, \dots , j$ integer constants, where $$f_N=\sum_{j=0}^{N}{N \brack j}q^{j^2}.$$ Here $(a;q)_n$ and $a \brack b$ denote the standard q-pochammer symbol and q binomial coefficient, resp. I turned to Mathematica for help in determining exactly what these integer constants must be. The bulk of my code is as follows:

f[n_] := Expand[Sum[FunctionExpand[QBinomial[n, j, q]]*q^(j^2), {j, 0, n}]]

lhs[n_] := Sum[(-1)^j*q^(j (5 j + 1)/2), {j, -n, n}]


The following implementation worked for 3 variables, but not for 6 variables:

rhs[n_, a_, b_, c_, d_, e_, f_] :=
Expand[Sum[f[r]*(-1)^{n - r}*q^{(n - r) (a*n + b*r + c)/2}*
FunctionExpand[QPochhammer[q, q, d*n + e*r + f]/QPochhammer[q, q, n - r]], {r,
0, n}]]

L[n_] :=
Do[
If[CoefficientList[lhs[n], q] == CoefficientList[rhs[n, a_, b_, c_, d_, e_, f_], q],
Print[{a, b, c}]], {a, t, -t + 1}, {b, t, -t + 1}, {c, t, -t + 1}, {d, 1, 3},
{e, 0, 3}, {f, 0, 3}]


It's also fairly sluggish in 3 variables.

A New Approach

Following George2079's suggestion, I now seek polynomials $g(r,n), h(r,n)$ in $q$ of fixed degree for which $$\sum_{j=-n}^{j=n}(-1)^jq^{j(5j+1)/2}=\sum_{r=0}^{n}f_{g(r,n)}\cdot h(r,n).$$ I am almost certain that such an $h(r,n)$ exists for all $n$ when $g(r,n)=ar+bn$. Then the following is a preliminary function on the assumption that $g(r,n)=r$, and finds something for $n\le8$:

Guess[n_] :=
Module[{maxr = Floor[Sqrt[n (5 n + 1)/2]], maxk = 20, mat, ns},
mat = Transpose[
PadRight[CoefficientList[#, q] & /@
Append[Flatten[Table[q^k*f[r], {r, 0, maxr}, {k, 0, maxk}]], lhs[n]]]];
ns = Select[NullSpace[mat], Last[#] =!= 0 &];
If[Length[ns] =!= 1, Throw["Error"]];
Return[((#.q^(Range[Length[#]] - 1)) & /@ Partition[Most[First[ns]], maxk + 1]).
Table[F[r], {r, 0, maxr}]];
];

Guess[5]
(-1 + q^2 + q^3 + q^8 + q^9 - 3 q^10 - 8 q^11 + q^12 + 12 q^13 +
7 q^14 - 5 q^15 - 12 q^16 + 10 q^17 + 16 q^18 + 3 q^19 -
5 q^20) F[0] -
13 q^20 F[1] + (19 q^19 + 11 q^20) F[
2] + (-8 q^16 - 17 q^17 - 17 q^18 - 9 q^19 + 5 q^20) F[
3] + (9 q^14 + 18 q^15 + 18 q^16 + 8 q^17 + 5 q^18 - 2 q^19 -
10 q^20) F[
4] + (-9 q^12 - 11 q^13 - 9 q^14 - 8 q^15 + q^17 + 3 q^18 +
4 q^19 + 6 q^20) F[
5] + (5 q^10 + 3 q^11 + 2 q^12 + q^13 + q^14 - q^15 - 2 q^16 -
2 q^17 - q^18 - q^19 - q^20) F[6] + (-q^8 - q^9 + q^11 + q^16) F[
7]

Expand[% /. F -> f]
-1 + q^2 + q^3 - q^9 - q^11 + q^21 + q^24 - q^38 - q^42 + q^60 + q^65

tab = {}; ii = 1;
While[True, AppendTo[tab, Guess[ii]]; ii++; Print[ii]]
2
3
4
5
6
7
8


And one gets Hold[Throw[Error]]. Can someone produce an implementation that caters well to guessing what the general $h(r,n)$ might look like/be?

-
You used Set (=) where you need Equal (==). Look at the two If statements in the code. – Szabolcs May 27 '14 at 17:39
Unrelated: instead of three nested Fors consider using a single Do with 3 iterators. Further improvement: instead of Print use Sow and Reap to collect results. – Szabolcs May 27 '14 at 17:41
Please review this post: mathematica.stackexchange.com/a/18395/12 – Szabolcs May 27 '14 at 17:43
@Szabolcs Thanks for this. I'm looking thru it, as it contains a lot of useful information. Can you clarify exactly how you mean to implement your Do loop with 3 iterators? Then I can work more efficiently to learn the relevant functions, as I am starting to see that there is a wealth of functions, and it is difficult to navigate which ones might be relevant to the problem at hand. – Jonny May 28 '14 at 0:35
Similar to this: Do[..., {i,1,10}, {j,1,5}, ...] – Szabolcs May 28 '14 at 1:31