with(Groebner):

sf_ord_lower := proc()
	local eqs, nterms, basis, realroots, posroots, i;

	eqs := [
	  -S + Y + Z, # + F0 - G0 + F1 - G1,
	  -E + Y + Z*x^2 + S*x^4,
	  -Y + P,
	  -Z,

	  -P + x + E*x
	];

	nterms := [ P, Z, Y, E, S ];

	basis := gbasis(eqs, lexdeg(nterms[1..-2], [S]));
	print("sf ord lower basis computed", k);

	realroots := fsolve(denom(solve(basis[1] = 0, S)) = 0, x);

	posroots := select(type, [realroots], positive);
	print(posroots);

	1/min(seq(posroots[i], i=1..nops(posroots)));
end:

sf_rpn_lower := proc()
        local eqs, nterms, basis, realroots, posroots, i;

        eqs := [
          -S + Y + Z,
          -E + Y + Z + S*x^2,
          -Y + P,
          -Z,

          -P + x + E*x^2
        ];

	nterms := [ P, Z, Y, E, S ];

	basis := gbasis(eqs, lexdeg(nterms[1..-2], [S]));
	print("sf rpn lower basis computed", k);

	realroots := fsolve(denom(solve(basis[1] = 0, S)) = 0, x);

        posroots := select(type, [realroots], positive);
	print(posroots);

        1/min(seq(posroots[i], i=1..nops(posroots)));
end:

sf_alph_lower := proc()
        local eqs, nterms, basis, realroots, posroots, i;

        eqs := [
          -S + Y + Z,
          -E + Y + Z + S,
          -Y + P,
          -Z,

          -P + x + E*x
        ];

	nterms := [ P, Z, Y, E, S ];
										
	basis := gbasis(eqs, lexdeg(nterms[1..-2], [S]));
	print("sf alph lower basis computed", k);

	realroots := fsolve(denom(solve(basis[1] = 0, S)) = 0, x);

        posroots := select(type, [realroots], positive);
	print(posroots);

        1/min(seq(posroots[i], i=1..nops(posroots)));
end: