source: git/factory/fac_berlekamp.cc @ 16f511

spielwiese
Last change on this file since 16f511 was 16f511, checked in by Oleksandr Motsak <motsak@…>, 11 years ago
Fixed the usage of "config.h" (if defined HAVE_CONFIG_H)
  • Property mode set to 100644
File size: 10.4 KB
Line 
1/* emacs edit mode for this file is -*- C++ -*- */
2
3#ifdef HAVE_CONFIG_H
4#include "config.h"
5#endif /* HAVE_CONFIG_H */
6
7#include "cf_assert.h"
8#include "debug.h"
9
10#include "cf_defs.h"
11#include "fac_berlekamp.h"
12#include "ffops.h"
13#include "gfops.h"
14#include "imm.h"
15#include "canonicalform.h"
16#include "cf_iter.h"
17#include "cf_generator.h"
18#include "fac_sqrfree.h"
19
20#ifdef DEBUGOUTPUT
21void QprintFF( int ** Q, int n )
22{
23    for ( int i = 0; i < n; i++ ) {
24        for ( int j = 0; j < n; j++ )
25            std::cerr << Q[i][j] << "  ";
26        std::cerr << std::endl;
27    }
28    std::cerr << std::endl;
29}
30#endif /* DEBUGOUTPUT */
31
32#ifdef DEBUGOUTPUT
33void QprintGF( int ** Q, int n )
34{
35    for ( int i = 0; i < n; i++ ) {
36        for ( int j = 0; j < n; j++ ) {
37            gf_print( std::cerr, Q[i][j] );
38            std::cerr << "  ";
39        }
40        std::cerr << std::endl;
41    }
42    std::cerr << std::endl;
43}
44#endif /* DEBUGOUTPUT */
45
46void QmatFF ( const CanonicalForm & f, int ** Q, int p )
47{
48    int n = degree( f ), nn = (n-1)*p + 1;
49    int i, m, rn;
50    int * a = new int [n];
51    int * r = new int [n];
52    int * q;
53
54    q = Q[0]; *q = r[0] = 1; a[0] = 0; q++;
55
56    for ( i = 1; i < n; i++, q++ )
57        *q = r[i] = a[i] = 0;
58    CFIterator I = f; I++;
59    while ( I.hasTerms() ) {
60        a[I.exp()] = I.coeff().intval();
61        I++;
62    }
63    for ( m = 1; m < nn; m++ ) {
64        rn = r[n-1];
65        for ( i = n-1; i > 0; i-- )
66            r[i] = ff_sub( r[i-1], ff_mul( rn, a[i] ) );
67        r[0] = ff_mul( ff_neg( rn ), a[0] );
68        if ( m % p == 0 ) {
69            q = Q[m/p];
70            for ( i = 0; i < n; i++, q++ )
71                *q = r[i];
72        }
73    }
74    for ( i = 0; i < n; i++ )
75        Q[i][i] = ff_sub( Q[i][i], 1 );
76
77    delete [] a;
78    delete [] r;
79}
80
81void QmatGF ( const CanonicalForm & f, int ** Q, int p )
82{
83    int n = degree( f ), nn = (n-1)*p + 1;
84    int i, m, rn;
85    int * a = new int [n];
86    int * r = new int [n];
87    int * q;
88
89    q = Q[0]; *q = r[0] = gf_one(); a[0] = gf_zero(); q++;
90
91    for ( i = 1; i < n; i++, q++ )
92        *q = r[i] = a[i] = gf_zero();
93    CFIterator I = f; I++;
94    while ( I.hasTerms() ) {
95        a[I.exp()] = imm2int( I.coeff().getval() );
96        I++;
97    }
98    for ( m = 1; m < nn; m++ ) {
99        rn = r[n-1];
100        for ( i = n-1; i > 0; i-- )
101            r[i] = gf_sub( r[i-1], gf_mul( rn, a[i] ) );
102        r[0] = gf_mul( gf_neg( rn ), a[0] );
103        if ( m % p == 0 ) {
104            q = Q[m/p];
105            for ( i = 0; i < n; i++, q++ )
106                *q = r[i];
107        }
108    }
109    for ( i = 0; i < n; i++ )
110        Q[i][i] = gf_sub( Q[i][i], gf_one() );
111
112    delete [] a;
113    delete [] r;
114}
115
116int nullSpaceFF ( int ** Q, int ** b, int n )
117{
118    int * c = new int[n];
119    int r, i, j, k, h, s, d;
120
121    r = 0;
122    for ( s = 0; s < n; s++ ) c[s] = -1;
123    for ( h = 0; h < n; h++ ) {
124        j = 0;
125        while ( j < n && ! ( Q[h][j] != 0 && c[j] < 0 ) ) j++;
126        if ( j < n ) {
127            d = ff_neg( ff_inv( Q[h][j] ) );
128            for ( s = 0; s < n; s++ )
129                Q[s][j] = ff_mul( d, Q[s][j] );
130            for ( i = 0; i < n; i++ ) {
131                if ( i != j ) {
132                    d = Q[h][i];
133                    for ( s = 0; s < n; s++ )
134                        Q[s][i] = ff_add( ff_mul( d, Q[s][j] ), Q[s][i] );
135                }
136            }
137            c[j] = h;
138        }
139        else {
140            b[r] = new int[n];
141            for ( j = 0; j < n; j++ ) {
142                if ( j == h )
143                    b[r][j] = 1;
144                else {
145                    k = 0;
146                    while ( k < n && c[k] != j ) k++;
147                    if ( k < n )
148                        b[r][j] = Q[h][k];
149                    else
150                        b[r][j] = 0;
151                }
152            }
153            r++;
154        }
155    }
156    delete [] c;
157    return r;
158}
159
160int nullSpaceGF ( int ** Q, int ** b, int n )
161{
162    int * c = new int[n];
163    int r, i, j, k, h, s, d;
164
165    r = 0;
166    for ( s = 0; s < n; s++ ) c[s] = -1;
167    for ( h = 0; h < n; h++ ) {
168        j = 0;
169        while ( j < n && ! ( ! gf_iszero( Q[h][j] ) && c[j] < 0 ) ) j++;
170        if ( j < n ) {
171            d = gf_neg( gf_inv( Q[h][j] ) );
172            for ( s = 0; s < n; s++ )
173                Q[s][j] = gf_mul( d, Q[s][j] );
174            for ( i = 0; i < n; i++ ) {
175                if ( i != j ) {
176                    d = Q[h][i];
177                    for ( s = 0; s < n; s++ )
178                        Q[s][i] = gf_add( gf_mul( d, Q[s][j] ), Q[s][i] );
179                }
180            }
181            c[j] = h;
182        }
183        else {
184            b[r] = new int[n];
185            for ( j = 0; j < n; j++ ) {
186                if ( j == h )
187                    b[r][j] = gf_one();
188                else {
189                    k = 0;
190                    while ( k < n && c[k] != j ) k++;
191                    if ( k < n )
192                        b[r][j] = Q[h][k];
193                    else
194                        b[r][j] = gf_zero();
195                }
196            }
197            r++;
198        }
199    }
200    delete [] c;
201    return r;
202}
203
204CanonicalForm cfFromIntVec( int * a, int n, const Variable & x )
205{
206    CanonicalForm result = power( x, n-1 ) * a[n-1];
207    for ( int i = n-2; i >= 0; i-- )
208        if ( a[i] != 0 )
209            result += power( x, i ) * a[i];
210    return result;
211}
212
213CanonicalForm cfFromGFVec( int * a, int n, const Variable & x )
214{
215    CanonicalForm result = power( x, n-1 ) * CanonicalForm( int2imm_gf( a[n-1] ) );
216    for ( int i = n-2; i >= 0; i-- )
217        if ( ! gf_iszero( a[i] ) )
218            result += power( x, i ) * CanonicalForm( int2imm_gf( a[i] ) );
219    return result;
220}
221
222typedef int * intptr;
223
224CFFList BerlekampFactorFF ( const CanonicalForm & f )
225{
226    CFFList F;
227    int p = getCharacteristic();
228    int r, s, len, i, k, n = degree( f );
229    Variable x = f.mvar();
230    CanonicalForm u, g;
231    intptr* Q = new intptr [n];
232    intptr* B = new intptr [n];
233    for ( i = 0; i < n; i++ )
234        Q[i] = new int[n];
235    QmatFF( f, Q, p );
236#ifdef DEBUGOUTPUT
237    DEBOUTLN( cerr, "Q = " );
238    QprintFF( Q, n );
239#endif /* DEBUGOUTPUT */
240    k = nullSpaceFF( Q, B, n );
241#ifdef DEBUGOUTPUT
242    DEBOUTLN( cerr, "Q = " );
243    QprintFF( Q, n );
244#endif /* DEBUGOUTPUT */
245    F.insert( CFFactor( f, 1 ) );
246    r = 1;
247    len = 1;
248    while ( len < k ) {
249        ASSERT( r < k, "fatal fatal" );
250        ListIterator<CFFactor> I = F;
251        while ( I.hasItem() && len < k ) {
252            u = I.getItem().factor();
253            for ( s = 0; s < p && len < k; s++ ) {
254                g = gcd( cfFromIntVec( B[r], n, x ) - s, u );
255                if ( degree( g ) > 0 && g != u ) {
256                    u /= g;
257                    I.append( CFFactor( g, 1 ) );
258                    I.append( CFFactor( u, 1 ) );
259                    I.remove( 1 );
260                    len++;
261                }
262            }
263            I++;
264        }
265        r++;
266    }
267    for ( i = 0; i < n; i++ )
268        delete [] Q[i];
269    for ( i = 0; i < r; i++ )
270        delete [] B[i];
271    delete [] B;
272    delete [] Q;
273    return F;
274}
275
276CFFList BerlekampFactorGF ( const CanonicalForm & f )
277{
278    CFFList F;
279    int r, len, i, k, n = degree( f );
280    Variable x = f.mvar();
281    CanonicalForm u, g;
282    intptr* Q = new intptr [n];
283    intptr* B = new intptr [n];
284    for ( i = 0; i < n; i++ )
285        Q[i] = new int[n];
286    QmatGF( f, Q, gf_q );
287#ifdef DEBUGOUTPUT
288    DEBOUTLN( cerr, "Q = " );
289    QprintGF( Q, n );
290#endif /* DEBUGOUTPUT */
291    k = nullSpaceGF( Q, B, n );
292#ifdef DEBUGOUTPUT
293    DEBOUTLN( cerr, "Q = " );
294    QprintFF( Q, n );
295#endif /* DEBUGOUTPUT */
296    F.insert( CFFactor( f, 1 ) );
297    r = 1;
298    len = 1;
299    GFGenerator s;
300    while ( len < k ) {
301        ASSERT( r < k, "fatal fatal" );
302        ListIterator<CFFactor> I = F;
303        while ( I.hasItem() && len < k ) {
304            u = I.getItem().factor();
305            for ( s.reset(); s.hasItems() && len < k; s++ ) {
306                g = gcd( cfFromGFVec( B[r], n, x ) - s.item(), u );
307                if ( degree( g ) > 0 && g != u ) {
308                    u /= g;
309                    I.append( CFFactor( g, 1 ) );
310                    I.append( CFFactor( u, 1 ) );
311                    I.remove( 1 );
312                    len++;
313                }
314            }
315            I++;
316        }
317        r++;
318    }
319    for ( i = 0; i < n; i++ )
320        delete [] Q[i];
321    for ( i = 0; i < r; i++ )
322        delete [] B[i];
323    delete [] B;
324    delete [] Q;
325    return F;
326}
327// {
328//     CFFList F;
329//     int p = getCharacteristic();
330//     int r, len, k, n = degree( f );
331//     Variable x = f.mvar();
332//     CanonicalForm u, g;
333//     intptr* Q = new intptr [n];
334//     for ( int i = 0; i < n; i++ )
335//         Q[i] = new int[n];
336//     QmatGF( f, Q, p );
337// //  Qprint( Q, n );
338//     k = nullSpaceGF( Q, n );
339// //  Qprint( Q, n );
340//     F.insert( CFFactor( f, 1 ) );
341//     r = 1;
342//     len = 1;
343//     GFIterator s;
344//     while ( len < k ) {
345//         ListIterator<CFFactor> I = F;
346//         while ( I.hasItem() && len < k ) {
347//             u = I.getItem().factor();
348//             for ( s.reset(); s.hasItems() && len < k; s++ ) {
349//                 g = gcd( cfFromGFVec( Q[r], n, x ) - s.item(), u );
350//                 if ( degree( g ) > 0 && g != u ) {
351//                     u /= g;
352//                     I.append( CFFactor( g, 1 ) );
353//                     I.append( CFFactor( u, 1 ) );
354//                     I.remove( 1 );
355//                     len++;
356//                 }
357//             }
358//             I++;
359//         }
360//         r++;
361//     }
362//     for ( i = 0; i < n; i++ )
363//         delete [] Q[i];
364//     return F;
365// }
366
367CFFList FpFactorizeUnivariateB( const CanonicalForm& f, bool issqrfree )
368{
369    CFFList F, G, H;
370    CanonicalForm fac;
371    ListIterator<CFFactor> i, k;
372    int d;
373    bool galoisfield = getGFDegree() > 1;
374
375    if ( LC( f ).isOne() )
376        if ( issqrfree )
377            F.append( CFFactor( f, 1 ) );
378        else
379            F = sqrFreeFp( f );
380    else {
381        H.append( LC( f ) );
382        if ( issqrfree )
383            F.append( CFFactor( f / LC( f ), 1 ) );
384        else
385            F = sqrFreeFp( f / LC( f ) );
386    }
387    for ( i = F; i.hasItem(); ++i ) {
388        d = i.getItem().exp();
389        fac = i.getItem().factor();
390        if ( galoisfield )
391            G = BerlekampFactorGF( fac / LC( fac ) );
392        else
393            G = BerlekampFactorFF( fac / LC( fac ) );
394        for ( k = G; k.hasItem(); ++k ) {
395            fac = k.getItem().factor();
396            H.append( CFFactor( fac / LC( fac ), d ) );
397        }
398    }
399    return H;
400}
Note: See TracBrowser for help on using the repository browser.