Skip to content

Fix modular reduction during PolyRatFun expansion #921

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
cbmarini wants to merge 1 commit into
base: master
Choose a base branch
from
Open

Fix modular reduction during PolyRatFun expansion #921

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
126 changes: 126 additions & 0 deletions check/fixes.frm
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5613,3 +5613,129 @@ Print;
assert succeeded?
assert result("test") =~ expr("0")
*--#] PullReq860_3 :
*--#[ PolyRatFunExpandModulus :
* Reduce expansion coefficients before sorting, including uncombined terms.
Modulus 13;
CFunction rat;
Symbol ep;
PolyRatFun rat(expand,ep,3);
Local F = rat(1,1+ep/2);
Local G = rat(1+ep,2+ep);
Local H = rat(1+ep,2);
Local N = rat(1,2);
Local S = rat(ep,2);
Local D = rat(1,2*ep);
Local L = rat(ep^2,2*ep+ep^2);
Local P2 = 6*rat(7);
* Continue beyond the quadratic coefficient that vanishes modulo 13.
Local R = rat(1,1+4*ep+3*ep^2);
Print +s;
.end
Modulus PlusMin,13;
CFunction rat;
Symbol ep;
PolyRatFun rat(expand,ep,3);
Local F = rat(1,1+ep/2);
Local G = rat(1+ep,2+ep);
Local H = rat(1+ep,2);
Local N = rat(1,2);
Local S = rat(ep,2);
Local D = rat(1,2*ep);
Local L = rat(ep^2,2*ep+ep^2);
Local P2 = 6*rat(7);
Local R = rat(1,1+4*ep+3*ep^2);
Print +s;
.end
assert succeeded?
modes = ["Modulus 13", "PlusMin 13"]
[
["F", "+ rat(1 + 6*ep + 10*ep^2 + 8*ep^3)", "+ rat(1 + 6*ep - 3*ep^2 - 5*ep^3)"],
["G", "+ rat(7 + 10*ep + 8*ep^2 + 9*ep^3)", "+ rat(-6 - 3*ep - 5*ep^2 - 4*ep^3)"],
["H", "+ rat(7 + 7*ep)", "+ rat(-6 - 6*ep)"],
["N", "+ rat(7)", "+ rat(-6)"],
["S", "+ rat(7*ep)", "+ rat(-6*ep)"],
["D", "+ rat(7*ep^-1)", "+ rat(-6*ep^-1)"],
["L", "+ rat(7*ep + 3*ep^2 + 5*ep^3 + 4*ep^4)", "+ rat(-6*ep + 3*ep^2 + 5*ep^3 + 4*ep^4)"],
["P2", "+ rat(3)", "+ rat(3)"],
["R", "+ rat(1 + 9*ep + 12*ep^3)", "+ rat(1 - 4*ep - ep^3)"]
].each do |name, *values|
values.each_with_index do |value, index|
assert result(name, index) =~ expr(value), "#{modes[index]}: #{name}"
end
end
*--#] PolyRatFunExpandModulus :
*--#[ PolyFunModulusCoefficient :
* Compact numeric arguments need the same reduction as general arguments.
Modulus 13;
CFunction pf;
PolyFun pf;
Local F = 6*pf(7);
Local G = 1/2*pf(7);
Local H = -6*pf(7);
Print +s;
.end
Modulus PlusMin,13;
CFunction pf;
PolyFun pf;
Local F = 6*pf(7);
Local G = 1/2*pf(7);
Local H = -6*pf(7);
Print +s;
.end
* A nonzero product can vanish when the modulus is composite.
Modulus 15;
CFunction pf;
PolyFun pf;
Local F = 3*pf(5);
Print;
.end
assert succeeded?
modes = ["Modulus 13", "PlusMin 13"]
[
["F", "+ pf(3)", "+ pf(3)"],
["G", "+ pf(10)", "+ pf(-3)"],
["H", "+ pf(10)", "+ pf(-3)"]
].each do |name, *values|
values.each_with_index do |value, index|
assert result(name, index) =~ expr(value), "#{modes[index]}: #{name}"
end
end
assert result("F", 2) =~ expr("0"), "Modulus 15: F"
*--#] PolyFunModulusCoefficient :
*--#[ DollarFactorizeModulusContent :
* Content extraction must reduce coefficients before sorting.
Modulus AlsoDollars,13;
Symbol x;
Local F = 1;
$p = 2*x+4*x^2;
FactDollar $p;
Print " N = %$;",$p[0];
Print " F1 = %$;",$p[1];
Print " F2 = %$;",$p[2];
Print " F3 = %$;",$p[3];
ModuleOption noparallel;
.end
Modulus PlusMin,AlsoDollars,13;
Symbol x;
Local F = 1;
$p = 2*x+4*x^2;
FactDollar $p;
Print " N = %$;",$p[0];
Print " F1 = %$;",$p[1];
Print " F2 = %$;",$p[2];
Print " F3 = %$;",$p[3];
ModuleOption noparallel;
.end
assert succeeded?
modes = ["Modulus 13", "PlusMin 13"]
[
["N", "3", "3"],
["F1", "x", "x"],
["F2", "x + 7", "x - 6"],
["F3", "4", "4"]
].each do |name, *values|
values.each_with_index do |value, index|
assert result(name, index) =~ expr(value), "#{modes[index]}: #{name}"
end
end
*--#] DollarFactorizeModulusContent :
12 changes: 10 additions & 2 deletions doc/manual/statements.tex
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3526,11 +3526,14 @@ \section{modulus}
The options are
\begin{description}
\item[NoFunctions] Modulus calculus is not performed inside function
arguments.
arguments, except for the active PolyFun (see \ref{substapolyfun}) and the
coefficients of the active PolyRatFun in expand mode (see
\ref{substapolyratfun}).
\item[AlsoFunctions] Modulus calculus is also performed inside function
arguments.
\item[CoefficientsOnly] Modulus calculus is neither performed inside function
arguments nor on powers of symbols.
arguments nor on powers of symbols. The exceptions listed under NoFunctions
also apply.
\item[PlusMin] The values of numbers are reduced to the range
$(-value+1)/2$ to $(value-1)/2$.
\item[Positive] The values of numbers are reduced to the range $0$ to
Expand Down Expand Up @@ -4632,6 +4635,11 @@ \section{polyratfun}
addition of two coefficients is a simple and straightforward operation that
does not need the expensive polynomial GCD computations.

\noindent When Modulus (see \ref{substamodulus}) is active, the coefficients
of this expansion are reduced modulo the specified value. This also applies
with the default NoFunctions option: the active PolyRatFun in expand mode
is exempt from this restriction, so AlsoFunctions is not required.

\noindent Of course one can program such expansions externally and maybe
better suited for the problem at hand, but using this option of the
polyratfun is much faster and gives fewer chances of mistakes.
Expand Down
2 changes: 2 additions & 0 deletions sources/proces.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5172,6 +5172,8 @@ int PrepPoly(PHEAD WORD *term,WORD par)
if ( m >= AT.WorkSpace && m < AT.WorkTop )
AT.WorkPointer = m;
if ( Normalize(BHEAD v) ) Terminate(-1);
/* Modulo reduction of the argument */
if ( AN.ncmod != 0 && Modulus(v) ) Terminate(-1);
AT.WorkPointer = oldworkpointer;
m = w;
if ( *v == 4 && v[2] == 1 && (v[1]&MAXPOSITIVE) == v[1] ) {
Expand Down
43 changes: 42 additions & 1 deletion sources/ratio.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1340,6 +1340,16 @@ WORD *MultiplyWithTerm(PHEAD WORD *in, WORD *term, WORD par)
*termout = tt - termout;
// in or term can contain non-symbols: we call this in TakeContent
Normalize(BHEAD termout);
/* Modulo reduction before sorting */
if ( AN.ncmod != 0 ) {
if ( *termout != 0 && Modulus(termout) ) {
LowerSortLevel();
AR.CompareRoutine = (COMPAREDUMMY)oldcompareroutine;
AR.SortType = oldsorttype;
goto CalledFrom;
}
if ( *termout == 0 ) { in += *in; continue; }
}
StoreTerm(BHEAD termout);
in += *in;
}
Expand Down Expand Up @@ -3221,6 +3231,10 @@ NormArg:;
rr = r + *r;
i = *r; rrr = outarg; NCOPY(rrr,r,i);
Normalize(BHEAD outarg);
if ( AN.ncmod != 0 && *outarg != 0 && Modulus(outarg) ) {
LowerSortLevel(); AT.TrimPower = 0;
error = 3; goto onerror;
}
if ( *outarg > 0 ) StoreTerm(BHEAD outarg);
}
r = fun+FUNHEAD+ARGHEAD;
Expand Down Expand Up @@ -3301,6 +3315,10 @@ NormArg:;
}
}
else { error = 2; goto onerror; } /* should not happen! */
if ( AN.ncmod != 0 && *rr != 0 ) {
if ( Modulus(rr) ) { error = 3; goto onerror; }
r = rr + *rr; *r = 0;
}
}
else { /* Multi-term numerator. */
m = arg1+ARGHEAD;
Expand Down Expand Up @@ -3332,6 +3350,12 @@ NormArg:;
m++; while ( m < mm ) *r++ = *m++;
*rrr = r-rrr;
Normalize(BHEAD rrr);
if ( AN.ncmod != 0 ) {
if ( *rrr != 0 && Modulus(rrr) ) {
LowerSortLevel(); error = 3; goto onerror;
}
if ( *rrr == 0 ) continue;
}
StoreTerm(BHEAD rrr);
}
EndSort(BHEAD rr,1);
Expand Down Expand Up @@ -3413,6 +3437,12 @@ NormArg:;
if ( r < AT.WorkTop && r >= AT.WorkSpace )
AT.WorkPointer = r;
Normalize(BHEAD rrr);
if ( AN.ncmod != 0 ) {
if ( *rrr != 0 && Modulus(rrr) ) {
LowerSortLevel(); error = 3; goto onerror;
}
if ( *rrr == 0 ) continue;
}
if ( ABS(rrr[*rrr-1]) == *rrr-1 ) {
if ( AR.PolyFunPow >= 0 ) {
StoreTerm(BHEAD rrr);
Expand All @@ -3425,6 +3455,11 @@ NormArg:;
}
EndSort(BHEAD rr,1);
}
/* Fast numerators bypass sorting, but still need reduction. */
if ( *arg1 < 0 && AN.ncmod != 0 && *rr != 0 ) {
if ( Modulus(rr) ) { error = 3; goto onerror; }
rr[*rr] = 0;
}
r = rr; while ( *r ) r += *r;
i = r-rr;
r = fun + FUNHEAD + ARGHEAD;
Expand Down Expand Up @@ -3590,6 +3625,9 @@ NormArg:;
*thecopy = rr - thecopy;
AT.WorkPointer = rr;
Normalize(BHEAD thecopy);
if ( AN.ncmod != 0 && *thecopy != 0 && Modulus(thecopy) ) {
LowerSortLevel(); error = 3; goto onerror;
}
if ( *thecopy > 0 ) StoreTerm(BHEAD thecopy);
AT.WorkPointer = thecopy;
}
Expand Down Expand Up @@ -3725,8 +3763,11 @@ int InvPoly(PHEAD WORD *inpoly, WORD maxpow, WORD sym)
}
}
/*
Copy c to the proper location
Reduce before reusing c, and represent modular zeroes by null pointers.
Copy c to the proper location.
*/
if ( AN.ncmod != 0 && lenc != 0 &&
TakeModulus((UWORD *)c,&lenc,AC.cmod,AN.ncmod,UNPACK) ) goto calcerror;
if ( lenc == 0 ) AT.pWorkSpace[outpointers+j] = 0;
else {
AT.pWorkSpace[outpointers+j] = w;
Expand Down
Loading