1from __future__ import print_function
2
3from .PyPolyBoRi import *
4from .easy_polynomials import (easy_linear_polynomials as
5    easy_linear_polynomials_func)
6from .statistics import used_vars_set
7from random import Random
8from warnings import warn
9import copy
10import os
11import sys
12
13
14class GeneratorLimitExceeded(Exception):
15    """docstring for GeneratorLimitExceeded"""
16
17    def __init__(self, strat):
18        #super(GeneratorLimitExceeded, self).__init__()
19        self.strat = strat
20
21#used_polynomials=list()
22
23def pkey(p):
24    return (p[0], len(p))
25mat_counter = 0
26
27
28def build_and_print_matrices(v, strat):
29    """old solution using PIL, the currently used implementation is done in C++
30    and plots the same matrices, as being calculated"""
31    treated = BooleSet()
32    v = list(v)
33    rows = 0
34    polys_in_mat = []
35    if len(v) == 0:
36        return
37    while(len(v) > 0):
38        rows = rows + 1
39        p = v[0]
40        v = v[1:]
41        #used_polynomials.append(p)
42        for m in list(p.terms()):
43            m = Monomial(m)
44            if not m in BooleSet(treated):
45                i = strat.select(m)
46                if i >= 0:
47                    p2 = strat[i]
48                    p2 = p2 * (m // p2.lead())
49                    v.append(p2)
50        polys_in_mat.append(p)
51        treated = treated.union(p.set())
52    m2i = dict([(v, k) for (k, v) in enumerate(list(Polynomial(BooleSet(
53        treated)).terms()))])
54    #print polys_in_mat
55    polys_in_mat.sort(key=Polynomial.lead, reverse=True)
56    polys_in_mat = [[m2i[t] for t in p.terms()] for p in polys_in_mat]
57
58    global mat_counter
59    mat_counter = mat_counter + 1
60    from PIL import Image
61    from PIL import ImageDraw
62
63    rows = len(polys_in_mat)
64    cols = len(m2i)
65    #print cols,rows
66    im = Image.new("1", (cols, rows), "white")
67    #im=Image.new("1",(,10000),"white")
68    for i in range(len(polys_in_mat)):
69        p = polys_in_mat[i]
70        for j in p:
71            assert i < rows
72            assert j < cols
73
74            im.putpixel((j, i), 0)
75
76    file_name = matrix_prefix + str(mat_counter) + ".png"
77    if os.path.exists(file_name):
78        os.remove(file_name)
79    im.save(file_name)
80    del im
81
82    print("MATRIX_SIZE:", rows, "x", cols)
83
84
85def multiply_polynomials(l, ring):
86    """
87    >>> r=Ring(1000)
88    >>> x=r.variable
89    >>> multiply_polynomials([x(3), x(2)+x(5)*x(6), x(0), x(0)+1], r)
90    0
91    """
92    l = [Polynomial(p) for p in l]
93
94    def sort_key(p):
95        return p.navigation().value()
96    l = sorted(l, key=sort_key)
97    res = Polynomial(ring.one())
98    for p in l:
99        res = p * res
100    return res
101
102
103def build_and_print_matrices_deg_colored(v, strat):
104    """old PIL solution using a different color for each degree"""
105    if len(v) == 0:
106        return
107
108    treated = BooleSet()
109    v = list(v)
110    rows = 0
111    polys_in_mat = []
112
113    while(len(v) > 0):
114        rows = rows + 1
115        p = v[0]
116        v = v[1:]
117        for m in list(p.terms()):
118            m = Monomial(m)
119            if not m in BooleSet(treated):
120                i = strat.select(m)
121                if i >= 0:
122                    p2 = strat[i]
123                    p2 = p2 * (m // p2.lead())
124                    v.append(p2)
125        polys_in_mat.append(p)
126        treated = treated.union(p.set())
127    m2i = dict([(v, k) for (k, v) in enumerate(BooleSet(treated))])
128    max_deg = max([m.deg() for m in BooleSet(treated)])
129    if max_deg == 0:
130        max_deg = 1
131    i2deg = dict([(m2i[m], m.deg()) for m in BooleSet(treated)])
132    polys_in_mat = [[m2i[t] for t in p.terms()] for p in polys_in_mat]
133    polys_in_mat.sort(key=pkey)
134    global mat_counter
135    mat_counter = mat_counter + 1
136    from PIL import Image
137    from PIL import ImageDraw
138    from PIL import ImageColor
139
140    rows = len(polys_in_mat)
141    cols = len(m2i)
142    im = Image.new("RGB", (cols, rows), "white")
143    for i in range(len(polys_in_mat)):
144        p = polys_in_mat[i]
145        for j in p:
146            assert i < rows
147            assert j < cols
148            im.putpixel((j, i), ImageColor.getrgb("hsl(" + str(270 - (270 *
149                i2deg[j]) / max_deg) + ",100%,50%)"))
150    file_name = matrix_prefix + str(mat_counter) + ".png"
151    if os.path.exists(file_name):
152        os.remove(file_name)
153    im.save(file_name)
154    del im
155
156    print("MATRIX_SIZE:", rows, "x", cols)
157
158
159def high_probability_polynomials_trick(p, strat):
160    lead_deg = p.lead_deg()
161    if lead_deg <= 4:
162        return
163
164    ring = p.ring()
165    factor = multiply_polynomials(easy_linear_factors(p), ring)
166    p = p / factor
167
168    #again, do it twice, it's cheap
169    lead_deg = p.lead_deg()
170    if lead_deg <= 3:
171        return
172
173    if lead_deg > 9:
174        return
175    uv = p.vars_as_monomial()
176
177    candidates = []
178
179    if uv.deg() <= 4:
180        return
181
182    if not uv.deg() <= lead_deg + 1:
183        return
184
185    space = uv.divisors()
186
187    lead = p.lead()
188    for v in lead.variables():
189        variable_selection = lead // v
190        vars_reversed = reversed(list(variable_selection.variables()))
191        #it's just a way to loop over the cartesian product
192        for assignment in variable_selection.divisors():
193            c_p = assignment
194            for v in vars_reversed:
195                if not assignment.reducible_by(v):
196                    c_p = (v + 1) * c_p
197
198            points = (c_p + 1).zeros_in(space)
199            if p.zeros_in(points).empty():
200                candidates.append(c_p * factor)
201        #there many more combinations depending on plugged in values
202    for c in candidates:
203        strat.add_as_you_wish(c)
204
205
206def symmGB_F2_python(G, deg_bound=1000000000000, over_deg_bound=0,
207    use_faugere=False, use_noro=False, opt_lazy=True, opt_red_tail=True,
208    max_growth=2.0, step_factor=1.0, implications=False, prot=False,
209    full_prot=False, selection_size=1000, opt_exchange=True,
210    opt_allow_recursion=False, ll=False,
211    opt_linear_algebra_in_last_block=True, max_generators=None,
212    red_tail_deg_growth=True, matrix_prefix='mat',
213    modified_linear_algebra=True, draw_matrices=False,
214    easy_linear_polynomials=True):
215    if use_noro and use_faugere:
216        raise ValueError('both use_noro and use_faugere specified')
217
218    def add_to_basis(strat, p):
219        if p.is_zero():
220            if prot:
221                print("-")
222        else:
223            if prot:
224                if full_prot:
225                    print(p)
226                print("Result: ", "deg:", p.deg(), "lm: ", p.lead(), "el: ", p
227                    .elength())
228            if easy_linear_polynomials and p.lead_deg() > 2:
229                lin = easy_linear_polynomials_func(p)
230                for q in lin:
231                    strat.add_generator_delayed(q)
232            old_len = len(strat)
233            strat.add_as_you_wish(p)
234            new_len = len(strat)
235            if new_len == 1 + old_len:
236                high_probability_polynomials_trick(p, strat)
237
238
239            if prot:
240                print("#Generators:", len(strat))
241
242    if (isinstance(G, list)):
243        if len(G) == 0:
244            return []
245        G = [Polynomial(g) for g in G]
246        current_ring = G[0].ring()
247        strat = GroebnerStrategy(current_ring)
248        strat.reduction_strategy.opt_red_tail = opt_red_tail
249        strat.opt_lazy = opt_lazy
250        strat.opt_exchange = opt_exchange
251        strat.opt_allow_recursion = opt_allow_recursion
252        strat.enabled_log = prot
253        strat.reduction_strategy.opt_ll = ll
254        strat.opt_modified_linear_algebra = modified_linear_algebra
255        strat.opt_linear_algebra_in_last_block = (
256            opt_linear_algebra_in_last_block)
257        strat.opt_red_by_reduced = False  # True
258        strat.reduction_strategy.opt_red_tail_deg_growth = red_tail_deg_growth
259
260        strat.opt_draw_matrices = draw_matrices
261        strat.matrix_prefix = matrix_prefix
262
263        for g in  G:
264            if not g.is_zero():
265                strat.add_generator_delayed(g)
266    else:
267        strat = G
268
269    if prot:
270        print("added delayed")
271    i = 0
272    try:
273        while strat.npairs() > 0:
274            if max_generators and len(strat) > max_generators:
275                raise GeneratorLimitExceeded(strat)
276            i = i + 1
277            if prot:
278                print("Current Degree:", strat.top_sugar())
279            if (strat.top_sugar() > deg_bound) and (over_deg_bound <= 0):
280                return strat
281            if (strat.top_sugar() > deg_bound):
282                ps = strat.some_spolys_in_next_degree(over_deg_bound)
283                over_deg_bound -= len(ps)
284            else:
285                ps = strat.some_spolys_in_next_degree(selection_size)
286
287            if ps and ps[0].ring().has_degree_order():
288                ps = [strat.reduction_strategy.cheap_reductions(p) for p in ps]
289                ps = [p for p in ps if not p.is_zero()]
290                if len(ps) > 0:
291                    min_deg = min((p.deg() for p in ps))
292                new_ps = []
293                for p in ps:
294                    if p.deg() <= min_deg:
295                        new_ps.append(p)
296                    else:
297                        strat.add_generator_delayed(p)
298                ps = new_ps
299
300            if prot:
301                print("(", strat.npairs(), ")")
302            if prot:
303                print("start reducing")
304                print("Chain Crit. : ", strat.chain_criterions, "VC:", strat.
305                    variable_chain_criterions, "EASYP", strat.
306                    easy_product_criterions, "EXTP", strat.
307                    extended_product_criterions)
308                print(len(ps), "spolys added")
309
310            if use_noro or use_faugere:
311                v = BoolePolynomialVector()
312
313                for p in ps:
314                    if not p.is_zero():
315                        v.append(p)
316                if use_noro:
317                    res = strat.noro_step(v)
318                else:
319                    res = strat.faugere_step_dense(v)
320
321            else:
322                v = BoolePolynomialVector()
323                for p in ps:
324                    p = Polynomial(
325                        mod_mon_set(
326                        BooleSet(p.set()), strat.reduction_strategy.monomials))
327                    if not p.is_zero():
328                        v.append(p)
329                if len(v) > 100:
330                    res = parallel_reduce(v, strat, int(step_factor * 10),
331                        max_growth)
332                else:
333                    if len(v) > 10:
334                        res = parallel_reduce(v, strat, int(step_factor * 30),
335                            max_growth)
336                    else:
337                        res = parallel_reduce(v, strat, int(step_factor * 100
338                            ), max_growth)
339
340            if prot:
341                print("end reducing")
342
343
344            if len(res) > 0 and res[0].ring().has_degree_order():
345                res_min_deg = min([p.deg() for p in res])
346                new_res = []
347                for p in res:
348                    if p.deg() == res_min_deg:
349                        new_res.append(p)
350                    else:
351                        strat.add_generator_delayed(p)
352                res = new_res
353
354            def sort_key(p):
355                return p.lead()
356            res_cp = sorted(res, key=sort_key)
357
358
359            for p in  res_cp:
360                old_len = len(strat)
361                add_to_basis(strat, p)
362                if implications and old_len == len(strat) - 1:
363                    strat.implications(len(strat) - 1)
364                if p.is_one():
365                    if prot:
366                        print("GB is 1")
367                    return strat
368                if prot:
369                    print("(", strat.npairs(), ")")
370
371            strat.clean_top_by_chain_criterion()
372        return strat
373    except KeyboardInterrupt:
374        #return strat
375        raise
376
377
378def GPS(G, vars_start, vars_end):
379    def step(strat, trace, var, val):
380        print("plugin: ", var, val)
381        print("npairs", strat.npairs())
382        strat = GroebnerStrategy(strat)
383        print("npairs", strat.npairs())
384        strat.add_generator_delayed(Polynomial(Monomial(Variable(var, strat.r)
385            ) + val))
386        strat = symmGB_F2_python(strat, prot=True, deg_bound=2,
387            over_deg_bound=10)
388        if var <= vars_start:
389            strat = symmGB_F2_python(strat, prot=True, opt_lazy=False,
390                redTail=False)
391        if strat.containsOne():
392            pass
393        else:
394            if var <= vars_start:
395                #bug: may contain Delayed polynomials
396                print("!!!!!!! SOLUTION", trace)
397                raise Exception
398                #yield trace
399            else:
400                branch(strat, trace + [(var, val)], var - 1)
401
402    def branch(strat, trace, var):
403        while(strat.variableHasValue(var)):
404            #remember to add value to trace
405            var = var - 1
406        step(strat, trace, var, 0)
407        step(strat, trace, var, 1)
408    if G:
409        strat = GroebnerStrategy(G[0].ring())
410        #strat.add_generator(G[0])
411        for g in G[:]:
412            strat.add_generator_delayed(g)
413        branch(strat, [], vars_end - 1)
414
415
416def GPS_with_proof_path(G, proof_path, deg_bound, over_deg_bound):
417    def step(strat, trace, proof_path, pos, val):
418        print(proof_path)
419        print("plugin: ", pos, val, proof_path[pos])
420        print("npairs", strat.npairs())
421        strat = GroebnerStrategy(strat)
422        print("npairs", strat.npairs())
423        print("npairs", strat.npairs())
424        plug_p = Polynomial(proof_path[pos]) + val
425        plug_p_lead = plug_p.lead()
426        if len(plug_p) == 2 and (plug_p + plug_p_lead).deg() == 0:
427            for v in plug_p_lead:
428                strat.add_generator_delayed(v + 1)
429        else:
430            strat.add_generator_delayed(plug_p)
431        print("npairs", strat.npairs())
432        print("pos:", pos)
433        strat = symmGB_F2_python(strat, deg_bound=deg_bound, opt_lazy=False,
434            over_deg_bound=over_deg_bound, prot=True)
435        print("npairs", strat.npairs())
436        pos = pos + 1
437        if pos >= len(proof_path):
438            print("OVER")
439            strat = symmGB_F2_python(strat, prot=True)
440        if strat.containsOne():
441            pass
442        else:
443            if pos >= len(proof_path):
444                print("npairs", strat.npairs())
445                print("minimized:")
446                for p in strat.minimalize_and_tail_reduce():
447                    print(p)
448                #bug: may contain Delayed polynomials
449                print("!!!!!!! SOLUTION", trace)
450                raise Exception
451                #yield trace
452            else:
453                branch(strat, trace + [(pos, val)], proof_path, pos)
454
455    def branch(strat, trace, proof_path, pos):
456
457        step(strat, trace, proof_path, pos, 0)
458        step(strat, trace, proof_path, pos, 1)
459    strat = GroebnerStrategy(G[0].ring())
460    strat.add_generator(Polynomial(G[0]))
461    for g in G[1:]:
462        strat.add_generator_delayed(Polynomial(g))
463    branch(strat, [], proof_path, 0)
464
465
466def GPS_with_suggestions(G, deg_bound, over_deg_bound, opt_lazy=True,
467    opt_red_tail=True, initial_bb=True):
468    def step(strat, trace, var, val):
469        print(trace)
470        plug_p = val + var
471        print("plugin: ", len(trace), plug_p)
472        trace = trace + [plug_p]
473        strat = GroebnerStrategy(strat)
474
475        strat.add_generator_delayed(plug_p)
476        print("npairs", strat.npairs())
477
478        strat = symmGB_F2_python(strat, deg_bound=deg_bound,
479            opt_lazy=opt_lazy, over_deg_bound=over_deg_bound, prot=True)
480
481        #pos=pos+1
482        if not strat.containsOne():
483            branch(strat, trace)
484
485    def branch(strat, trace):
486        print("branching")
487        index = strat.suggestPluginVariable()
488
489        if index < 0:
490            uv = set(used_vars_set(strat))
491            lv = set([next(iter(p.lead())).index() for p in strat if p.
492                lead_deg() == 1])
493            candidates = uv.difference(lv)
494            if len(candidates) > 0:
495                index = next(iter(candidates)).index()
496        if index >= 0:
497            print("chosen index:", index)
498            step(strat, trace, Polynomial(Monomial(Variable(index))), 0)
499            step(strat, trace, Polynomial(Monomial(Variable(index))), 1)
500        else:
501            print("FINAL!!!", index)
502            strat = symmGB_F2_python(strat, prot=True)
503            if not strat.containsOne():
504                print("TRACE", trace)
505                print("SOLUTION")
506                for p in strat.minimalize_and_tail_reduce():
507                    print(p)
508                raise Exception
509
510    def sort_crit(p):
511        #return (p.deg(),p.lead(),p.elength())
512        return (p.lead(), p.deg(), p.elength())
513    if not G:
514        return
515    strat = GroebnerStrategy(G[0].ring())
516    strat.reduction_strategy.opt_red_tail = opt_red_tail  # True
517    strat.opt_exchange = False
518    strat.opt_allow_recursion = False
519    #strat.opt_red_tail_deg_growth=False
520    strat.opt_lazy = opt_lazy
521    #strat.opt_lazy=True
522    first_deg_bound = 1
523    G = [Polynomial(p) for p in G]
524    G.sort(key=sort_crit)
525    higher_deg = {}
526    if initial_bb:
527        for g in G:
528            print(g)
529            index = strat.select(g.lead())
530            if p.deg() == 1:  # (index<0):
531                strat.add_as_you_wish(g)
532            else:
533                first_deg_bound = max(first_deg_bound, g.deg())
534                strat.add_generator_delayed(g)
535            print(g, len(strat))
536    else:
537        for g in G:
538            strat.add_as_you_wish(g)
539    if initial_bb:
540        strat = symmGB_F2_python(strat, deg_bound=max(deg_bound,
541            first_deg_bound), opt_lazy=opt_lazy, over_deg_bound=0, prot=True)
542    strat.opt_lazy = opt_lazy
543    print("INITIALIZED")
544    branch(strat, [])
545
546
547def GPS_with_non_binary_proof_path(G, proof_path, deg_bound, over_deg_bound):
548    def step(strat, trace, proof_path, pos, choice):
549        print(proof_path)
550        print("plugin: ", pos)
551        print("npairs", strat.npairs())
552        strat = GroebnerStrategy(strat)
553        print("npairs", strat.npairs())
554        print("npairs", strat.npairs())
555        for p in proof_path[pos][choice]:
556            print(p)
557            strat.add_generator_delayed(Polynomial(p))
558
559        print("npairs", strat.npairs())
560        print("pos:", pos)
561        strat = symmGB_F2_python(strat, deg_bound=deg_bound,
562            over_deg_bound=over_deg_bound, prot=True)
563        print("npairs", strat.npairs())
564        pos = pos + 1
565        if pos >= len(proof_path):
566            print("OVER")
567            strat = symmGB_F2_python(strat)
568        if strat.containsOne():
569            pass
570        else:
571            if pos >= len(proof_path):
572                print("npairs", strat.npairs())
573                #strat.to_std_out()
574                l = [p for p in strat]
575                strat2 = symmGB_F2_python(l)
576                #strat2.to_std_out()
577                #bug: may contain Delayed polynomials
578                print("!!!!!!! SOLUTION", trace)
579                raise Exception
580                #yield trace
581            else:
582                branch(strat, trace + [(pos, choice)], proof_path, pos)
583                #workaround because of stack depth
584                #step(strat,trace+[(var,val)],var-1, 0)
585                #step(strat,trace+[(var,val)],var-1, 1)
586
587    def branch(strat, trace, proof_path, pos):
588        for i in range(len(proof_path[pos])):
589            step(strat, trace, proof_path, pos, i)
590
591    strat = GroebnerStrategy(G[0].ring())
592    strat.add_generator(G[0])
593    for g in G[1:]:
594        strat.add_generator_delayed(g)
595    branch(strat, [], proof_path, 0)
596
597
598def symmGB_F2_C(G, opt_exchange=True,
599    deg_bound=1000000000000, opt_lazy=False,
600    over_deg_bound=0, opt_red_tail=True,
601    max_growth=2.0, step_factor=1.0,
602    implications=False, prot=False,
603    full_prot=False, selection_size=1000,
604    opt_allow_recursion=False, use_noro=False, use_faugere=False,
605    ll=False, opt_linear_algebra_in_last_block=True,
606    max_generators=None, red_tail_deg_growth=True,
607    modified_linear_algebra=True, matrix_prefix="",
608    draw_matrices=False):
609    #print implications
610    if use_noro:
611        raise NotImplementedError("noro not implemented for symmgb")
612    if (isinstance(G, list)):
613        if len(G) == 0:
614            return []
615
616        G = [Polynomial(g) for g in G]
617        strat = GroebnerStrategy(G[0].ring())
618        strat.reduction_strategy.opt_red_tail = opt_red_tail
619        strat.enabled_log = prot
620        strat.opt_lazy = opt_lazy
621        strat.opt_exchange = opt_exchange
622        strat.reduction_strategy.opt_ll = ll
623        strat.opt_allow_recursion = opt_allow_recursion
624        strat.opt_linear_algebra_in_last_block = (
625            opt_linear_algebra_in_last_block)
626        strat.enabled_log = prot
627        strat.opt_modified_linear_algebra = modified_linear_algebra
628        strat.matrix_prefix = matrix_prefix
629        strat.opt_draw_matrices = draw_matrices
630        strat.reduction_strategy.opt_red_tail_deg_growth = red_tail_deg_growth
631        #strat.add_generator(G[0])
632
633        strat.redByReduced = False  # True
634
635        #if PROT:
636        #    print "added first"
637        for g in G:  # [1:]:
638            if not g.is_zero():
639                strat.add_generator_delayed(g)
640    strat.symmGB_F2()
641    return strat
642
643
644def normal_form(poly, ideal, reduced=True):
645    """ Simple normal form computation of a polynomial  against an ideal.
646    >>> from brial import declare_ring, normal_form
647    >>> r=declare_ring(['x','y'], globals())
648    >>> normal_form(x+y, [y],reduced=True)
649    x
650    >>> normal_form(x+y,[x,y])
651    0
652    """
653    ring = poly.ring()
654    strat = ReductionStrategy(ring)
655    strat.opt_red_tail = reduced
656    ideal = [Polynomial(p) for p in ideal if p != 0]
657    ideal = sorted(ideal, key=Polynomial.lead)
658    last = None
659    for p in ideal:
660        if p.lead() != last:
661            strat.add_generator(p)
662        else:
663            warn("%s will not used for reductions" % p)
664        last = p.lead()
665    return strat.nf(poly)
666
667
668def _test():
669    import doctest
670    doctest.testmod()
671
672if __name__ == "__main__":
673    _test()
674