source: git/factory/fac_multivar.cc @ 125195

spielwiese
Last change on this file since 125195 was 125195, checked in by Rüdiger Stobbe <stobbe@…>, 28 years ago
Version in work. git-svn-id: file:///usr/local/Singular/svn/trunk@47 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 8.4 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2// $Id: fac_multivar.cc,v 1.1 1996-07-23 09:18:35 stobbe Exp $
3
4/*
5$Log: not supported by cvs2svn $
6Revision 1.0  1996/05/17 10:59:45  stobbe
7Initial revision
8
9*/
10
11#define TIMING
12#undef DEBUGOUTPUT
13
14#include "timing.h"
15
16#include "assert.h"
17#include "cf_defs.h"
18#include "fac_multivar.h"
19#include "fac_univar.h"
20#include "cf_reval.h"
21#include "cf_map.h"
22#include "fac_util.h"
23#include "cf_binom.h"
24#include "cf_iter.h"
25
26
27TIMING_DEFINE_PRINT(fac_content);
28TIMING_DEFINE_PRINT(fac_findeval);
29TIMING_DEFINE_PRINT(fac_distrib);
30
31
32static CFArray
33conv_to_factor_array( const CFFList & L )
34{
35    int n;
36    CFFListIterator I = L;
37    bool negate = false;
38
39    if ( ! I.hasItem() )
40        n = 0;
41    else  if ( I.getItem().factor().inBaseDomain() ) {
42        negate = I.getItem().factor().sign() < 0;
43        I++;
44        n = L.length();
45    }
46    else
47        n = L.length() + 1;
48    CFFListIterator J = I;
49    while ( J.hasItem() ) {
50        n += J.getItem().exp() - 1;
51        J++;
52    }
53    CFArray result( 1, n-1 );
54    int i, j, k;
55    i = 1;
56    while ( I.hasItem() ) {
57        k = I.getItem().exp();
58        for ( j = 1; j <= k; j++ ) {
59            result[i] = I.getItem().factor();
60            i++;
61        }
62        I++;
63    }
64    if ( negate )
65        result[1] = -result[1];
66    return result;
67}
68
69static modpk
70coeffBound ( const CanonicalForm & f, int p )
71{
72    int * degs = degrees( f );
73    int M = 0, i, k = f.level();
74    for ( i = 1; i <= k; i++ )
75        M += degs[i];
76    CanonicalForm b = 2 * maxCoeff( f ) * power( CanonicalForm( 3 ), M );
77    CanonicalForm B = p;
78    k = 1;
79    while ( B < b ) {
80        B *= p;
81        k++;
82    }
83    return modpk( p, k );
84}
85
86static bool
87nonDivisors ( CanonicalForm omega, CanonicalForm delta, const CFArray & F, CFArray & d )
88{
89    DEBOUTLN( cerr, "nondivisors omega = ", omega );
90    DEBOUTLN( cerr, "nondivisors delta = ", delta );
91    DEBOUTLN( cerr, "nondivisors F = ", F );
92    CanonicalForm q, r;
93    int k = F.size();
94    d = CFArray( 0, k );
95    d[0] = delta * omega;
96    for ( int i = 1; i <= k; i++ ) {
97        q = abs(F[i]);
98        for ( int j = i-1; j >= 0; j-- ) {
99            r = d[j];
100            do {
101                r = gcd( r, q );
102                q = q / r;
103            } while ( r != 1 );
104            if ( q == 1 )
105                return false;
106        }
107        d[i] = q;
108    }
109    return true;
110}
111
112static void
113findEvaluation ( const CanonicalForm & U, const CanonicalForm & V, const CanonicalForm & omega, const CFFList & F, Evaluation & A, CanonicalForm & U0, CanonicalForm & delta, CFArray & D, int r )
114{
115    CanonicalForm Vn;
116    CFFListIterator I;
117    int j;
118    bool found = false;
119    CFArray FF = CFArray( 1, F.length() );
120    if ( r > 0 )
121        A.nextpoint();
122    while ( ! found ) {
123        Vn = A( V );
124        if ( Vn != 0 ) {
125            U0 = A( U );
126            DEBOUTLN( cerr, "findev U0 = ", U0 );
127            if ( isSqrFree( U0 ) ) {
128                delta = content( U0 );
129                DEBOUTLN( cerr, "findev content( U0 ) = ", delta );
130                for ( I = F, j = 1; I.hasItem(); I++, j++ )
131                    FF[j] = A( I.getItem().factor() );
132                found = nonDivisors( omega, delta, FF, D );
133            }
134            else {
135                DEBOUTLN( cerr, "findev not sqrfree :", sqrFree( U0 ) );
136            }
137        }
138        if ( ! found )
139            A.nextpoint();
140    }
141}
142
143static CFArray
144ZFactorizeMulti ( const CanonicalForm & arg )
145{
146    DEBOUTLN( cerr, "MULTIFACTOR START ----------------------------------- level = ", arg.level() );
147    CFMap M;
148    CanonicalForm UU, U = compress( arg, M );
149    CanonicalForm delta, omega, V = LC( U, 1 );
150    int t = U.level();
151    CFFList F = factorize( V );
152    CFFListIterator I, J;
153    CFArray G, lcG, D;
154    int i, j, k, m, r, maxdeg, h;
155    REvaluation A( 2, t, IntRandom( 100 ) );
156    CanonicalForm U0;
157    CanonicalForm ft, ut, gt, d;
158    modpk b;
159    bool negate = false;
160
161#ifdef DEBUGOUTPUT
162    cerr << "fac U = " << U << endl;
163    for ( i = 1; i <= level( U ); i++ )
164        cerr << "fac deg(U," << Variable(i) << ") = "
165             << degree( U, Variable(i) ) << endl
166             << "fac lc(U," << Variable(i) << ") = "
167             << LC( U, Variable(i) ) << endl;
168#endif
169    maxdeg = 0;
170    for ( i = 2; i <= t; i++ ) {
171        j = U.degree( Variable( i ) );
172        if ( j > maxdeg ) maxdeg = j;
173    }
174
175    if ( F.getFirst().factor().inCoeffDomain() ) {
176        omega = F.getFirst().factor();
177        F.removeFirst();
178        if ( omega < 0 ) {
179            negate = true;
180            omega = -omega;
181            U = -U;
182        }
183    }
184    else
185        omega = 1;
186
187    bool goodeval = false;
188    r = 0;
189
190//    for ( i = 0; i < 10*t; i++ )
191//      A.nextpoint();
192
193    DEBOUTLN( cerr, "---------------------------------------", ' ' );
194
195    while ( ! goodeval ) {
196        TIMING_START(fac_findeval);
197        findEvaluation( U, V, omega, F, A, U0, delta, D, r );
198        TIMING_END(fac_findeval);
199        DEBOUTLN( cerr, "fac evaluation = ", A );
200        G = conv_to_factor_array( factorize( U0, false ) );
201        DEBOUTLN( cerr, "fac fac(U0) = ", G );
202        b = coeffBound( U, getZFacModulus().getp() );
203#ifdef DEBUGOUTPUT
204        cerr << "p^k(" << U.level() << ") = " << b.getp() << "^" << b.getk() << endl;
205        cerr << "(fac: U  = " << U << endl
206             << "      U0 = " << U0 << endl
207             << "      V  = " << V << endl
208             << "      F  = " << F << endl
209             << "      a  = " << A << endl
210             << "      d  = " << delta << endl
211             << "      b  = " << b.getp() << "^" << b.getk() << endl
212             << "      ub = " << b.getp() << "^" << getZFacModulus().getk() << endl
213             << "      G  = " << G << endl
214             << "      D  = " << D << " )" << endl;
215#endif
216        r = G.size();
217        DEBOUTLN( cerr, "fac SIZE OF UNIFAC = ", r );
218        lcG = CFArray( 1, r );
219        for ( j = 1; j <= r; j ++ )
220            lcG[j] = 1;
221
222        TIMING_START(fac_distrib);
223        goodeval = true;
224        for ( I = F; goodeval && I.hasItem(); I++ ) {
225            ft = A( I.getItem().factor() );
226            m = I.getItem().exp();
227            j = 1;
228            while ( m > 0 && j <= r ) {
229                ut = lc( G[j] ) * delta;
230                while ( m > 0 && divides( ft, ut ) ) {
231                    m--; ut /= ft;
232                    lcG[j] *= I.getItem().factor();
233                }
234                j++;
235            }
236            goodeval = (m == 0);
237        }
238        TIMING_END(fac_distrib);
239        if ( goodeval ) {
240            if ( delta != 1 ) {
241                for ( j = 1; j <= r; j++ ) {
242                    gt = A( lcG[j] );
243                    d = gcd( gt, lc( G[j] ) );
244                    lcG[j] *= lc( G[j] ) / d;
245                    gt /= d;
246                    G[j] *= gt;
247                    delta /= gt;
248                }
249                DEBOUTLN( cerr, "fac delta = ", delta );
250                if ( delta != 1 ) {
251                    for ( j = 1; j <= r; j++ ) {
252                        G[j] *= delta;
253                        lcG[j] *= delta;
254                    }
255                    UU = U * power( delta, r-1 );
256                }
257                else
258                    UU = U;
259            }
260            else
261                UU = U;
262            DEBOUTLN( cerr, "fac good evaluation, lcG = ", lcG );
263            DEBOUTLN( cerr, "fac                    G = ", G );
264            DEBOUTLN( cerr, "fac delta = ", delta );
265            DEBOUTLN( cerr, "fac omega = ", omega );
266            for ( j = 1; j <= r; j++ ) {
267                gt = A( lcG[j] );
268                if ( gt != lc( G[j] ) ) {
269                    gt = lc( G[j] ) / gt;
270                    lcG[j] *= gt;
271                    omega /= gt;
272                }
273            }
274            DEBOUTLN( cerr, "fac good evaluation, lcG = ", lcG );
275            DEBOUTLN( cerr, "fac                    G = ", G );
276            DEBOUTLN( cerr, "fac delta = ", delta );
277            DEBOUTLN( cerr, "fac omega = ", omega );
278            goodeval = Hensel( UU, G, lcG, A, b, Variable(1) );
279        }
280    }
281    for ( i = 1; i <= r; i++ )
282        G[i] /= icontent( G[i] );
283    // negate noch beachten !
284    if ( negate )
285        G[1] = -G[1];
286    DEBOUTLN( cerr, "MULTIFACTOR END   ----------------------------------- level = ", arg.level() );
287    return G;
288}
289
290CFFList
291ZFactorizeMultivariate ( const CanonicalForm & f, bool issqrfree )
292{
293    CFFList G, F, R;
294    CFArray GG;
295    CFFListIterator i, j;
296    CFMap M;
297    CanonicalForm g, cont;
298    Variable v1, vm;
299    int k, m, n;
300
301    v1 = Variable(1);
302    if ( issqrfree )
303        F = CFFactor( f, 1 );
304    else
305        F = sqrFree( f );
306
307    for ( i = F; i.hasItem(); i++ ) {
308        if ( i.getItem().factor().inCoeffDomain() ) {
309            if ( ! i.getItem().factor().isOne() )
310                R.append( CFFactor( i.getItem().factor(), i.getItem().exp() ) );
311        }
312        else {
313            TIMING_START(fac_content);
314            g = compress( i.getItem().factor(), M );
315            // now after compress g contains Variable(1)
316            vm = g.mvar();
317            g = swapvar( g, v1, vm );
318            cont = content( g );
319            g = swapvar( g / cont, v1, vm );
320            cont = swapvar( cont, v1, vm );
321            n = i.getItem().exp();
322            TIMING_END(fac_content);
323            DEBOUTLN( cerr, "now after content ...", ' ' );
324            if ( g.isUnivariate() ) {
325                G = factorize( g, true );
326                for ( j = G; j.hasItem(); j++ )
327                    if ( ! j.getItem().factor().isOne() )
328                        R.append( CFFactor( M( j.getItem().factor() ), n ) );
329            }
330            else {
331                GG = ZFactorizeMulti( g );
332                m = GG.max();
333                for ( k = GG.min(); k <= m; k++ )
334                    if ( ! GG[k].isOne() )
335                        R.append( CFFactor( M( GG[k] ), n ) );
336            }
337            G = factorize( cont, true );
338            for ( j = G; j.hasItem(); j++ )
339                if ( ! j.getItem().factor().isOne() )
340                    R.append( CFFactor( M( j.getItem().factor() ), n ) );
341        }
342    }
343    return R;
344}
Note: See TracBrowser for help on using the repository browser.