FindQ := proc(FL, AL)
	global L,H,y,U;
	local i, A, F, C, P, CY;
	A := mul(1+y*exp(-i), i = AL);
	F := mul(1+y*exp(-i), i = FL);
	C := F/((1+y)*A);
	P := convert(FL,`*`);
	CY := P * eval(A/F, y = -1);
	PiStar(factor(normal(C * CY, expanded)), P);
	map(factor, convert(subs(exp(L) = 1/U, %),parfrac,U) );
end:

# Evaluate C at H = point.  Use L'Hopital if necessary.
EVAL_AT := proc( C,  point )
	global H;
	local F,G;
	F := normal(C);
	F := [numer(F), denom(F)];
	G := expand( eval(F, H = point) );
	while testeq(G[2]) do
		if G[1] <> 0 then
			error "non-zero divided by zero"
		fi;
		F := [diff(F[1],H), diff(F[2],H)];
		G := expand( eval(F, H = point) )
	od;
	factor(G[1]/G[2])
end:

# The input an is a list where each entry is of the form [co, i, point]
# The program evaluates the i'th derivative of F at the point H = point,
# multiplies that by co, and adds this up for all members of the list an.
Apply_AN := proc(an, F)
	local i;
	add(`if`(i[1]=0, 0, i[1]*EVAL_AT(DIFF(F,i[2]),i[3])), i=an )
end:

# i'th derivative of F w.r.t. H
DIFF := proc(F,i)
	global H;
	options remember;
	if i=0 then
		F
	else
		normal( diff(procname(F,i-1), H) )
	fi
end:

# The PiStar(C, P) program returns Apply_AN(an, C) for a suitable list an.
#
# See Apply_AN for notations. For each [co, i, point] in the list an, we can read
# off the possible i's and points from the list P.  Initially, we'll use unknowns
# for the co's. Next, we apply Apply_AN(an, H^i) for a number of i's, and compare
# with what the result should have been.  This gives a set of equations, from which
# the unknown co's can be computed.
PiStar := proc(C, P)
	global H;
	local A, V, an, i,j, co, eq;
	if type(P,list) then
		A := convert(P,`*`)
	else
		A := P
	fi;
	V := factors(A * H^(degree(A,H)-ldegree(A,H)) )[2];
	an := [seq(seq([co[i,j], j, solve(V[i][1],H)], j=0..V[i][2]-1), i=1..nops(V))];
	eq := {seq(Apply_AN(an, H^i) + residue(H^i/A, H=infinity), i=0..nops(an))};
	an := subs(solve(eq, {seq(i[1],i=an)}), an);
	Apply_AN(an,C)
end:

# These are the Q's listed in the paper:
#
Q_D5 := FindQ([H, H+L, H+L, H+L],   [2*H + 2*L,  2*H + 2*L]);
Q_E6 := FindQ([H, H+L, H+L], [3*H + 3*L]);
Q_E7 := FindQ([H, H+L, H+2*L, H+2*L], [2*H + 2*L, 2*H + 4*L]);
Q_E8 := FindQ([H, H + 2*L, H + 3*L], [3*H + 6*L]);

# The program can do the same for other cases, for example:
Q_E7star := FindQ([H,H,H+L], [3*H+2*L]);