Skip to content

Commit

Permalink
Allow more rings to be used with libsingular
Browse files Browse the repository at this point in the history
  • Loading branch information
@user202729
user202729 committed Dec 3, 2024
1 parent 39ebbe4 commit 2e88a31
Showing 1 changed file with 144 additions and 79 deletions.
223 changes: 144 additions & 79 deletions src/sage/libs/singular/ring.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ AUTHORS:
from sage.cpython.string cimport str_to_bytes, bytes_to_str

from sage.libs.gmp.types cimport __mpz_struct
from sage.libs.gmp.mpz cimport mpz_init_set_ui
from sage.libs.gmp.mpz cimport mpz_init_set

from sage.libs.singular.decl cimport ring, currRing
from sage.libs.singular.decl cimport rChangeCurrRing, rComplete, rDelete, idInit
Expand All @@ -31,6 +31,7 @@ from sage.libs.singular.decl cimport n_coeffType
from sage.libs.singular.decl cimport rDefault, GFInfo, ZnmInfo, nInitChar, AlgExtInfo, TransExtInfo


from sage.rings.integer cimport Integer
from sage.rings.integer_ring cimport IntegerRing_class
from sage.rings.integer_ring import ZZ
import sage.rings.abc
Expand Down Expand Up @@ -80,7 +81,7 @@ if bytes_to_str(rSimpleOrdStr(ringorder_ip)) == "rp":

#############################################################################
cdef ring *singular_ring_new(base_ring, n, names, term_order) except NULL:
"""
r"""
Create a new Singular ring over the ``base_ring`` in ``n``
variables with the names ``names`` and the term order
``term_order``.
Expand Down Expand Up @@ -159,17 +160,118 @@ cdef ring *singular_ring_new(base_ring, n, names, term_order) except NULL:
sage: R.<x,y,z> = F[]
sage: from sage.libs.singular.function import singular_function
sage: sing_print = singular_function('print')
sage: sing_print(R)
'polynomial ring, over a field, global ordering\n// coefficients: ZZ/7(a, b)\n// number of vars : 3\n// block 1 : ordering dp\n// : names x y z\n// block 2 : ordering C'
sage: print(sing_print(R))
polynomial ring, over a field, global ordering
// coefficients: ZZ/7(a, b)
// number of vars : 3
// block 1 : ordering dp
// : names x y z
// block 2 : ordering C
::
sage: F = PolynomialRing(QQ, 's,t').fraction_field()
sage: R.<x,y,z> = F[]
sage: from sage.libs.singular.function import singular_function
sage: sing_print = singular_function('print')
sage: sing_print(R)
'polynomial ring, over a field, global ordering\n// coefficients: QQ(s, t)\n// number of vars : 3\n// block 1 : ordering dp\n// : names x y z\n// block 2 : ordering C'
sage: print(sing_print(R))
polynomial ring, over a field, global ordering
// coefficients: QQ(s, t)
// number of vars : 3
// block 1 : ordering dp
// : names x y z
// block 2 : ordering C
Small primes::
sage: R = PolynomialRing(GF(2), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a field, global ordering
// coefficients: ZZ/2
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
sage: R = PolynomialRing(GF(3), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a field, global ordering
// coefficients: ZZ/3
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
sage: R = PolynomialRing(GF(1000000007), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a field, global ordering
// coefficients: ZZ/1000000007
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Large prime (note that the print is wrong, the field in fact doesn't have zero-divisors)::
sage: R = PolynomialRing(GF(2^128+51), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a ring (with zero-divisors), global ordering
// coefficients: ZZ/(bigint(340282366920938463463374607431768211507)^1)
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Finite field with large degree (note that if stack size is too small and the exponent is too large
a stack overflow may happen inside libsingular)::
sage: R = PolynomialRing(GF(2^160), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a field, global ordering
// coefficients: ZZ/2[z160]/(z160^160+z160^159+z160^155+z160^154+z160^153+z160^152+z160^151+z160^149+z160^148+z160^147+z160^146+z160^145+z160^144+z160^143+z160^141+z160^139+z160^137+z160^131+z160^129+z160^128+z160^127+z160^126+z160^123+z160^122+z160^121+z160^117+z160^116+z160^115+z160^113+z160^111+z160^110+z160^108+z160^106+z160^102+z160^100+z160^99+z160^97+z160^96+z160^95+z160^94+z160^93+z160^92+z160^91+z160^87+z160^86+z160^82+z160^80+z160^79+z160^78+z160^74+z160^73+z160^72+z160^71+z160^70+z160^67+z160^66+z160^65+z160^62+z160^59+z160^58+z160^57+z160^55+z160^54+z160^53+z160^52+z160^51+z160^49+z160^47+z160^44+z160^40+z160^35+z160^32+z160^30+z160^28+z160^27+z160^26+z160^24+z160^23+z160^21+z160^20+z160^18+z160^16+z160^11+z160^10+z160^8+z160^7+1)
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Integer modulo small power of 2::
sage: R = PolynomialRing(Zmod(2^32), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a ring (with zero-divisors), global ordering
// coefficients: ZZ/(2^32)
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Integer modulo large power of 2::
sage: R = PolynomialRing(Zmod(2^1000), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a ring (with zero-divisors), global ordering
// coefficients: ZZ/(bigint(2)^1000)
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Integer modulo large power of odd prime::
sage: R = PolynomialRing(Zmod(3^300), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a ring (with zero-divisors), global ordering
// coefficients: ZZ/(bigint(3)^300)
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Integer modulo non-prime::
sage: R = PolynomialRing(Zmod(15^20), ("a", "b"), implementation="singular"); print(sing_print(R))
polynomial ring, over a ring (with zero-divisors), global ordering
// coefficients: ZZ/bigint(332525673007965087890625)
// number of vars : 2
// block 1 : ordering dp
// : names a b
// block 2 : ordering C
Non-prime finite field with large characteristic (not supported, see :issue:`33319`)::
sage: PolynomialRing(GF((2^31+11)^2), ("a", "b"), implementation="singular")
Traceback (most recent call last):
...
TypeError: characteristic must be <= 2147483647.
"""
cdef long cexponent
cdef GFInfo* _param
Expand All @@ -182,7 +284,7 @@ cdef ring *singular_ring_new(base_ring, n, names, term_order) except NULL:
cdef int offset
cdef int nvars
cdef int characteristic
cdef int modbase
cdef Integer ch, modbase
cdef int ringorder_column_pos
cdef int ringorder_column_asc

Expand Down Expand Up @@ -377,39 +479,56 @@ cdef ring *singular_ring_new(base_ring, n, names, term_order) except NULL:
_cf = nInitChar( n_Z, NULL) # integer coefficient ring
_ring = rDefault (_cf, nvars, _names, nblcks, _order, _block0, _block1, _wvhdl)

elif (isinstance(base_ring, FiniteField_generic) and base_ring.is_prime_field()):
if base_ring.characteristic() <= 2147483647:
elif isinstance(base_ring, sage.rings.abc.IntegerModRing):

ch = base_ring.characteristic()
if ch < 2:
raise NotImplementedError(f"polynomials over {base_ring} are not supported in Singular")

isprime = ch.is_prime()

if isprime and ch <= 2147483647:
assert isinstance(base_ring, FiniteField_generic)
characteristic = base_ring.characteristic()

# example for simpler ring creation interface without monomial orderings:
#_ring = rDefault(characteristic, nvars, _names)

_ring = rDefault(characteristic, nvars, _names, nblcks, _order, _block0, _block1, _wvhdl)

else:
raise TypeError("Characteristic p must be <= 2147483647.")
modbase, cexponent = ch.perfect_power()

# example for simpler ring creation interface without monomial orderings:
#_ring = rDefault(characteristic, nvars, _names)
if modbase == 2:
_cf = nInitChar(n_Z2m, <void *>cexponent)

_ring = rDefault(characteristic, nvars, _names, nblcks, _order, _block0, _block1, _wvhdl)
elif modbase.is_prime():
_info.base = <__mpz_struct*>omAlloc(sizeof(__mpz_struct))
mpz_init_set(_info.base, modbase.value)
_info.exp = cexponent
_cf = nInitChar( n_Znm, <void *>&_info )

else:
_info.base = <__mpz_struct*>omAlloc(sizeof(__mpz_struct))
mpz_init_set(_info.base, ch.value)
_info.exp = 1
_cf = nInitChar( n_Zn, <void *>&_info )
_ring = rDefault(_cf, nvars, _names, nblcks, _order, _block0, _block1, _wvhdl)

elif isinstance(base_ring, FiniteField_generic):
if base_ring.characteristic() <= 2147483647:
characteristic = -base_ring.characteristic() # note the negative characteristic
else:
assert not base_ring.is_prime_field() # would have been handled above
if base_ring.characteristic() > 2147483647:
raise TypeError("characteristic must be <= 2147483647.")

# TODO: This is lazy, it should only call Singular stuff not PolynomialRing()
k = PolynomialRing(base_ring.prime_subfield(),
name=base_ring.variable_name(), order='lex', implementation='singular')
minpoly = base_ring.polynomial()(k.gen())

ch = base_ring.characteristic()
F = ch.factor()
assert(len(F)==1)

modbase = F[0][0]
cexponent = F[0][1]

_ext_names = <char**>omAlloc0(sizeof(char*))
_name = str_to_bytes(k._names[0])
_ext_names[0] = omStrDup(_name)
_cfr = rDefault( modbase, 1, _ext_names )
_cfr = rDefault( <int>base_ring.characteristic(), 1, _ext_names )

_cfr.qideal = idInit(1,1)
rComplete(_cfr, 1)
Expand All @@ -422,60 +541,6 @@ cdef ring *singular_ring_new(base_ring, n, names, term_order) except NULL:

_ring = rDefault (_cf, nvars, _names, nblcks, _order, _block0, _block1, _wvhdl)

elif isinstance(base_ring, sage.rings.abc.IntegerModRing):

ch = base_ring.characteristic()
if ch < 2:
raise NotImplementedError(f"polynomials over {base_ring} are not supported in Singular")

isprime = ch.is_prime()

if not isprime and ch.is_power_of(2):
exponent = ch.nbits() -1
cexponent = exponent

if exponent <= 30:
ringtype = n_Z2m
else:
ringtype = n_Znm

if ringtype == n_Znm:
F = ch.factor()

modbase = F[0][0]
cexponent = F[0][1]

_info.base = <__mpz_struct*>omAlloc(sizeof(__mpz_struct))
mpz_init_set_ui(_info.base, modbase)
_info.exp = cexponent
_cf = nInitChar(ringtype, <void *>&_info)
else: # ringtype == n_Z2m
_cf = nInitChar(ringtype, <void *>cexponent)

elif not isprime and ch.is_prime_power() and ch < ZZ(2)**160:
F = ch.factor()
assert(len(F)==1)

modbase = F[0][0]
cexponent = F[0][1]

_info.base = <__mpz_struct*>omAlloc(sizeof(__mpz_struct))
mpz_init_set_ui(_info.base, modbase)
_info.exp = cexponent
_cf = nInitChar( n_Znm, <void *>&_info )

else:
try:
characteristic = ch
except OverflowError:
raise NotImplementedError("Characteristic %d too big." % ch)

_info.base = <__mpz_struct*>omAlloc(sizeof(__mpz_struct))
mpz_init_set_ui(_info.base, characteristic)
_info.exp = 1
_cf = nInitChar( n_Zn, <void *>&_info )
_ring = rDefault(_cf, nvars, _names, nblcks, _order, _block0, _block1, _wvhdl)

else:
raise NotImplementedError(f"polynomials over {base_ring} are not supported in Singular")

Expand Down

0 comments on commit 2e88a31

Please sign in to comment.