code: def dup_zz_zassenhaus(f, K):
      """Factor primitive square-free polynomials in `Z[x]`."""
      n = dup_degree(f)
  
      if n == 1:
          return [f]
  
      from sympy.ntheory import isprime
  
      fc = f[-1]
      A = dup_max_norm(f, K)
      b = dup_LC(f, K)
      B = int(abs(K.sqrt(K(n + 1)) * 2**n * A * b))
      C = int((n + 1) ** (2 * n) * A ** (2 * n - 1))
      gamma = int(_ceil(2 * _log(C, 2)))
      bound = int(2 * gamma * _log(gamma))
      a = []
      # choose a prime number `p` such that `f` be square free in Z_p
      # if there are many factors in Z_p, choose among a few different `p`
      # the one with fewer factors
      for px in range(3, bound + 1):
          if not isprime(px) or b % px == 0:
              continue
  
          px = K.convert(px)
  
          F = gf_from_int_poly(f, px)
  
          if not gf_sqf_p(F, px, K):
              continue
          fsqfx = gf_factor_sqf(F, px, K)[1]
          a.append((px, fsqfx))
          if len(fsqfx) < 15 or len(a) > 4:
              break
      p, fsqf = min(a, key=lambda x: len(x[1]))
  
      l = int(_ceil(_log(2 * B + 1, p)))
  
      modular = [gf_to_int_poly(ff, p) for ff in fsqf]
  
      g = dup_zz_hensel_lift(p, f, modular, l, K)
  
      sorted_T = range(len(g))
      T = set(sorted_T)
      factors, s = [], 1
      pl = p**l
  
      while 2 * s <= len(T):
          for S in subsets(sorted_T, s):
              # lift the constant coefficient of the product `G` of the factors
              # in the subset `S`; if it is does not divide `fc`, `G` does
              # not divide the input polynomial
  
              if b == 1:
                  q = 1
                  for i in S:
                      q = q * g[i][-1]
                  q = q % pl
                  if not _test_pl(fc, q, pl):
                      continue
              else:
                  G = [b]
                  for i in S:
                      G = dup_mul(G, g[i], K)
                  G = dup_trunc(G, pl, K)
                  G = dup_primitive(G, K)[1]
                  q = G[-1]
                  if q and fc % q != 0:
                      continue
  
              H = [b]
              S = set(S)
              T_S = T - S
  
              if b == 1:
                  G = [b]
                  for i in S:
                      G = dup_mul(G, g[i], K)
                  G = dup_trunc(G, pl, K)
  
              for i in T_S:
                  H = dup_mul(H, g[i], K)
  
              H = dup_trunc(H, pl, K)
  
              G_norm = dup_l1_norm(G, K)
              H_norm = dup_l1_norm(H, K)
  
              if G_norm * H_norm <= B:
                  T = T_S
                  sorted_T = [i for i in sorted_T if i not in S]
  
                  G = dup_primitive(G, K)[1]
                  f = dup_primitive(H, K)[1]
  
                  factors.append(G)
                  b = dup_LC(f, K)
  
                  break
          else:
              s += 1
  
      return factors + [f]
  
imports: from sympy.external.gmpy import GROUND_TYPES
  from sympy.core.random import _randint
  from sympy.polys.galoistools import gf_from_int_poly, gf_to_int_poly, gf_lshift, gf_add_mul, gf_mul, gf_div, gf_rem, gf_gcdex, gf_sqf_p, gf_factor_sqf, gf_factor
  from sympy.polys.densebasic import dup_LC, dmp_LC, dmp_ground_LC, dup_TC, dup_convert, dmp_convert, dup_degree, dmp_degree, dmp_degree_in, dmp_degree_list, dmp_from_dict, dmp_zero_p, dmp_one, dmp_nest, dmp_raise, dup_strip, dmp_ground, dup_inflate, dmp_exclude, dmp_include, dmp_inject, dmp_eject, dup_terms_gcd, dmp_terms_gcd
  from sympy.polys.densearith import dup_neg, dmp_neg, dup_add, dmp_add, dup_sub, dmp_sub, dup_mul, dmp_mul, dup_sqr, dmp_pow, dup_div, dmp_div, dup_quo, dmp_quo, dmp_expand, dmp_add_mul, dup_sub_mul, dmp_sub_mul, dup_lshift, dup_max_norm, dmp_max_norm, dup_l1_norm, dup_mul_ground, dmp_mul_ground, dup_quo_ground, dmp_quo_ground
  from sympy.polys.densetools import dup_clear_denoms, dmp_clear_denoms, dup_trunc, dmp_ground_trunc, dup_content, dup_monic, dmp_ground_monic, dup_primitive, dmp_ground_primitive, dmp_eval_tail, dmp_eval_in, dmp_diff_eval_in, dmp_compose, dup_shift, dup_mirror
  from sympy.polys.euclidtools import dmp_primitive, dup_inner_gcd, dmp_inner_gcd
  from sympy.polys.sqfreetools import dup_sqf_p, dup_sqf_norm, dmp_sqf_norm, dup_sqf_part, dmp_sqf_part
  from sympy.polys.polyutils import _sort_factors
  from sympy.polys.polyconfig import query
  from sympy.polys.polyerrors import ExtraneousFactors, DomainError, CoercionFailed, EvaluationFailed
  from sympy.utilities import subsets
  from math import ceil, log
label_desc: Use image node_id 8 for calling a global function with example usage: dup_zz_zassenhaus(f, K) and returns: unknown, list