rrspace for nodal curves

parent b97bcb23
...@@ -61,19 +61,22 @@ vanish at singular points of the curve (this case cannot occur if the curve is ...@@ -61,19 +61,22 @@ vanish at singular points of the curve (this case cannot occur if the curve is
smooth). In this case, please use another model of the curve. smooth). In this case, please use another model of the curve.
1. `rrspace`: the program `rrspace` computes the Riemann-Roch space associated 1. `rrspace`: the program `rrspace` computes the Riemann-Roch space associated
to a given divisor D on a curve q. It reads its input on the standard input, to a given divisor D on a smooth or nodal plane curve described by a
and it returns a basis of the Riemann-Roch space on the standard output. bivariate polynomial q. It reads its input on the standard input, and it
returns a basis of the Riemann-Roch space on the standard output.
The format for its input (see next section for more details) is The format for its input (see next section for more details) is
``` ```
<PRIME> <PRIME>
<CURVE> <CURVE>
<EFFECTIVEDIVISOR>
<DIVISOR> <DIVISOR>
``` ```
where `<DIVISOR>` represents a divisor whose support does not contain any where `<DIVISOR>` represents a divisor whose support does not contain any
singular point of the curve. singular point of the curve, and `<EFFECTIVEDIVISOR>` represents the effective
divisor of all the nodes of the curve.
It outputs the following data: It outputs the following data:
* the dimension D of the Riemann-Roch space; * the dimension D of the Riemann-Roch space;
...@@ -84,20 +87,22 @@ smooth). In this case, please use another model of the curve. ...@@ -84,20 +87,22 @@ smooth). In this case, please use another model of the curve.
Riemann-Roch space. Riemann-Roch space.
2. `class_grp_arith`: this program takes as input two elements of the divisor 2. `class_grp_arith`: this program takes as input two elements of the divisor
class group of a curve `C` of genus `g` and it returns a representation of the class group of a nodal curve `C` of genus `g` with `r` nodes and it returns
sum. Here, we are given a curve, together with a degree 1 divisor on `C`. A a representation of the sum. Here, we are given a curve, together with a degree
class of divisors (modulo principal divisors) is represented by an effective 1 divisor on `C`. A class of divisors (modulo principal divisors) is
divisor `D1` of degree `g`: it represents the class of divisors equivalent to represented by an effective divisor `D1` of degree `g` such that its support
the degree-0 divisor `D1 - g*O`. contains no singular point of the curve: it represents the class of divisors
equivalent to the degree-0 divisor `D1-g*O`.
`class_grp_arith` reads its input on the standard input, and it returns an `class_grp_arith` reads its input on the standard input, and it returns an
effective divisor of degree at most `g` on the standard output. effective divisor of degree `g+r` on the standard output.
The format for its input (see next section for more details) is The format for its input (see next section for more details) is
``` ```
<PRIME> <PRIME>
<CURVE> <CURVE>
<EFFECTIVEDIVISOR>
<GENUS> <GENUS>
<DIVISOR> <DIVISOR>
<EFFECTIVEDIVISOR> <EFFECTIVEDIVISOR>
...@@ -105,12 +110,15 @@ smooth). In this case, please use another model of the curve. ...@@ -105,12 +110,15 @@ smooth). In this case, please use another model of the curve.
``` ```
The divisor provided on the 4th line must have degree 1. The two effective The divisor provided on the 4th line must have degree 1. The two effective
divisors `D1` and `D2` provided on the last two lines must have degree `<GENUS>`. divisors `D1` and `D2` provided on the last two lines must have degree
The supports of all the divisors must not contain any singular point of the `<GENUS>`. The first effective divisor describes all the nodes of the curve.
curve. The supports of the three last divisors must not contain any singular point of
the curve.
The program returns a representation of a degree-g effective divisor `D3` such that
`(D1-g*O) + (D2-g*O) = (D3-g*O)` in the divisor class group of `C`. The program returns a representation of a degree-`(g+r)` effective divisor
`D3 ≥ E` such that `(D1-g*O) + (D2-g*O) = (D3-E-g*O)` in the divisor class
group of `C`, where `E` is the degree-`r` effective divisor representing the
nodes of the curve.
3. `tests`: this software does not take any input. It just runs two tests: one 3. `tests`: this software does not take any input. It just runs two tests: one
for the computation of Riemann-Roch spaces, and one for the arithmetic in the for the computation of Riemann-Roch spaces, and one for the arithmetic in the
......
...@@ -118,7 +118,8 @@ Interpolate(const EffectiveDivisor& D) { ...@@ -118,7 +118,8 @@ Interpolate(const EffectiveDivisor& D) {
} }
EffectiveDivisor EffectiveDivisor
PrincipalDivisor(const Curve& C, const BivPol& h) { PrincipalDivisor(const Curve& C, const BivPol& h, const EffectiveDivisor& E) {
assert(E.curve() == C);
size_t deg_res = h.degree() * C.degree(); // Bézout bound size_t deg_res = h.degree() * C.degree(); // Bézout bound
vec_ZZ_p eval_pts; vec_ZZ_p eval_pts;
for (long i = 0; i <= (long)deg_res; ++i) for (long i = 0; i <= (long)deg_res; ++i)
...@@ -150,7 +151,7 @@ PrincipalDivisor(const Curve& C, const BivPol& h) { ...@@ -150,7 +151,7 @@ PrincipalDivisor(const Curve& C, const BivPol& h) {
ZZ_pX subres0 = interpolate(eval_pts, eval_subres0); ZZ_pX subres0 = interpolate(eval_pts, eval_subres0);
ZZ_pX subres1 = interpolate(eval_pts, eval_subres1); ZZ_pX subres1 = interpolate(eval_pts, eval_subres1);
return EffectiveDivisor(C, res, -MulMod(InvMod(subres1, res), subres0, res)); return PositiveDifference(EffectiveDivisor(C, res, -MulMod(InvMod(subres1, res), subres0, res)), E);
} }
EffectiveDivisor EffectiveDivisor
...@@ -161,6 +162,29 @@ PositiveDifference(const EffectiveDivisor& D1, const EffectiveDivisor& D2) { ...@@ -161,6 +162,29 @@ PositiveDifference(const EffectiveDivisor& D1, const EffectiveDivisor& D2) {
return EffectiveDivisor(D1.curve(), f/LeadCoeff(f), D1.get_g() % f); return EffectiveDivisor(D1.curve(), f/LeadCoeff(f), D1.get_g() % f);
} }
// This simple sum only works if D1 and D2 do not have any point with the same
// x-coordinate in their support
EffectiveDivisor
SimpleSum(const EffectiveDivisor& D1, const EffectiveDivisor& D2) {
assert(D1.curve() == D2.curve());
if (deg(D1.get_f()) == 0)
return D2;
if (deg(D2.get_f()) == 0)
return D1;
NTL::ZZ_pX gcd, a1, a2;
XGCD(gcd, a1, a2, D1.get_f(), D2.get_f());
assert(deg(gcd) == 0);
ZZ_pX new_f = D1.get_f()*D2.get_f();
ZZ_pX new_g = (D1.get_g()*a2*D2.get_f() +
D2.get_g()*a1*D1.get_f()) % new_f;
assert(new_g % D1.get_f() == D1.get_g());
assert(new_g % D2.get_f() == D2.get_g());
return EffectiveDivisor(D1.curve(), new_f, new_g);
}
// This adding function works only if all points in the support are
// non-singular
EffectiveDivisor EffectiveDivisor
Sum(const EffectiveDivisor& D1, const EffectiveDivisor& D2) { Sum(const EffectiveDivisor& D1, const EffectiveDivisor& D2) {
assert(D1.curve() == D2.curve()); assert(D1.curve() == D2.curve());
...@@ -175,6 +199,7 @@ Sum(const EffectiveDivisor& D1, const EffectiveDivisor& D2) { ...@@ -175,6 +199,7 @@ Sum(const EffectiveDivisor& D1, const EffectiveDivisor& D2) {
ZZ_pX new_f = D1.get_f()*D2.get_f(); ZZ_pX new_f = D1.get_f()*D2.get_f();
ZZ_pX new_g = (D1.get_g()*a2*(D2.get_f()/gcd) + ZZ_pX new_g = (D1.get_g()*a2*(D2.get_f()/gcd) +
D2.get_g()*a1*(D1.get_f()/gcd)) % (new_f/gcd); D2.get_g()*a1*(D1.get_f()/gcd)) % (new_f/gcd);
new_g = NewtonHenselStep(*D1.curve().get_pdefpol(), new_g, new_f/gcd); new_g = NewtonHenselStep(*D1.curve().get_pdefpol(), new_g, new_f/gcd);
new_g = new_g % new_f; new_g = new_g % new_f;
assert(new_g % D1.get_f() == D1.get_g()); assert(new_g % D1.get_f() == D1.get_g());
...@@ -226,7 +251,7 @@ EffectiveBasisRRinDegree(const EffectiveDivisor& D, size_t deg) { ...@@ -226,7 +251,7 @@ EffectiveBasisRRinDegree(const EffectiveDivisor& D, size_t deg) {
ZZ_pX(0)); ZZ_pX(0));
for (size_t i = 0; i < vec_mons.size(); ++i) { for (size_t i = 0; i < vec_mons.size(); ++i) {
pair<size_t, size_t> p = vec_mons[i]; pair<size_t, size_t> p = vec_mons[i];
interp_pol[p.second] += kerM.get(k, i)* (ZZ_pX(1) << p.first); interp_pol[p.second] += kerM.get(k, i)*(ZZ_pX(1) << p.first);
} }
while (IsZero(interp_pol.back())) while (IsZero(interp_pol.back()))
interp_pol.pop_back(); interp_pol.pop_back();
...@@ -237,12 +262,14 @@ EffectiveBasisRRinDegree(const EffectiveDivisor& D, size_t deg) { ...@@ -237,12 +262,14 @@ EffectiveBasisRRinDegree(const EffectiveDivisor& D, size_t deg) {
} }
RRspace RRspace
RiemannRochBasis(const Divisor& D) { RiemannRochBasis(const Divisor& D, const EffectiveDivisor& E) {
assert(D.curve().get_pdefpol()->mod_eval(D.get_pos().get_g(), D.get_pos().get_f()) == 0); assert(D.curve().get_pdefpol()->mod_eval(D.get_pos().get_g(), D.get_pos().get_f()) == 0);
assert(D.curve().get_pdefpol()->mod_eval(D.get_neg().get_g(), D.get_neg().get_f()) == 0); assert(D.curve().get_pdefpol()->mod_eval(D.get_neg().get_g(), D.get_neg().get_f()) == 0);
BivPol h = Interpolate(D.get_pos()); assert(D.curve().get_pdefpol()->mod_eval(E.get_g(), E.get_f()) == 0);
EffectiveDivisor Dp = PrincipalDivisor(D.curve(), h); BivPol h = Interpolate(SimpleSum(D.get_pos(), E));
EffectiveDivisor Dp = PrincipalDivisor(D.curve(), h, E);
EffectiveDivisor Dp2 = PositiveDifference(Dp, D.get_pos()); EffectiveDivisor Dp2 = PositiveDifference(Dp, D.get_pos());
EffectiveDivisor newDneg = Sum(Dp2, D.get_neg()); EffectiveDivisor Dp3 = PositiveDifference(Dp2, E);
EffectiveDivisor newDneg = SimpleSum(Sum(Dp3, D.get_neg()), E);
return RRspace(h, EffectiveBasisRRinDegree(newDneg, h.degree()), newDneg); return RRspace(h, EffectiveBasisRRinDegree(newDneg, h.degree()), newDneg);
} }
...@@ -39,7 +39,7 @@ BuildInterpolationMatrixInDegree(const EffectiveDivisor& D, std::size_t deg); ...@@ -39,7 +39,7 @@ BuildInterpolationMatrixInDegree(const EffectiveDivisor& D, std::size_t deg);
// Given a form h on the curve, compute the principal effective divisor (h) // Given a form h on the curve, compute the principal effective divisor (h)
EffectiveDivisor EffectiveDivisor
PrincipalDivisor(const Curve& C, const BivPol& h); PrincipalDivisor(const Curve& C, const BivPol& h, const EffectiveDivisor& E);
EffectiveDivisor EffectiveDivisor
PositiveDifference(const EffectiveDivisor& D1, const EffectiveDivisor& D2); PositiveDifference(const EffectiveDivisor& D1, const EffectiveDivisor& D2);
...@@ -55,6 +55,9 @@ MultiplyByInt(const EffectiveDivisor& D, unsigned int k); ...@@ -55,6 +55,9 @@ MultiplyByInt(const EffectiveDivisor& D, unsigned int k);
Divisor Divisor
Sum(const Divisor& D1, const Divisor& D2); Sum(const Divisor& D1, const Divisor& D2);
Divisor
SimpleSum(const Divisor& D1, const Divisor& D2);
Divisor Divisor
MultiplyByInt(const Divisor& D, unsigned int k); MultiplyByInt(const Divisor& D, unsigned int k);
...@@ -65,7 +68,7 @@ struct RRspace { ...@@ -65,7 +68,7 @@ struct RRspace {
BivPol denom; BivPol denom;
std::vector<BivPol> num_basis; std::vector<BivPol> num_basis;
// Residual divisor added with zeroes: (h) - D_+ + D_- // Residual divisor: (h) - D_+ + D_-
// Useful for implementing the arithmetic in the Jacobian with reduced // Useful for implementing the arithmetic in the Jacobian with reduced
// representations of divisors. // representations of divisors.
EffectiveDivisor Dnum; EffectiveDivisor Dnum;
...@@ -75,6 +78,6 @@ struct RRspace { ...@@ -75,6 +78,6 @@ struct RRspace {
}; };
RRspace RRspace
RiemannRochBasis(const Divisor& D); RiemannRochBasis(const Divisor& D, const EffectiveDivisor& E);
#endif #endif
...@@ -32,17 +32,24 @@ int main() { ...@@ -32,17 +32,24 @@ int main() {
ZZ_p::init(ZZ(p)); ZZ_p::init(ZZ(p));
cin >> eq_curve; cin >> eq_curve;
Curve C(&eq_curve); Curve C(&eq_curve);
// Singular divisor for nodal curves
EffectiveDivisor E(C);
cin >> E;
size_t g; size_t g;
cin >> g; cin >> g;
Divisor D(C); Divisor O(C);
cin >> D; cin >> O;
EffectiveDivisor D1(C), D2(C); EffectiveDivisor D1(C), D2(C);
cin >> D1; cin >> D1;
cin >> D2; cin >> D2;
EffectiveDivisor R = AddReduced(g, D, D1, D2); // The following assertion should be true for a nodal curve
assert(g == C.degree() - E.degree());
EffectiveDivisor R = AddReduced(g, O, E, D1, D2);
cout << R << endl; cout << R << endl;
return 0; return 0;
......
...@@ -66,6 +66,8 @@ class Divisor { ...@@ -66,6 +66,8 @@ class Divisor {
public: public:
Divisor() = delete; Divisor() = delete;
Divisor(const Curve& C) : Dpos(C), Dneg(C) {} Divisor(const Curve& C) : Dpos(C), Dneg(C) {}
Divisor(const EffectiveDivisor& _Dpos)
: Dpos(_Dpos), Dneg(EffectiveDivisor(_Dpos.curve())) {}
Divisor(const EffectiveDivisor& _Dpos, const EffectiveDivisor& _Dneg) Divisor(const EffectiveDivisor& _Dpos, const EffectiveDivisor& _Dneg)
: Dpos(_Dpos), Dneg(_Dneg) { : Dpos(_Dpos), Dneg(_Dneg) {
assert(Dpos.curve() == Dneg.curve()); assert(Dpos.curve() == Dneg.curve());
......
...@@ -21,30 +21,36 @@ Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA ...@@ -21,30 +21,36 @@ Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
#include "divisor_class_group.h" #include "divisor_class_group.h"
EffectiveDivisor EffectiveDivisor
DivisorReduction(const Divisor& D1, const Divisor& D2, std::size_t g) { DivisorReduction(const Divisor& D, const Divisor& O, const EffectiveDivisor& E, std::size_t g) {
assert(D1.degree() == 0); assert(D.degree() == 0);
assert(D2.degree() == 1); assert(O.degree() == 1);
assert(D1.curve() == D2.curve()); assert(D.curve() == O.curve());
assert(D.curve() == E.curve());
RRspace E = RiemannRochBasis(Sum(D1, MultiplyByInt(D2, g)));
assert(E.num_basis.size() > 0); RRspace rr = RiemannRochBasis(Sum(D, MultiplyByInt(O, g)), E);
return PositiveDifference(PrincipalDivisor(D1.curve(), E.num_basis[0]), E.Dnum); assert(rr.num_basis.size() > 0);
return PositiveDifference(PrincipalDivisor(D.curve(), rr.num_basis[0], E), rr.Dnum);
} }
// D degree 1 divisor // D degree 1 divisor
// E divisor representing the nodes
EffectiveDivisor EffectiveDivisor
AddReduced(std::size_t g, AddReduced(std::size_t g,
const Divisor& O, const Divisor& O,
const EffectiveDivisor& E,
const EffectiveDivisor& D1, const EffectiveDivisor& D1,
const EffectiveDivisor& D2) { const EffectiveDivisor& D2) {
assert(O.curve() == D1.curve()); assert( O.curve() == D1.curve());
assert(O.curve() == D2.curve()); assert( O.curve() == D2.curve());
assert(O.degree() == 1); assert( O.curve() == E.curve());
assert( O.degree() == 1);
assert(D1.degree() <= g); assert(D1.degree() <= g);
assert(D2.degree() <= g); assert(D2.degree() <= g);
Divisor gO = MultiplyByInt(O, g); Divisor gO = MultiplyByInt(O, g);
RRspace E = RiemannRochBasis(Divisor(Sum(Sum(D1, D2), gO.get_neg()), gO.get_pos())); RRspace rr = RiemannRochBasis(Divisor(Sum(Sum(D1, D2), gO.get_neg()), gO.get_pos()), E);
assert(E.num_basis.size() > 0); assert(rr.num_basis.size() > 0);
return PositiveDifference(PrincipalDivisor(O.curve(), E.num_basis[0]), E.Dnum); EffectiveDivisor res = PositiveDifference(PrincipalDivisor(O.curve(), rr.num_basis[0], E), rr.Dnum);
assert(res.degree() == g + E.degree());
return res;
} }
...@@ -40,7 +40,8 @@ DivisorReduction(const Divisor& D1, const Divisor& D2, std::size_t g); ...@@ -40,7 +40,8 @@ DivisorReduction(const Divisor& D1, const Divisor& D2, std::size_t g);
// D degree 1 divisor // D degree 1 divisor
EffectiveDivisor EffectiveDivisor
AddReduced(std::size_t g, AddReduced(std::size_t g,
const Divisor& D, const Divisor& O,
const EffectiveDivisor& E,
const EffectiveDivisor& D1, const EffectiveDivisor& D1,
const EffectiveDivisor& D2); const EffectiveDivisor& D2);
......
1009 1009
[ 4 [810 944 972 589 71] [790 206 709 905] [236 969 139] [311 861] [1] ] [ 4 [810 944 972 589 71] [790 206 709 905] [236 969 139] [311 861] [1] ]
< [1] [] >
3 3
{ < [507 1] [832] > < [1] [] > } { < [507 1] [832] > < [1] [] > }
< [374 74 387 1] [78 130 82] > < [374 74 387 1] [78 130 82] >
......
1009 1009
[ 4 [0 0 1 1008] [] [1008] [] [1] ] [ 4 [0 0 1 1008] [] [1008] [] [1] ]
{ < [0 768 952 607 467 1] [0 377 661 970 690] > < [0 0 1] [0 1] > } < [0 1] [] >
{ < [33 702 755 1] [194 780 875] > < [247 1] [62] > }
...@@ -33,14 +33,14 @@ int main() { ...@@ -33,14 +33,14 @@ int main() {
cin >> eq_curve; cin >> eq_curve;
Curve C(&eq_curve); Curve C(&eq_curve);
// // Singular divisor for nodal curves // Singular divisor for nodal curves
// EffectiveDivisor E(C); EffectiveDivisor E(C);
// cin >> E; cin >> E;
Divisor D(C); Divisor D(C);
cin >> D; cin >> D;
RRspace basisRR = RiemannRochBasis(D); RRspace basisRR = RiemannRochBasis(D, E);
cout << "Dimension: " << basisRR.num_basis.size() << endl; cout << "Dimension: " << basisRR.num_basis.size() << endl;
cout << "Denominator: " << endl; cout << "Denominator: " << endl;
......
...@@ -50,7 +50,7 @@ class EllipticCurveTest : public CppUnit::TestFixture { ...@@ -50,7 +50,7 @@ class EllipticCurveTest : public CppUnit::TestFixture {
} }
void testGroupLaw() { void testGroupLaw() {
NTL::ZZ_pX fO, gO, f1, g1, f2, g2; NTL::ZZ_pX fO, gO, fE, gE, f1, g1, f2, g2;
// O = (781 : 643 : 1) // O = (781 : 643 : 1)
std::stringstream("[228 1]") >> fO; std::stringstream("[228 1]") >> fO;
...@@ -59,6 +59,11 @@ class EllipticCurveTest : public CppUnit::TestFixture { ...@@ -59,6 +59,11 @@ class EllipticCurveTest : public CppUnit::TestFixture {
EffectiveDivisor(*C, fO, gO), EffectiveDivisor(*C, fO, gO),
EffectiveDivisor(*C)); EffectiveDivisor(*C));
// E is the zero divisor
std::stringstream("[1]") >> fE;
std::stringstream("[]") >> gE;
EffectiveDivisor E(*C, fE, gE);
// P1 = (381 : 777 : 1) // P1 = (381 : 777 : 1)
std::stringstream("[628 1]") >> f1; std::stringstream("[628 1]") >> f1;
std::stringstream("[777]") >> g1; std::stringstream("[777]") >> g1;
...@@ -70,7 +75,7 @@ class EllipticCurveTest : public CppUnit::TestFixture { ...@@ -70,7 +75,7 @@ class EllipticCurveTest : public CppUnit::TestFixture {
EffectiveDivisor P2(*C, f2, g2); EffectiveDivisor P2(*C, f2, g2);
EffectiveDivisor P3 = EffectiveDivisor P3 =
AddReduced(1, O, P1, P2); AddReduced(1, O, E, P1, P2);
NTL::ZZ_pX fR, gR; NTL::ZZ_pX fR, gR;
...@@ -107,7 +112,7 @@ class SmoothQuarticTest : public CppUnit::TestFixture { ...@@ -107,7 +112,7 @@ class SmoothQuarticTest : public CppUnit::TestFixture {
} }
void testGroupLaw() { void testGroupLaw() {
NTL::ZZ_pX fO, gO, f1, g1, f2, g2; NTL::ZZ_pX fO, gO, fE, gE, f1, g1, f2, g2;
// O = (502 : 832 : 1) // O = (502 : 832 : 1)
std::stringstream("[507 1]") >> fO; std::stringstream("[507 1]") >> fO;
...@@ -116,6 +121,11 @@ class SmoothQuarticTest : public CppUnit::TestFixture { ...@@ -116,6 +121,11 @@ class SmoothQuarticTest : public CppUnit::TestFixture {
EffectiveDivisor(*C, fO, gO), EffectiveDivisor(*C, fO, gO),
EffectiveDivisor(*C)); EffectiveDivisor(*C));
// E is the zero divisor
std::stringstream("[1]") >> fE;
std::stringstream("[]") >> gE;
EffectiveDivisor E(*C, fE, gE);
// D1 // D1
std::stringstream("[374 74 387 1]") >> f1; std::stringstream("[374 74 387 1]") >> f1;
std::stringstream("[78 130 82]") >> g1; std::stringstream("[78 130 82]") >> g1;
...@@ -127,7 +137,7 @@ class SmoothQuarticTest : public CppUnit::TestFixture { ...@@ -127,7 +137,7 @@ class SmoothQuarticTest : public CppUnit::TestFixture {
EffectiveDivisor D2(*C, f2, g2); EffectiveDivisor D2(*C, f2, g2);
EffectiveDivisor D3 = EffectiveDivisor D3 =
AddReduced(3, O, D1, D2); AddReduced(3, O, E, D1, D2);
NTL::ZZ_pX fR, gR; NTL::ZZ_pX fR, gR;
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment