module dipoly1; % Distributive polynomial algorithms.
% Authors: R. Gebauer, A. C. Hearn, H. Kredel.
% Modification for REDUCE > 3.3: H. Melenk.
fluid '(dipvars!* dipzero);
symbolic procedure dipconst!? p;
not dipzero!? p
and dipzero!? dipmred p
and evzero!? dipevlmon p;
symbolic procedure terprit n;
for i:=1:n do terpri();
symbolic procedure dfcprint pl;
% h polynomial factor list of distributive polynomials print.
for each p in pl do dfcprintin p;
symbolic procedure dfcprintin p;
% factor with exponent print.
( if cdr p neq 1 then << prin2 " ( "; dipprint1(p1,nil); prin2 " )** ";
prin2 cdr p; terprit 2 >> else << prin2 " "; dipprint p1>> )
where p1:= dipmonic a2dip prepf car p;
symbolic procedure dfcprin p;
% print content, factors and exponents of factorized polynomial p.
<< terpri(); prin2 " content of factorized polynomials = ";
prin2 car p;
terprit 2; dfcprint cdr p >>;
symbolic procedure diplcm p;
% Distributive polynomial least common multiple of denominators.
% p is a distributive rational polynomial. diplcm(p) calculates
% the least common multiple of the denominators and returns a
% base coefficient of the form 1/lcm(denom bc1,.....,denom bci).
if dipzero!? p then mkbc(1,1)
else bclcmd(diplbc p, diplcm dipmred p);
symbolic procedure diprectoint(p,u);
% Distributive polynomial conversion rational to integral.
% p is a distributive rational polynomial, u is a base coefficient
% ( 1/lcm denom p ). diprectoint(p,u) returns a primitive
% associate pseudo integral ( denominators are 1 ) distributive
% polynomial.
if bczero!? u then dipzero
else if bcone!? u then p
else diprectoint1(p,u);
symbolic procedure diprectoint1(p,u);
% Distributive polynomial conversion rational to integral internal 1.
% diprectoint1 is used in diprectoint.
if dipzero!? p then dipzero
else dipmoncomp(bclcmdprod(u,diplbc p),dipevlmon p,
diprectoint1(dipmred p,u));
symbolic procedure dipresul(p1,p2);
% test for fast downwards calculation
% p1 and p2 are two bivariate distributive polynomials.
% dipresul(p1,p2) returns the resultant of p1 and p2 with respect
% respect to the first variable of the variable list (car dipvars!*).
begin scalar q1,q2,q,ct;
q1:=dip2a diprectoint(p1,diplcm p1);
q2:=dip2a diprectoint(p2,diplcm p2);
ct := time();
q:= a2dip prepsq simpresultant list(q1,q2,car dipvars!*);
ct := time() - ct;
prin2 " resultant : "; dipprint dipmonic q; terpri();
prin2 " time resultant : "; prin2 ct; terpri();
end;
symbolic procedure dipbcprod (p,a);
% /* Distributive polynomial base coefficient product.
% p is a distributive polynomial, a is a base coefficient.
% dipbcprod(p,a) computes p*a, a distributive polynomial. */
if bczero!? a then dipzero
else if bcone!? a then p
else dipbcprodin(p,a);
symbolic procedure dipbcprodin (p,a);
% /* Distributive polynomial base coefficient product internal.
% p is a distributive polynomial, a is a base coefficient,
% where a is not equal 0 and not equal 1.
% dipbcprodin(p,a) computes p*a, a distributive polynomial. */
if dipzero!? p then dipzero
else dipmoncomp(bcprod(a, diplbc p),
dipevlmon p,
dipbcprodin(dipmred p, a));
symbolic procedure dipdif (p1,p2);
% /* Distributive polynomial difference. p1 and p2 are distributive
% polynomials. dipdif(p1,p2) calculates the difference of the
% two distributive polynomials p1 and p2, a distributive polynomial*/
if dipzero!? p1 then dipneg p2
else if dipzero!? p2 then p1
else ( if sl = 1 then dipmoncomp(diplbc p1,
ep1,
dipdif(dipmred p1, p2) )
else if sl = -1 then dipmoncomp(bcneg diplbc p2,
ep2,
dipdif(p1,dipmred p2))
else ( if bczero!? al
then dipdif(dipmred p1, dipmred p2)
else dipmoncomp(al,
ep1,
dipdif(dipmred p1,
dipmred p2) )
) where al = bcdif(diplbc p1, diplbc p2)
) where sl = evcomp(ep1, ep2)
where ep1 = dipevlmon p1, ep2 = dipevlmon p2;
symbolic procedure diplength p;
% /* Distributive polynomial length. p is a distributive
% polynomial. diplength(p) returns the number of terms
% of the distributive polynomial p, a digit.*/
if dipzero!? p then 0 else 1 + diplength dipmred p;
symbolic procedure diplistsum pl;
% /* Distributive polynomial list sum. pl is a list of distributive
% polynomials. diplistsum(pl) calculates the sum of all polynomials
% and returns a list of one distributive polynomial. */
if null pl or null cdr pl then pl
else diplistsum(dipsum(car pl, cadr pl) . diplistsum cddr pl);
symbolic procedure diplmerge (pl1,pl2);
% /* Distributive polynomial list merge. pl1 and pl2 are lists
% of distributive polynomials where pl1 and pl2 are in non
% decreasing order. diplmerge(pl1,pl2) returns the merged
% distributive polynomial list of pl1 and pl2. */
if null pl1 then pl2
else if null pl2 then pl1
else ( if sl >= 0 then cpl1 . diplmerge(cdr pl1, pl2)
else cpl2 . diplmerge(cdr pl2, pl1)
) where sl = evcomp(ep1, ep2)
where ep1 = dipevlmon cpl1, ep2 = dipevlmon cpl2
where cpl1 = car pl1, cpl2 = car pl2;
symbolic procedure diplsort pl;
% /* Distributive polynomial list sort. pl is a list of
% distributive polynomials. diplsort(pl) returns the
% sorted distributive polynomial list of pl.
sort(pl, function dipevlcomp);
symbolic procedure dipevlcomp (p1,p2);
% /* Distributive polynomial exponent vector leading monomial
% compare. p1 and p2 are distributive polynomials.
% dipevlcomp(p1,p2) returns a boolean expression true if the
% distributive polynomial p1 is smaller or equal the distributive
% polynomial p2 else false. */
not evcompless!?(dipevlmon p1, dipevlmon p2);
symbolic procedure dipmonic p;
% /* Distributive polynomial monic. p is a distributive
% polynomial. dipmonic(p) computes p/lbc(p) if p is
% not equal dipzero and returns a distributive
% polynomial, else dipmonic(p) returns dipzero. */
if dipzero!? p then p
else dipbcprod(p, bcinv diplbc p);
symbolic procedure dipneg p;
% /* Distributive polynomial negative. p is a distributive
% polynomial. dipneg(p) returns the negative of the distributive
% polynomial p, a distributive polynomial. */
if dipzero!? p then p
else dipmoncomp ( bcneg diplbc p,
dipevlmon p,
dipneg dipmred p );
symbolic procedure dipone!? p;
% /* Distributive polynomial one. p is a distributive polynomial.
% dipone!?(p) returns a boolean value. If p is the distributive
% polynomial one then true else false. */
not dipzero!? p
and dipzero!? dipmred p
and evzero!? dipevlmon p
and bcone!? diplbc p;
symbolic procedure dippairsort pl;
% /* Distributive polynomial list pair merge sort. pl is a list
% of distributive polynomials. dippairsort(pl) returns the
% list of merged and in non decreasing order sorted
% distributive polynomials. */
if null pl or null cdr pl then pl
else diplmerge(diplmerge( car(pl) . nil, cadr(pl) . nil ),
dippairsort cddr pl);
symbolic procedure dipprod (p1,p2);
% /* Distributive polynomial product. p1 and p2 are distributive
% polynomials. dipprod(p1,p2) calculates the product of the
% two distributive polynomials p1 and p2, a distributive polynomial*/
if diplength p1 <= diplength p2 then dipprodin(p1, p2)
else dipprodin(p2, p1);
symbolic procedure dipprodin (p1,p2);
% /* Distributive polynomial product internal. p1 and p2 are distrib
% polynomials. dipprod(p1,p2) calculates the product of the
% two distributive polynomials p1 and p2, a distributive polynomial*/
if dipzero!? p1 or dipzero!? p2 then dipzero
else ( dipmoncomp(bcprod(bp1, diplbc p2),
evsum(ep1, dipevlmon p2),
dipsum(dipprodin(dipfmon(bp1, ep1),
dipmred p2),
dipprodin(dipmred p1, p2) ) )
) where bp1 = diplbc p1,
ep1 = dipevlmon p1;
symbolic procedure dipprodls (p1,p2);
% /* Distributive polynomial product. p1 and p2 are distributive
% polynomials. dipprod(p1,p2) calculates the product of the
% two distributive polynomials p1 and p2, a distributive polynomial
% using distributive polynomials list sum (diplistsum). */
if dipzero!? p1 or dipzero!? p2 then dipzero
else car diplistsum if diplength p1 <= diplength p2
then dipprodlsin(p1, p2)
else dipprodlsin(p2, p1);
symbolic procedure dipprodlsin (p1,p2);
% /* Distributive polynomial product. p1 and p2 are distributive
% polynomials. dipprod(p1,p2) calculates the product of the
% two distributive polynomials p1 and p2, a distributive polynomial
% using distributive polynomials list sum (diplistsum). */
if dipzero!? p1 or dipzero!? p2 then nil
else ( dipmoncomp(bcprod(bp1, diplbc p2),
evsum(ep1, dipevlmon p2),
car dipprodlsin(dipfmon(bp1, ep1),
dipmred p2))
. dipprodlsin(dipmred p1, p2)
) where bp1 = diplbc p1,
ep1 = dipevlmon p1;
%symbolic procedure dipsum (p1,p2);
% /* Distributive polynomial sum. p1 and p2 are distributive
% polynomials. dipsum(p1,p2) calculates the sum of the
% two distributive polynomials p1 and p2, a distributive polynomial*/
%
% if dipzero!? p1 then p2
% else if dipzero!? p2 then p1
% else ( if sl = 1 then dipmoncomp(diplbc p1,
% ep1,
% dipsum(dipmred p1, p2) )
% else if sl = -1 then dipmoncomp(diplbc p2,
% ep2,
% dipsum(p1,dipmred p2))
% else ( if bczero!? al then dipsum(dipmred p1,
% dipmred p2)
% else dipmoncomp(al,
% ep1,
% dipsum(dipmred p1,
% dipmred p2) )
% ) where al = bcsum(diplbc p1, diplbc p2)
% ) where sl = evcomp(ep1, ep2)
% where ep1 = dipevlmon p1, ep2 = dipevlmon p2;
symbolic procedure dipsum (p1,p2);
% Distributive polynomial sum. p1 and p2 are distributive
% polynomials. dipsum(p1,p2) calculates the sum of the
% two distributive polynomials p1 and p2.
% Iterative version, better suited for very long polynomials.
% Warning: this routine uses "dipmred" == "cdr cdr" for a
% destructive concatenation.
begin scalar w,rw,sl,ep1,ep2,nt,al,done;
while not done do
<<if dipzero!? p1 then <<nt:=p2; done:=t>>
else
if dipzero!? p2 then <<nt:=p1; done:=t>>
else
<<ep1 := dipevlmon p1;
ep2 := dipevlmon p2;
sl := evcomp(ep1, ep2);
% Compute the next term.
if sl #= 1 then
<<nt := dipmoncomp(diplbc p1, ep1, nil);
p1 := dipmred p1>>
else
if sl #= -1 then
<<nt := dipmoncomp(diplbc p2, ep2, nil);
p2 := dipmred p2>>
else
<<al := bcsum(diplbc p1, diplbc p2);
nt := if not bczero!? al then dipmoncomp(al,ep1,nil);
p1 := dipmred p1; p2 := dipmred p2>>;
>>;
% Append the term to the sum polynomial.
if nt then
if null w then w:=rw:=nt
else <<cdr cdr rw := nt; rw := nt>>;
>>;
return w;
end;
endmodule;
end;