+static umod rem_xq(int q, const umod& b)
+{
+ cl_univpoly_modint_ring UPR = b.ring();
+ cl_modint_ring R = UPR->basering();
+
+ int n = degree(b);
+ if ( n > q ) {
+ umod c = UPR->create(q);
+ c.set_coeff(q, R->one());
+ c.finalize();
+ return c;
+ }
+
+ mvec c(n+1, R->zero());
+ int k = q-n;
+ c[n] = R->one();
+ DCOUTVAR(k);
+
+ int ofs = 0;
+ do {
+ cl_MI qk = div(c[n-ofs], coeff(b, n));
+ if ( !zerop(qk) ) {
+ for ( int i=1; i<=n; ++i ) {
+ c[n-i+ofs] = c[n-i] - qk * coeff(b, n-i);
+ }
+ ofs = ofs ? 0 : 1;
+ DCOUTVAR(ofs);
+ DCOUTVAR(c);
+ }
+ } while ( k-- );
+
+ if ( ofs ) {
+ c.pop_back();
+ }
+ else {
+ c.erase(c.begin());
+ }
+ umod res = umod_from_modvec(c);
+ return res;
+}
+
+static void distinct_degree_factor(const umod& a_, umodvec& result)
+{
+ umod a = COPY(a, a_);
+
+ DCOUT(distinct_degree_factor);
+ DCOUTVAR(a);
+
+ cl_univpoly_modint_ring UPR = a.ring();
+ cl_modint_ring R = UPR->basering();
+ int q = cl_I_to_int(R->modulus);
+ int n = degree(a);
+ size_t nhalf = n/2;
+
+
+ size_t i = 1;
+ umod w = UPR->create(1);
+ w.set_coeff(1, R->one());
+ w.finalize();
+ umod x = COPY(x, w);
+
+ umodvec ai;
+
+ while ( i <= nhalf ) {
+ w = expt_pos(w, q);
+ w = rem(w, a);
+
+ ai.push_back(gcd(a, w-x));
+
+ if ( ai.back() != UPR->one() ) {
+ a = div(a, ai.back());
+ w = rem(w, a);
+ }
+
+ ++i;
+ }
+
+ result = ai;
+ DCOUTVAR(result);
+ DCOUT(END distinct_degree_factor);
+}
+
+static void same_degree_factor(const umod& a, umodvec& result)
+{
+ DCOUT(same_degree_factor);
+
+ cl_univpoly_modint_ring UPR = a.ring();
+ cl_modint_ring R = UPR->basering();
+ int deg = degree(a);
+
+ umodvec buf;
+ distinct_degree_factor(a, buf);
+ int degsum = 0;
+
+ for ( size_t i=0; i<buf.size(); ++i ) {
+ if ( buf[i] != UPR->one() ) {
+ degsum += degree(buf[i]);
+ umodvec upv;
+ berlekamp(buf[i], upv);
+ for ( size_t j=0; j<upv.size(); ++j ) {
+ result.push_back(upv[j]);
+ }
+ }
+ }
+
+ if ( degsum < deg ) {
+ result.push_back(a);
+ }
+
+ DCOUTVAR(result);
+ DCOUT(END same_degree_factor);
+}
+
+static void distinct_degree_factor_BSGS(const umod& a, umodvec& result)
+{
+ DCOUT(distinct_degree_factor_BSGS);
+ DCOUTVAR(a);
+
+ cl_univpoly_modint_ring UPR = a.ring();
+ cl_modint_ring R = UPR->basering();
+ int q = cl_I_to_int(R->modulus);
+ int n = degree(a);
+
+ cl_N pm = 0.3;
+ int l = cl_I_to_int(ceiling1(the<cl_F>(expt(n, pm))));
+ DCOUTVAR(l);
+ umodvec h(l+1, UPR->create(-1));
+ umod qk = UPR->create(1);
+ qk.set_coeff(1, R->one());
+ qk.finalize();
+ h[0] = qk;
+ DCOUTVAR(h[0]);
+ for ( int i=1; i<=l; ++i ) {
+ qk = expt_pos(h[i-1], q);
+ h[i] = rem(qk, a);
+ DCOUTVAR(i);
+ DCOUTVAR(h[i]);
+ }
+
+ int m = std::ceil(((double)n)/2/l);
+ DCOUTVAR(m);
+ umodvec H(m, UPR->create(-1));
+ int ql = std::pow(q, l);
+ H[0] = COPY(H[0], h[l]);
+ DCOUTVAR(H[0]);
+ for ( int i=1; i<m; ++i ) {
+ qk = expt_pos(H[i-1], ql);
+ H[i] = rem(qk, a);
+ DCOUTVAR(i);
+ DCOUTVAR(H[i]);
+ }
+
+ umodvec I(m, UPR->create(-1));
+ for ( int i=0; i<m; ++i ) {
+ I[i] = UPR->one();
+ for ( int j=0; j<l; ++j ) {
+ I[i] = I[i] * (H[i] - h[j]);
+ }
+ DCOUTVAR(i);
+ DCOUTVAR(I[i]);
+ I[i] = rem(I[i], a);
+ DCOUTVAR(I[i]);
+ }
+
+ umodvec F(m, UPR->one());
+ umod f = COPY(f, a);
+ for ( int i=0; i<m; ++i ) {
+ DCOUTVAR(i);
+ umod g = gcd(f, I[i]);
+ if ( g == UPR->one() ) continue;
+ F[i] = g;
+ f = div(f, g);
+ DCOUTVAR(F[i]);
+ }
+
+ result.resize(n, UPR->one());
+ if ( f != UPR->one() ) {
+ result[n] = f;
+ }
+ for ( int i=0; i<m; ++i ) {
+ DCOUTVAR(i);
+ umod f = COPY(f, F[i]);
+ for ( int j=l-1; j>=0; --j ) {
+ umod g = gcd(f, H[i]-h[j]);
+ result[l*(i+1)-j-1] = g;
+ f = div(f, g);
+ }
+ }
+
+ DCOUTVAR(result);
+ DCOUT(END distinct_degree_factor_BSGS);
+}
+
+static void cantor_zassenhaus(const umod& a, umodvec& result)
+{
+}
+
git://www.ginac.de / ginac.git / blobdiff