From abf5da67450b7f2b55c001e97d6a7d352eb816ce Mon Sep 17 00:00:00 2001 From: Alexander Demin Date: Tue, 11 Jul 2023 22:44:20 +0300 Subject: [PATCH] merge master --- src/Groebner.jl | 1 - src/f4/basis.jl | 6 +- src/f4/f4.jl | 2 +- src/gb/correctness.jl | 216 ---------------------------- src/gb/groebner.jl | 253 --------------------------------- src/interface.jl | 118 ++++++++-------- src/monomials/packedtuples.jl | 259 ---------------------------------- src/utils/types.jl | 6 +- test/array_normalform.jl | 33 ----- test/groebner_maxpairs.jl | 60 -------- test/isgroebner_stress.jl | 90 ------------ test/monoms/monom_orders.jl | 6 +- test/monoms/packedtuples.jl | 66 +-------- test/runtests.jl | 3 +- 14 files changed, 73 insertions(+), 1046 deletions(-) delete mode 100644 src/gb/correctness.jl delete mode 100644 src/gb/groebner.jl delete mode 100644 test/array_normalform.jl delete mode 100644 test/groebner_maxpairs.jl delete mode 100644 test/isgroebner_stress.jl diff --git a/src/Groebner.jl b/src/Groebner.jl index a130bff1..f815f5a0 100644 --- a/src/Groebner.jl +++ b/src/Groebner.jl @@ -129,7 +129,6 @@ export kbase export Lex, DegLex, DegRevLex, InputOrdering, WeightedOrdering, BlockOrdering, MatrixOrdering -export NotPacked, Packed, best_monom_representation @doc read(joinpath(dirname(@__DIR__), "README.md"), String) Groebner diff --git a/src/f4/basis.jl b/src/f4/basis.jl index 4b53d4ba..65729d54 100644 --- a/src/f4/basis.jl +++ b/src/f4/basis.jl @@ -250,12 +250,12 @@ function update_pairset!( !is_gcd_const(ht.monoms[basis.monoms[i][1]], ht.monoms[new_lead]) lcms[i] = get_lcm(basis.monoms[i][1], new_lead, ht, update_ht) deg = update_ht.hashdata[lcms[i]].deg - ps[newidx] = SPair(i, idx, lcms[i], pr(deg)) + ps[newidx] = SPair(Int32(i), Int32(idx), lcms[i], pr(deg)) else # lcm == 0 will mark redundancy of an S-pair lcms[i] = MonomIdx(0) # ps[newidx] = SPair(i, idx, MonomIdx(0), pr(deg)) - ps[newidx] = SPair{pr}(i, idx, MonomIdx(0), typemax(pr)) + ps[newidx] = SPair{pr}(Int32(i), Int32(idx), MonomIdx(0), typemax(pr)) end end @@ -385,7 +385,7 @@ function is_redundant!( # add new S-pair corresponding to Spoly(i, idx) lcm_new = get_lcm(lead_i, lead_new, ht, ht) psidx = pairset.load + 1 - ps[psidx] = SPair{pt}(i, idx, lcm_new, pt(ht.hashdata[lcm_new].deg)) + ps[psidx] = SPair{pt}(Int32(i), Int32(idx), lcm_new, pt(ht.hashdata[lcm_new].deg)) # mark redundant basis.isredundant[idx] = true diff --git a/src/f4/f4.jl b/src/f4/f4.jl index 617fc288..5d00b008 100644 --- a/src/f4/f4.jl +++ b/src/f4/f4.jl @@ -592,7 +592,7 @@ function select_normal!( end end - @info "Selected $(npairs) pairs" + @log level=-4 "Selected $(npairs) pairs in select_normal!" reinitialize_matrix!(matrix, npairs) diff --git a/src/gb/correctness.jl b/src/gb/correctness.jl deleted file mode 100644 index 563270c2..00000000 --- a/src/gb/correctness.jl +++ /dev/null @@ -1,216 +0,0 @@ - -# If one of the primes in the modular computation -# was unlucky all along -function _monte_carlo_error(msg) - @warn msg - throw(RecoverableException(msg)) -end - -function _not_a_basis_error(basis, msg) - throw(DomainError(basis, msg)) -end - -# Checks if basis is groebner basis with a randomized algorithm -function _check_isgroebner(basis) - !isgroebner(basis, certify=false) && - _not_a_basis_error(basis, "Input does not look like a groebner basis.") -end - -# Checks that all computed bases are of the same shape, -# if that is not the case, throws -function basis_shape_control(gb_coeffs1::AbstractVector, gb_coeffs2::AbstractVector) - if length(gb_coeffs1) != length(gb_coeffs2) - _monte_carlo_error("Unlucky reduction to zero in probabilistic linear algebra.") - end - @inbounds for p in 1:length(gb_coeffs1) - if length(gb_coeffs1[p]) != length(gb_coeffs2[p]) - _monte_carlo_error("Unlucky cancellation of basis coefficients modulo a prime.") - end - end - nothing -end - -# Heuristic check for correctness of the reconstructed basis over rationals -function rational_correctness_bound(modulo::BigInt) - setprecision(2 * Base.GMP.MPZ.sizeinbase(modulo, 2)) do - ceil(BigInt, BigFloat(modulo)^(1 / 10)) - end -end - -function basis_coeff_control(primetracker, coeffaccum) - modulo = primetracker.modulo - bound = rational_correctness_bound(modulo) - for c in coeffaccum.gb_coeffs_qq - for cc in c - if cc > bound - return nothing - end - end - end - _monte_carlo_error("Unlucky cancellation of basis coefficients modulo a prime.") - nothing -end - -# Check that the basis is reconstructed correctly. -# There are 3 levels of checks: -# - heuristic check (discards obviously bad cases), -# - randomized check, -# - certification. -# -# By default, only the first two are active, which gives the correct basis -# with a high probability -function correctness_check!( - coeffaccum, - coeffbuffer, - primetracker, - meta, - ring, - exps, - coeffs, - coeffs_zz, - gens_temp_ff, - gb_ff, - ht, - maxpairs -) - - # first we check coefficients with a heuristic check only - if meta.heuristiccheck - if !heuristic_correctness_check(coeffaccum.gb_coeffs_qq, primetracker.modulo) - @info "Heuristic check failed." - return false - end - @info "Heuristic check passed!" - end - - # then check that a basis is also a basis modulo a prime - if meta.randomizedcheck - if !randomized_correctness_check!( - coeffbuffer, - coeffaccum, - ring, - coeffs_zz, - gens_temp_ff, - gb_ff, - primetracker, - ht, - maxpairs - ) - @info "Randomized check failed." - return false - end - @info "Randomized check passed!" - end - - if meta.guaranteedcheck - return guaranteed_correctness_check!( - ring, - gb_ff.monoms, - coeffaccum.gb_coeffs_qq, - exps, - coeffs, - gens_temp_ff, - ht, - maxpairs - ) - end - - return true -end - -threshold_in_heuristic(sznum, szden, szmod) = 1.30 * (sznum + sznum) >= szmod - -# Checks that -# ln(num) + ln(den) < c ln(modulo) -# for all coefficients num/den -function heuristic_correctness_check(gbcoeffs_qq, modulo) - lnm = Base.GMP.MPZ.sizeinbase(modulo, 2) - @inbounds for i in 1:length(gbcoeffs_qq) - for j in 1:length(gbcoeffs_qq[i]) - n = numerator(gbcoeffs_qq[i][j]) - d = denominator(gbcoeffs_qq[i][j]) - if threshold_in_heuristic( - Base.GMP.MPZ.sizeinbase(n, 2), - Base.GMP.MPZ.sizeinbase(d, 2), - lnm - ) - return false - end - end - end - return true -end - -function randomized_correctness_check!( - coeffbuffer, - coeffaccum, - ring, - coeffs_zz, - gens_temp_ff, - gb_ff, - primetracker, - ht, - maxpairs -) - goodprime = nextgoodprime!(primetracker) - - reduce_modulo!(coeffbuffer, coeffs_zz, gens_temp_ff.coeffs, goodprime) - gens_ff_copy = deepcopy_basis(gens_temp_ff) - cleanup_basis!(ring, gens_ff_copy, goodprime) - - gb_coeffs_zz = scale_denominators(coeffbuffer, coeffaccum.gb_coeffs_qq) - gb_ff_copy = deepcopy_basis(gb_ff) - reduce_modulo!(coeffbuffer, gb_coeffs_zz, gb_ff_copy.coeffs, goodprime) - cleanup_basis!(ring, gb_ff_copy, goodprime) - - # check that initial ideal contains in the computed groebner basis modulo goodprime - normal_form_f4!(ring, gb_ff_copy, ht, gens_ff_copy) - for i in 1:(gens_ff_copy.ndone) - # meaning that something is not reduced - if !iszero_coeffvector(gens_ff_copy.coeffs[i]) - return false - end - end - - # check that the basis is a groebner basis modulo goodprime - if !isgroebner_f4!(ring, gb_ff_copy, ht) - return false - end - - return true -end - -function guaranteed_correctness_check!( - ring, - gbexps, - gb_coeffs_qq, - exps, - coeffs, - gens_tmp_ff, - ht, - maxpairs -) - @info "Setting parameter certify=true in groebner is not recommended." - - gens_qq, _ = initialize_structures(ring, gens_tmp_ff.monoms[1:(gens_tmp_ff.ntotal)], coeffs, ht) - gb_qq, _ = initialize_structures(ring, gbexps, gb_coeffs_qq, ht) - - normalize_basis!(ring, gb_qq) - normalize_basis!(ring, gens_qq) - - gens_qq_copy = deepcopy_basis(gens_qq) - - normal_form_f4!(ring, gb_qq, ht, gens_qq_copy) - for i in 1:(gens_qq_copy.ndone) - # meaning that it is not reduced - if !iszero_coeffvector(gens_qq_copy.coeffs[i]) - return false - end - end - - if !isgroebner_f4!(ring, gb_qq, ht) - return false - end - - return true -end diff --git a/src/gb/groebner.jl b/src/gb/groebner.jl deleted file mode 100644 index f72263c4..00000000 --- a/src/gb/groebner.jl +++ /dev/null @@ -1,253 +0,0 @@ - -function groebner( - polynomials::AbstractVector, - representation::RepresentationStyle, - reduced::Bool, - ordering::AbstractMonomialOrdering, - certify::Bool, - linalg::Symbol, - rng, - maxpairs::Int -) - #= extract ring information, exponents and coefficients - from input polynomials =# - # Copies input, so that polynomials would not be changed itself. - ring, exps, coeffs = convert_to_internal(representation, polynomials, ordering) - - #= check and set algorithm parameters =# - metainfo = set_metaparameters(ring, ordering, certify, linalg, rng, maxpairs) - # now ring stores computation ordering - # metainfo is now a struct to store target ordering - - iszerobasis = remove_zeros_from_input!(ring, exps, coeffs) - iszerobasis && (return convert_to_output(ring, polynomials, exps, coeffs, metainfo)) - - #= change input ordering if needed =# - newring = assure_ordering!(ring, exps, coeffs, metainfo.targetord) - - #= compute the groebner basis =# - bexps, bcoeffs = groebner(newring, exps, coeffs, reduced, metainfo, maxpairs) - - # ring contains ordering of computation, it is the requested ordering - #= convert result back to representation of input =# - convert_to_output(newring, polynomials, bexps, bcoeffs, metainfo) -end - -#------------------------------------------------------------------------------ -# Finite field groebner - -# groebner over integers modulo a prime is simple: -# just initialize and call f4 modulo a prime. - -function groebner( - ring::PolyRing, - exps::Vector{Vector{M}}, - coeffs::Vector{Vector{C}}, - reduced::Bool, - meta::GroebnerMetainfo{Rng}, - maxpairs -) where {M <: Monom, C <: CoeffFF, Rng <: Random.AbstractRNG} - - # select hashtable size - tablesize = select_tablesize(ring, exps) - - # initialize basis and hashtable structures - basis, ht = initialize_structures(ring, exps, coeffs, meta.rng, tablesize) - - f4!(ring, basis, ht, reduced, meta.linalg, meta.rng, maxpairs=maxpairs) - - # extract exponents from hashtable - gbexps = hash_to_exponents(basis, ht) - - gbexps, basis.coeffs -end - -#------------------------------------------------------------------------------ -# Rational numbers groebner - -# groebner over rationals is implemented roughly in the following way: -# -# k = 1 -# while !(correctly reconstructed) -# k = k*2 -# select a batch of small prime numbers p1..pk -# compute a batch of finite field groebner bases gb1..gbk -# reconstruct gb1..gbk to gb_zz modulo prod(p1..pk) with CRT -# reconstruct gb_zz to gb_qq with rational reconstruction -# if the basis gb_qq is correct, then break -# end -# return gb_qq -# - -initial_batchsize() = 1 -initial_gaps() = (1, 1, 1, 1, 1) -batchsize_multiplier() = 2 - -function groebner( - ring::PolyRing, - exps::Vector{Vector{M}}, - coeffs::Vector{Vector{C}}, - reduced::Bool, - meta::GroebnerMetainfo, - maxpairs -) where {M <: Monom, C <: CoeffQQ} - - # we can mutate coeffs and exps here - - # select hashtable size - tablesize = select_tablesize(ring, exps) - @info "Selected tablesize $tablesize" - - # initialize hashtable and finite field basis structs - gens_temp_ff, ht = initialize_structures_ff(ring, exps, coeffs, meta.rng, tablesize) - gens_ff = deepcopy_basis(gens_temp_ff) - - # now hashtable is filled correctly, - # and gens_temp_ff exponents are correct and in correct order. - # gens_temp_ff coefficients are filled with random stuff and - # gens_temp_ff.ch is 0 - - # to store integer and rational coefficients of groebner basis - coeffaccum = CoeffAccum{BigInt, Rational{BigInt}}() - # to store BigInt buffers and reduce overall memory usage - coeffbuffer = CoeffBuffer() - - # scale coefficients of input to integers - coeffs_zz = scale_denominators(coeffbuffer, coeffs) - - # keeps track of used prime numbers - primetracker = PrimeTracker(coeffs_zz) - - # exps, coeffs and coeffs_zz **must be not changed** during whole computation - i = 1 - - # copy basis so that we initial exponents dont get lost - gens_ff = deepcopy_basis(gens_temp_ff) - - prime = nextluckyprime!(primetracker) - @info "$i: selected lucky prime $prime" - - # perform reduction and store result in gens_ff - reduce_modulo!(coeffbuffer, coeffs_zz, gens_ff.coeffs, prime) - - # do some things to ensure generators are correct - cleanup_basis!(ring, gens_ff, prime) - - # compute groebner basis in finite field. - # Need to make sure input invariants in f4! are satisfied, see f4/f4.jl for details - - tracer = Tracer() - pairset = initialize_pairset(powertype(M)) - - f4!(ring, gens_ff, tracer, pairset, ht, reduced, meta.linalg, meta.rng, maxpairs) - - # reconstruct into integers - @info "CRT modulo ($(primetracker.modulo), $(prime))" - - reconstruct_crt!(coeffbuffer, coeffaccum, primetracker, gens_ff.coeffs, prime) - - # reconstruct into rationals - @info "Reconstructing modulo $(primetracker.modulo)" - reconstruct_modulo!(coeffbuffer, coeffaccum, primetracker) - - correct = false - if correctness_check!( - coeffaccum, - coeffbuffer, - primetracker, - meta, - ring, - exps, - coeffs, - coeffs_zz, - gens_temp_ff, - gens_ff, - ht, - maxpairs - ) - @info "Reconstructed successfuly" - correct = true - end - - # initial batch size - batchsize = initial_batchsize() - gaps = initial_gaps() - # multiplier for the size of the batch of prime numbers - multiplier = batchsize_multiplier() - # if first prime was not successfull - while !correct - if i <= length(gaps) - batchsize = gaps[i] - else - batchsize = batchsize * multiplier - end - - for j in 1:batchsize - i += 1 - - # copy basis so that initial exponents dont get lost - gens_ff = deepcopy_basis(gens_temp_ff) - prime = nextluckyprime!(primetracker) - @info "$i: selected lucky prime $prime" - # perform reduction and store result in gens_ff - reduce_modulo!(coeffbuffer, coeffs_zz, gens_ff.coeffs, prime) - # do some things to ensure generators are correct - cleanup_basis!(ring, gens_ff, prime) - # compute groebner basis in finite field - # Need to make sure input invariants in f4! are satisfied, see f4/f4.jl for details - - f4!( - ring, - gens_ff, - tracer, - pairset, - ht, - reduced, - meta.linalg, - meta.rng, - maxpairs - ) - # reconstruct to integers - @info "CRT modulo ($(primetracker.modulo), $(prime))" - - batchsize > 1 && basis_shape_control(coeffaccum.gb_coeffs_qq, gens_ff.coeffs) - - reconstruct_crt!(coeffbuffer, coeffaccum, primetracker, gens_ff.coeffs, prime) - end - - # reconstruct to rationals - @info "Reconstructing modulo $(primetracker.modulo)" - success = reconstruct_modulo!(coeffbuffer, coeffaccum, primetracker) - !success && continue - - if correctness_check!( - coeffaccum, - coeffbuffer, - primetracker, - meta, - ring, - exps, - coeffs, - coeffs_zz, - gens_temp_ff, - gens_ff, - ht, - maxpairs - ) - @info "Success!" - correct = true - end - - # i > 2^10 && basis_coeff_control(primetracker, coeffaccum) - end - - #= - # TODO - if meta.usefglm - fglm_f4!(ring, gens_ff, ht) - end - =# - - gb_exps = hash_to_exponents(gens_ff, ht) - gb_exps, coeffaccum.gb_coeffs_qq -end diff --git a/src/interface.jl b/src/interface.jl index 3173dbbf..9e037f21 100644 --- a/src/interface.jl +++ b/src/interface.jl @@ -5,13 +5,12 @@ Computes a Groebner basis of the ideal generated by array `polynomials`. ## Possible Options -TODO(alex): default values are not clear. The `groebner` routine takes the following options: -- `reduced`: If the returned basis must be reduced and unique *(default)*. +- `reduced`: If the returned basis must be reduced and unique (default is `true`). - `ordering`: Specifies the monomial ordering. Available monomial orderings are: - `InputOrdering()` for inferring the ordering from the given `polynomials` - *(default)*, + (default), - `Lex()` for lexicographic, - `DegLex()` for degree lexicographic, - `DegRevLex()` for degree reverse lexicographic, @@ -21,22 +20,21 @@ The `groebner` routine takes the following options: For details and examples see the corresponding documentation page. - `certify`: Certify the obtained basis. When this option is `false`, the algorithm is randomized, and the result is correct with high probability - *(default)*. -- `linalg`: Linear algebra backend. Available options are: `:exact` for - deterministic sparse linear algebra, `:prob` for probabilistic sparse linear - algebra *(default)*. + (default is `false`). +- `linalg`: Linear algebra backend. Available options are: `:deterministic` for + deterministic sparse linear algebra, `:randomized` for probabilistic sparse linear + algebra (default is `:randomized`). - `monoms`: Monomial representation used during the computations. The algorithm - automatically tries to choose the best monomial representation. Otherwise, + automatically tries to choose the most suitable monomial representation. Otherwise, one can set `monoms` to one of the following: - - `default_monom_representation()`, - - `NotPacked{<:Unsigned}`, e.g., `NotPacked{UInt32}()`, for not packed - representation with `32` bits per exponent, - - `Packed{<:Unsigned}`, e.g., `Packed{UInt8}()`, for packed representation - with `8` bits per exponent. + - `:default` for the automatic choice (which is default), + - `:packed` for packed representation, + - `:sparse` for sparse representation. - `seed`: The seed for randomization. Default value is `42`. - Random number generator is -- `loglevel`: Logging level. Default value is `Logging.Warn`, - so that only warnings are produced. + The random number generator used in Groebner is `Random.Xoshiro`. +- `loglevel`: Logging level, an integer. Higher values mean less output from the program. + Default value is `0`, which means that only warnings are printed. + Set `loglevel` to negative values, e.g., `-1` for debugging. - `maxpairs`: The maximum number of critical pairs used at once in matrix construction. By default, this is not specified. Tweak this option to lower the amount of RAM consumed. @@ -77,6 +75,16 @@ See also `groebner_apply!`. ## Example +```jldoctest +using Groebner, AbstractAlgebra +R, (x, y) = GF(2^31-1)["x", "y"] + +# Learn +graph, gb_1 = groebner_learn([x*y^2 + x, y*x^2 + y]) + +# Apply +flag, gb_2 = groebner_apply!(graph, [2x*y^2 + 3x, 4y*x^2 + 5y]) +``` """ function groebner_learn(polynomials::AbstractVector; options...) @@ -86,8 +94,22 @@ end """ groebner_apply!(graph, polynomials; options...) +Computes a Groebner basis of `polynomials` using the given computation graph. + +See also `groebner_learn`. + ## Example +```jldoctest +using Groebner, AbstractAlgebra +R, (x, y) = GF(2^31-1)["x", "y"] + +# Learn +graph, gb_1 = groebner_learn([x*y^2 + x, y*x^2 + y]) + +# Apply +flag, gb_2 = groebner_apply!(graph, [2x*y^2 + 3x, 4y*x^2 + 5y]) +``` """ function groebner_apply!(graph, polynomials::AbstractVector; options...) @@ -104,7 +126,7 @@ Checks if the given array `polynomials` forms a Groebner basis. The `isgroebner` routine takes the following options: - `ordering`: Specifies the monomial ordering. Available monomial orderings are: - `InputOrdering()` for inferring the ordering from the given `polynomials` - *(default)*, + (default), - `Lex()` for lexicographic, - `DegLex()` for degree lexicographic, - `DegRevLex()` for degree reverse lexicographic, @@ -112,12 +134,12 @@ The `isgroebner` routine takes the following options: - `BlockOrdering(args...)` for block ordering, - `MatrixOrdering(matrix)` for matrix ordering. For details and examples see the corresponding documentation page. -- `certify`: Certify the obtained basis. When this option is `false`, the - algorithm is randomized, and the result is correct with high probability - *(default)*. +- `certify`: Use deterministic algorithm (default is `false`). - `seed`: The seed for randomization. Default value is `42`. -- `loglevel`: Logging level. Default value is `Logging.Warn`, - so that only warnings are produced. + The random number generator used in Groebner is `Random.Xoshiro`. +- `loglevel`: Logging level, an integer. Higher values mean less output from the program. + Default value is `0`, which means that only warnings are printed. + Set `loglevel` to negative values, e.g., `-1` for debugging. ## Example @@ -150,10 +172,11 @@ Computes the normal form of polynomials `tobereduced` w.r.t `basis`. ## Possible Options -The `isgroebner` routine takes the following options: +The `normalform` routine takes the following options: +- `check`: Check if the given array `basis` forms a Groebner basis (default is `true`). - `ordering`: Specifies the monomial ordering. Available monomial orderings are: - `InputOrdering()` for inferring the ordering from the given `polynomials` - *(default)*, + (default), - `Lex()` for lexicographic, - `DegLex()` for degree lexicographic, - `DegRevLex()` for degree reverse lexicographic, @@ -161,41 +184,12 @@ The `isgroebner` routine takes the following options: - `BlockOrdering(args...)` for block ordering, - `MatrixOrdering(matrix)` for matrix ordering. For details and examples see the corresponding documentation page. -- `certify`: Certify the obtained basis. When this option is `false`, the - algorithm is randomized, and the result is correct with high probability - *(default)*. -- `rng`: Random number generator object (must be `<: Random.AbstractRNG`) - used in the computations. - Default RNG is `Random.Xoshiro(42)`. -- `loglevel`: Logging level. Default value is `Logging.Warn`, - so that only warnings are produced. +- `loglevel`: Logging level, an integer. Higher values mean less output from the program. + Default value is `0`, which means that only warnings are printed. + Set `loglevel` to negative values, e.g., `-1` for debugging. ## Example -Using `DynamicPolynomials`: - -```jldoctest -using Groebner, DynamicPolynomials -@polyvar x y; -isgroebner([x*y^2 + x, y*x^2 + y]) -``` - -Using `AbstractAlgebra`: - -```jldoctest -using Groebner, AbstractAlgebra -R, (x, y) = QQ["x", "y"] -isgroebner([x*y^2 + x, y*x^2 + y]) - -The `basis` is assumed to be a Groebner basis. -If `check=true`, we use randomized algorithm to quickly check -if `basis` is indeed a Groebner basis (default). - -Uses the ordering on `basis` by default. -If `ordering` is explicitly specified, it takes precedence. - -# Example - ```jldoctest julia> using Groebner, DynamicPolynomials julia> @polyvar x y; @@ -215,11 +209,15 @@ normalform(basis::AbstractVector, tobereduced; options...) = Computes the basis of polynomial ring modulo the zero-dimensional ideal generated by `basis` as a basis of vector space. -The `basis` is assumed to be a Groebner basis. -If `check=true`, we use randomized algorithm to quickly check -if `basis` is indeed a Groebner basis (default). +## Possible Options + +The `kbase` routine takes the following options: +- `check`: Check if the given array `basis` forms a Groebner basis (default is `true`). +- `loglevel`: Logging level, an integer. Higher values mean less output from the program. + Default value is `0`, which means that only warnings are printed. + Set `loglevel` to negative values, e.g., `-1` for debugging. -# Example +## Example ```jldoctest julia> using Groebner, DynamicPolynomials diff --git a/src/monomials/packedtuples.jl b/src/monomials/packedtuples.jl index f6cf2adb..c30775cd 100644 --- a/src/monomials/packedtuples.jl +++ b/src/monomials/packedtuples.jl @@ -158,39 +158,6 @@ function construct_monom(::Type{PackedTuple3{T, B}}, ev::Vector{U}) where {T, B, a1 |= s << (indent * 8) PackedTuple3{T, B}(a1, a2, a3) end -function make_ev(::Type{PackedPair4{T, B}}, ev::Vector{U}) where {T, B, U} - n = length(ev) - epc = elperchunk(T, B) - @assert n < 4 * epc - if n < 3 * epc - small = make_ev(PackedPair3{T, B}, ev) - return PackedPair4{T, B}(small.a1, small.a2, small.a3, zero(T)) - end - indent = sizeof(T) - degsize(T, B, n) - a1, a2, a3, a4 = zero(T), zero(T), zero(T), zero(T) - s = zero(T) - @inbounds for i in n:-1:1 - _overflow_check(ev[i], B) - d = T(ev[i]) - if div(i - 1, epc) == 3 - a1 = a1 << (sizeof(B) * 8) - a1 = a1 | d - elseif div(i - 1, epc) == 2 - a2 = a2 << (sizeof(B) * 8) - a2 = a2 | d - elseif div(i - 1, epc) == 1 - a3 = a3 << (sizeof(B) * 8) - a3 = a3 | d - else - a4 = a4 << (sizeof(B) * 8) - a4 = a4 | d - end - _overflow_check(s, B) - s += d - end - a1 |= s << (indent * 8) - PackedPair4{T, B}(a1, a2, a3, a4) -end # Creates a zero exponent vector of the given type of length n function construct_const_monom(::Type{PackedTuple1{T, B}}, n::Integer) where {T, B} @@ -202,9 +169,6 @@ end function construct_const_monom(::Type{PackedTuple3{T, B}}, n::Integer) where {T, B} PackedTuple3{T, B}(zero(T), zero(T), zero(T)) end -function make_zero_ev(::Type{PackedPair4{T, B}}, n::Integer) where {T, B} - PackedPair4{T, B}(zero(T), zero(T), zero(T), zero(T)) -end # Hash of exponent vector `x` # Must be linear in x: @@ -238,20 +202,6 @@ function monom_hash(x::PackedTuple3{T, B}, b::Vector{MH}) where {T, B, MH} ) mod(h, MonomHash) end -function Base.hash(x::PackedPair4{T, B}, b::Vector{MH}) where {T, B, MH} - epc = elperchunk(T, B) - h = packeddot(x.a4, b, B, 0) - h = h + packeddot(x.a3, view(b, (epc + 1):(2 * epc)), B, 0) - h = h + packeddot(x.a2, view(b, (2epc + 1):(3 * epc)), B, 0) - h = - h + packeddot( - x.a1, - view(b, (3 * epc + 1):length(b)), - B, - epc - max(epc - 1, length(b) - 3 * epc) - ) - mod(h, MonomHash) -end # Creates a regular vector from an exponent vector `pv` # and writes the answer to `tmp` @@ -281,20 +231,6 @@ function monom_to_dense_vector!(tmp::Vector{I}, pv::PackedTuple3{T, B}) where {I packedunpack!(view(tmp, (2 * epc + 1):length(tmp)), pv.a1, B, indent) tmp end -function make_dense!(tmp::Vector{I}, pv::PackedPair4{T, B}) where {I, T, B} - epc = elperchunk(T, B) - (length(tmp) < 3 * epc) && - return make_dense!(tmp, PackedPair3{T, B}(pv.a1, pv.a2, pv.a3)) - indent = 0 - packedunpack!(tmp, pv.a4, B, indent) - indent = 0 - packedunpack!(view(tmp, (epc + 1):(2 * epc)), pv.a3, B, indent) - indent = 0 - packedunpack!(view(tmp, (2 * epc + 1):(3 * epc)), pv.a2, B, indent) - indent = epc - min(epc - 1, length(tmp) - 3 * epc) - packedunpack!(view(tmp, (3 * epc + 1):length(tmp)), pv.a1, B, indent) - tmp -end #------------------------------------------------------------------------------ # Monomial orderings implementations @@ -334,17 +270,6 @@ function monom_isless( b = construct_monom(ExponentVector{T}, monom_to_dense_vector!(tmp2, eb)) monom_isless(a, b, ord) end -function monom_isless( - ea::PackedPair4{T, B}, - eb::PackedPair4{T, B}, - ord::Ord -) where {T, B, Ord <: AbstractMonomialOrdering} - s = 4 * div(sizeof(T), sizeof(B)) - 1 - tmp1, tmp2 = Vector{T}(undef, s), Vector{T}(undef, s) - a = make_ev(PowerVector{T}, make_dense!(tmp1, ea)) - b = make_ev(PowerVector{T}, make_dense!(tmp2, eb)) - monom_isless(a, b, ord) -end function monom_isless( ea::PackedTuple1{T, B}, @@ -410,34 +335,6 @@ function monom_isless( end end -function monom_isless( - ea::PackedPair4{T, B}, - eb::PackedPair4{T, B}, - ::DegRevLex -) where {T, B} - da, db = totaldeg(ea), totaldeg(eb) - if da < db - return true - end - if da > db - return false - end - - if ea.a1 == eb.a1 - if ea.a2 == eb.a2 - if ea.a3 == eb.a3 - return !(ea.a4 <= eb.a4) - else - return !(ea.a3 <= eb.a3) - end - else - return !(ea.a2 <= eb.a2) - end - else - return !(ea.a1 <= eb.a1) - end -end - #------------------------------------------------------------------------------ # Monomial-Monomial arithmetic. @@ -477,20 +374,6 @@ function monom_lcm!( _monom_overflow_check(ans) ans end -function monom_lcm!( - ec::PackedPair4{T, B}, - ea::PackedPair4{T, B}, - eb::PackedPair4{T, B} -) where {T, B} - x1, si1 = packedmax(ea.a1, eb.a1, B, Val(1)) - x2, si2 = packedmax(ea.a2, eb.a2, B, Val(0)) - x3, si3 = packedmax(ea.a3, eb.a3, B, Val(0)) - x4, si4 = packedmax(ea.a4, eb.a4, B, Val(0)) - x1 = x1 + ((si1 + si2 + si3 + si4) << ((sizeof(T) - sizeof(B)) * 8)) - ans = PackedPair4{T, B}(x1, x2, x3, x4) - _overflow_check(ans) - ans -end function is_gcd_const(ea::PackedTuple1{T, B}, eb::PackedTuple1{T, B}) where {T, B} if !packedorth(ea.a1, eb.a1, B, Val(1)) @@ -519,21 +402,6 @@ function is_gcd_const(ea::PackedTuple3{T, B}, eb::PackedTuple3{T, B}) where {T, end return true end -function is_gcd_const(ea::PackedPair4{T, B}, eb::PackedPair4{T, B}) where {T, B} - if !packedorth(ea.a1, eb.a1, B, Val(1)) - return false - end - if !packedorth(ea.a2, eb.a2, B, Val(0)) - return false - end - if !packedorth(ea.a3, eb.a3, B, Val(0)) - return false - end - if !packedorth(ea.a4, eb.a4, B, Val(0)) - return false - end - return true -end function monom_product!( ec::PackedTuple1{T, B}, @@ -568,19 +436,6 @@ function monom_product!( _monom_overflow_check(ans) ans end -function monom_product!( - ec::PackedPair4{T, B}, - ea::PackedPair4{T, B}, - eb::PackedPair4{T, B} -) where {T, B} - x1 = ea.a1 + eb.a1 - x2 = ea.a2 + eb.a2 - x3 = ea.a3 + eb.a3 - x4 = ea.a4 + eb.a4 - ans = PackedPair4{T, B}(x1, x2, x3, x4) - _overflow_check(ans) - ans -end function monom_division!( ec::PackedTuple1{T, B}, @@ -612,18 +467,6 @@ function monom_division!( ans = PackedTuple3{T, B}(x1, x2, x3) ans end -function monom_division!( - ec::PackedPair4{T, B}, - ea::PackedPair4{T, B}, - eb::PackedPair4{T, B} -) where {T, B} - x1 = ea.a1 - eb.a1 - x2 = ea.a2 - eb.a2 - x3 = ea.a3 - eb.a3 - x4 = ea.a4 - eb.a4 - ans = PackedPair4{T, B}(x1, x2, x3, x4) - ans -end function is_monom_divisible(ea::PackedTuple1{T, B}, eb::PackedTuple1{T, B}) where {T, B} if !packedge(ea.a1, eb.a1, B, Val(1)) @@ -652,21 +495,6 @@ function is_monom_divisible(ea::PackedTuple3{T, B}, eb::PackedTuple3{T, B}) wher end return true end -function is_monom_divisible(ea::PackedPair4{T, B}, eb::PackedPair4{T, B}) where {T, B} - if !packedge(ea.a1, eb.a1, B, Val(1)) - return false - end - if !packedge(ea.a2, eb.a2, B, Val(0)) - return false - end - if !packedge(ea.a3, eb.a3, B, Val(0)) - return false - end - if !packedge(ea.a4, eb.a4, B, Val(0)) - return false - end - return true -end function is_monom_divisible!( ec::PackedTuple1{T, B}, @@ -698,16 +526,6 @@ function is_monom_divisible!( ans && (e = monom_division!(ec, ea, eb)) ans, e end -function is_monom_divisible!( - ec::PackedPair4{T, B}, - ea::PackedPair4{T, B}, - eb::PackedPair4{T, B} -) where {T, B} - ans = is_monom_divisible(ea, eb) - e = ec - ans && (e = monom_division!(ec, ea, eb)) - ans, e -end function is_monom_elementwise_eq(ea::PackedTuple1{T, B}, eb::PackedTuple1{T, B}) where {T, B} ea.a1 == eb.a1 @@ -718,9 +536,6 @@ end function is_monom_elementwise_eq(ea::PackedTuple3{T, B}, eb::PackedTuple3{T, B}) where {T, B} ea.a1 == eb.a1 && ea.a2 == eb.a2 && ea.a3 == eb.a3 end -function is_monom_elementwise_eq(ea::PackedPair4{T, B}, eb::PackedPair4{T, B}) where {T, B} - ea.a1 == eb.a1 && ea.a2 == eb.a2 && ea.a3 == eb.a3 && ea.a4 == eb.a4 -end #------------------------------------------------------------------------------ # Monomial division masks. @@ -850,77 +665,3 @@ function monom_divmask( res end - -function monom_divmask( - e::PackedPair4{T, B}, - DM::Type{Mask}, - ndivvars, - divmap, - ndivbits -) where {T, B, Mask} - epc = elperchunk(T, B) - - if ndivvars < 3 * epc - return monom_divmask( - PackedPair3{T, B}(e.a1, e.a2, e.a3), - DM, - ndivvars, - divmap, - ndivbits - ) - end - - ctr = one(Mask) - res = zero(Mask) - o = one(Mask) - - a4 = e.a4 - for i in 1:epc - ei = mod(a3, B) - a4 = a4 >> (sizeof(B) * 8) - for j in 1:ndivbits - @inbounds if ei >= divmap[ctr] - res |= o << (ctr - 1) - end - ctr += o - end - end - - a3 = e.a3 - for i in (epc + 1):(2 * epc) - ei = mod(a3, B) - a3 = a3 >> (sizeof(B) * 8) - for j in 1:ndivbits - @inbounds if ei >= divmap[ctr] - res |= o << (ctr - 1) - end - ctr += o - end - end - - a2 = e.a2 - for i in (2 * epc + 1):(3 * epc) - ei = mod(a2, B) - a2 = a2 >> (sizeof(B) * 8) - for j in 1:ndivbits - @inbounds if ei >= divmap[ctr] - res |= o << (ctr - 1) - end - ctr += o - end - end - - a1 = e.a1 - for i in (3 * epc + 1):min(4 * epc - 1, ndivvars) - ei = mod(a1, B) - a1 = a1 >> (sizeof(B) * 8) - for j in 1:ndivbits - @inbounds if ei >= divmap[ctr] - res |= o << (ctr - 1) - end - ctr += o - end - end - - res -end diff --git a/src/utils/types.jl b/src/utils/types.jl index 5ee9227d..02950117 100644 --- a/src/utils/types.jl +++ b/src/utils/types.jl @@ -98,10 +98,10 @@ const DivisionMask = UInt32 const ColumnIdx = Int32 # TODO -# f4 may fail in some cases and throw a RecoverableException. -# If we catch a RecoverableException, there is a hope to recover the program. +# f4 may fail in some cases and throw a MonomialDegreeOverflow. +# If we catch a MonomialDegreeOverflow, there is a hope to recover the program. # -# Currently, RecoverableException can be caused by one of the following: +# Currently, MonomialDegreeOverflow can be caused by one of the following: # - overflow in monomial operations (monoms/packedpairs.jl), # - bad choce of a prime during rational computation (gb/groebner.jl), # - fail of randomized sparse linear algebra (f4/matrix.jl), diff --git a/test/array_normalform.jl b/test/array_normalform.jl deleted file mode 100644 index 5e337ad4..00000000 --- a/test/array_normalform.jl +++ /dev/null @@ -1,33 +0,0 @@ - -using AbstractAlgebra - -@testset "normalform of an array" begin - for field in [GF(17), GF(2^31 - 1), ZZ, QQ] - R, (x, y) = PolynomialRing(field, ["x", "y"]) - gb = [x, y] - - @test Groebner.normalform(gb, [R(0)]) == [R(0)] - @test Groebner.normalform(gb, [R(0), R(1), R(0), R(0)]) == [R(0), R(1), R(0), R(0)] - @test Groebner.normalform(gb, [R(0), R(0), R(0), R(0)]) == [R(0), R(0), R(0), R(0)] - @test Groebner.normalform(gb, [R(0), x, R(0), x + 1, y, R(0)]) == [R(0), R(0), R(0), R(1), R(0), R(0)] - - @test Groebner.normalform(gb, [x, y + 1]) == [R(0), R(1)] - @test Groebner.normalform(gb, [y + 1, x]) == [R(1), R(0)] - @test Groebner.normalform( - gb, - [x, 3y, 4x + 5y, 6x * y, 7x + 1, y + 2, x^2 + y^2 + 3] - ) == [R(0), R(0), R(0), R(0), R(1), R(2), R(3)] - - gb = [x] - @test Groebner.normalform( - gb, - [x, 3y, 4x + 5y, 6x * y, 7x + 1, y + 2, x^2 + y^2 + 3] - ) == [R(0), 3y, 5y, R(0), R(1), y + 2, y^2 + 3] - - gb = [x^2 + y^2, y^3 - y] - @test Groebner.normalform( - gb, - [R(1), R(5), R(16), 3x, y, 4x * y + 2, y^2, x * y^2 - 8x + y] - ) == [R(1), R(5), R(16), 3x, y, 4x * y + 2, y^2, x * y^2 - 8x + y] - end -end diff --git a/test/groebner_maxpairs.jl b/test/groebner_maxpairs.jl deleted file mode 100644 index e5764901..00000000 --- a/test/groebner_maxpairs.jl +++ /dev/null @@ -1,60 +0,0 @@ - -using Logging - -@testset "groebner maxpairs" begin - s = Groebner.noonn(5, ordering=:degrevlex, ground=GF(2^31 - 1)) - gb = Groebner.groebner(s) - @test gb == Groebner.groebner(s, maxpairs=100) - @test gb == Groebner.groebner(s, maxpairs=10) - @test gb == Groebner.groebner(s, maxpairs=2) - @test gb == Groebner.groebner(s, maxpairs=1) - @test_throws AssertionError Groebner.groebner(s, maxpairs=0) - - s = Groebner.katsuran(5, ordering=:degrevlex, ground=QQ) - gb = Groebner.groebner(s) - @test gb == Groebner.groebner(s, maxpairs=100) - @test gb == Groebner.groebner(s, maxpairs=10) - @test gb == Groebner.groebner(s, maxpairs=2) - @test gb == Groebner.groebner(s, maxpairs=1) - @test_throws AssertionError Groebner.groebner(s, maxpairs=0) - - s = Groebner.katsuran(6, ordering=:deglex, ground=QQ) - gb = Groebner.groebner(s) - @test gb == Groebner.groebner(s, maxpairs=100) - @test gb == Groebner.groebner(s, maxpairs=10) - @test gb == Groebner.groebner(s, maxpairs=2) - @test gb == Groebner.groebner(s, maxpairs=1) - @test_throws AssertionError Groebner.groebner(s, maxpairs=0) - - s = Groebner.cyclicn(5, ordering=:deglex, ground=QQ) - gb = Groebner.groebner(s) - @test gb == Groebner.groebner(s, maxpairs=100, linalg=:exact) - @test gb == Groebner.groebner(s, maxpairs=10, linalg=:exact) - @test gb == Groebner.groebner(s, maxpairs=2, linalg=:exact) - @test gb == Groebner.groebner(s, maxpairs=2, linalg=:prob) - @test gb == Groebner.groebner(s, maxpairs=1, linalg=:exact) - @test gb == Groebner.groebner(s, maxpairs=1, linalg=:prob) - @test_throws AssertionError Groebner.groebner(s, maxpairs=0) - - s = Groebner.noonn(7, ordering=:degrevlex, ground=GF(2^31 - 1)) - gb = Groebner.groebner(s) - @test gb == Groebner.groebner(s, maxpairs=100) - @test gb == Groebner.groebner(s, maxpairs=10) - @test gb == Groebner.groebner(s, maxpairs=1) - @test gb == Groebner.groebner(s, maxpairs=1) - - R, (x, y, z) = QQ["x", "y", "z"] - for i in 1:100 - s = [x * y + z * (1 // BigInt(2)^i), x * z + y * (1 // BigInt(3)^i), y * z + x * (1 // BigInt(5)^i)] - gblex = Groebner.groebner(s, ordering=Groebner.Lex()) - gbdlex = Groebner.groebner(s, ordering=Groebner.DegLex()) - gbdrevlex = Groebner.groebner(s, ordering=Groebner.DegRevLex()) - for maxpairs in [2^10, 10, 5, 1] - @test gbdlex == - Groebner.groebner(s, ordering=Groebner.DegLex(), maxpairs=maxpairs) - @test gblex == Groebner.groebner(s, ordering=Groebner.Lex(), maxpairs=maxpairs) - @test gbdrevlex == - Groebner.groebner(s, ordering=Groebner.DegRevLex(), maxpairs=maxpairs) - end - end -end diff --git a/test/isgroebner_stress.jl b/test/isgroebner_stress.jl deleted file mode 100644 index 40c65581..00000000 --- a/test/isgroebner_stress.jl +++ /dev/null @@ -1,90 +0,0 @@ - -import Random -using AbstractAlgebra - -function test_params_isgb( - rng, - nvariables, - exps, - nterms, - npolys, - grounds, - coeffssize, - orderings -) - without(x, k) = x[1:end .!= k] - - for n in nvariables - for e in exps - for nt in nterms - for np in npolys - for gr in grounds - for ord in orderings - for csz in coeffssize - set = Groebner.generate_set(n, e, nt, np, csz, rng, gr, ord) - isempty(set) && continue - gb = Groebner.groebner(set) - - if Groebner._isgroebner_reference(gb) - @test Groebner.isgroebner(gb) - end - end - end - end - end - end - end - end -end - -@testset "isgroebner random stress tests" begin - rng = Random.MersenneTwister(5) - - # :degrevlex finite case tests - nvariables = [2, 3, 4] - exps = [1:2, 2:4, 4:5] - nterms = [1:1, 1:2, 4:5] - npolys = [1:1, 1:3, 4:5] - grounds = [GF(1031), GF(2^31 - 1)] - coeffssize = [3, 1000, 2^31 - 1] - orderings = [:degrevlex] - p = prod(map(length, (nvariables, exps, nterms, npolys, grounds, orderings, coeffssize))) - @info "producing $p $(orderings[1]) tests for isgroebner" - test_params_isgb(rng, nvariables, exps, nterms, npolys, grounds, coeffssize, orderings) - - # :lex finite case tests - nvariables = [2, 3] - exps = [1:2, 2:4, 2:3] - nterms = [1:1, 1:2, 2:3] - npolys = [1:1, 1:3, 2:3] - grounds = [GF(1031), GF(2^31 - 1)] - coeffssize = [3, 1000, 2^31 - 1] - orderings = [:lex] - p = prod(map(length, (nvariables, exps, nterms, npolys, grounds, orderings, coeffssize))) - @info "producing $p $(orderings[1]) tests for isgroebner" - test_params_isgb(rng, nvariables, exps, nterms, npolys, grounds, coeffssize, orderings) - - # :deglex finite case tests - nvariables = [2, 3] - exps = [1:2, 2:4, 2:3] - nterms = [1:1, 1:2, 2:3] - npolys = [1:1, 1:3, 2:3] - grounds = [GF(1031), GF(2^31 - 1)] - coeffssize = [3, 1000, 2^31 - 1] - orderings = [:deglex] - p = prod(map(length, (nvariables, exps, nterms, npolys, grounds, orderings, coeffssize))) - @info "producing $p $(orderings[1]) tests for isgroebner" - test_params_isgb(rng, nvariables, exps, nterms, npolys, grounds, coeffssize, orderings) - - # all other test cases - nvariables = [2, 3] - exps = [1:2, 2:4] - nterms = [1:1, 1:2] - npolys = [1:1, 1:3] - grounds = [QQ] - coeffssize = [3, 1000, 2^31 - 1] - orderings = [:deglex, :degrevlex, :lex] - p = prod(map(length, (nvariables, exps, nterms, npolys, grounds, orderings, coeffssize))) - @info "producing $p tests for isgroebner" - test_params_isgb(rng, nvariables, exps, nterms, npolys, grounds, coeffssize, orderings) -end diff --git a/test/monoms/monom_orders.jl b/test/monoms/monom_orders.jl index 8f85a4cf..6be2844c 100644 --- a/test/monoms/monom_orders.jl +++ b/test/monoms/monom_orders.jl @@ -87,9 +87,9 @@ implementations_to_test = [ @test !dl(a, b) @test !drl(a, b) - 30 > Groebner.capacity(EV{T}) && continue - a = make_ev(EV{T}, ones(UInt, 30)) - b = make_ev(EV{T}, ones(UInt, 30)) + 30 > Groebner.max_vars_in_monom(EV{T}) && continue + a = construct_monom(EV{T}, ones(UInt, 30)) + b = construct_monom(EV{T}, ones(UInt, 30)) @test !lex(a, b) @test !dl(a, b) @test !drl(a, b) diff --git a/test/monoms/packedtuples.jl b/test/monoms/packedtuples.jl index d0b51854..897efa0d 100644 --- a/test/monoms/packedtuples.jl +++ b/test/monoms/packedtuples.jl @@ -38,12 +38,12 @@ end EV{T, U} = Groebner.PackedTuple2{T, U} where {T, U} x = [i for i in 1:15] - ev = Groebner.make_ev(EV{UInt64, UInt8}, x) + ev = Groebner.construct_monom(EV{UInt64, UInt8}, x) @test ev.a1 == 0x780f0e0d0c0b0a09 @test ev.a2 == 0x0807060504030201 @test Groebner.totaldeg(ev) == sum(x) tmp = similar(x) - @test Groebner.make_dense!(tmp, ev) == x + @test Groebner.monom_to_dense_vector!(tmp, ev) == x x = [1, 2, 3, 0, 4] ev = Groebner.construct_monom(EV{UInt64, UInt8}, x) @@ -86,13 +86,13 @@ end EV{T, U} = Groebner.PackedTuple3{T, U} where {T, U} x = [div(i, 5) for i in 1:23] - ev = Groebner.make_ev(EV{UInt64, UInt8}, x) + ev = Groebner.construct_monom(EV{UInt64, UInt8}, x) @test ev.a1 == 0x2e04040404030303 @test ev.a2 == 0x0303020202020201 @test ev.a3 == 0x0101010100000000 @test Groebner.totaldeg(ev) == sum(x) tmp = similar(x) - @test Groebner.make_dense!(tmp, ev) == x + @test Groebner.monom_to_dense_vector!(tmp, ev) == x x = [1, 2, 3, 0, 4] ev = Groebner.construct_monom(EV{UInt64, UInt8}, x) @@ -101,64 +101,6 @@ end @test ev.a3 == 0x0000000000000000 @test typeof(Groebner.totaldeg(ev)) === UInt64 @test Groebner.totaldeg(ev) == UInt64(10) - @test Groebner.powertype(ev) === Groebner.MonomHash - tmp = similar(x) - @test Groebner.make_dense!(tmp, ev) == [1, 2, 3, 0, 4] - - @test Groebner.capacity(Groebner.PackedPair2{UInt64, UInt8}) == 15 - @test Groebner.capacity(Groebner.PackedPair2{UInt64, UInt16}) == 7 - - ev = Groebner.make_ev(EV{UInt64, UInt16}, x) - @test (ev.a1, ev.a2) == (0x000a000000000004, 0x0000000300020001) - @test typeof(Groebner.totaldeg(ev)) === UInt64 - @test Groebner.totaldeg(ev) == UInt16(10) - - y = [10, 20, 30, 40, 50] - @test_throws Groebner.RecoverableException Groebner.make_ev(EV{UInt64, UInt8}, y) - - y = [1, 2, 3, 4, 5] - ev = Groebner.make_ev(EV{UInt64, UInt32}, y) - @test ev.a1 == 0x0000000f00000005 - @test ev.a2 == 0x0000000400000003 - @test ev.a3 == 0x0000000200000001 - - for B in (UInt8, UInt16) - for T in (UInt64, UInt32) - for k in 1:10 - maxvars = 3 * div(sizeof(T), sizeof(B)) - 1 - n = rand(1:maxvars) - @assert maxvars > 0 - x = rand(0:div(typemax(B), 4 * n), n) - tmp = similar(x) - ev = Groebner.make_ev(EV{T, B}, x) - @test x == Groebner.make_dense!(tmp, ev) - @test x == tmp - end - end - end -end - -@testset "packed exponent pair-4" begin - EV{T, U} = Groebner.PackedPair4{T, U} where {T, U} - - x = [div(i, 4) for i in 1:31] - ev = Groebner.make_ev(EV{UInt64, UInt8}, x) - @test ev.a1 == 0x7007070707060606 - @test ev.a2 == 0x0605050505040404 - @test ev.a3 == 0x0403030303020202 - @test ev.a4 == 0x0201010101000000 - @test Groebner.totaldeg(ev) == sum(x) - tmp = similar(x) - @test Groebner.make_dense!(tmp, ev) == x - - x = [1, 2, 3, 0, 4] - ev = Groebner.make_ev(EV{UInt64, UInt8}, x) - @test ev.a1 == 0x0a00000400030201 - @test ev.a2 == 0x0000000000000000 - @test ev.a3 == 0x0000000000000000 - @test ev.a4 == 0x0000000000000000 - @test typeof(Groebner.totaldeg(ev)) === UInt64 - @test Groebner.totaldeg(ev) == UInt64(10) @test Groebner.entrytype(ev) === Groebner.MonomHash tmp = similar(x) @test Groebner.monom_to_dense_vector!(tmp, ev) == [1, 2, 3, 0, 4] diff --git a/test/runtests.jl b/test/runtests.jl index f9b9a09c..520a0bb9 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -52,9 +52,8 @@ end @includetests [ "normalform/normalform", "normalform/normalform_stress", - "normalform/handling_zeros" ] - @includetests ["fglm/fglm", "fglm/kbase"] + @includetests ["fglm/kbase"] if try_import(:DynamicPolynomials) @includetests ["input-output/DynamicPolynomials"]