source: git/factory/fac_multivar.cc @ b23c90

spielwiese
Last change on this file since b23c90 was c6caf1, checked in by Jens Schmidt <schmidt@…>, 26 years ago
* fac_util.cc (maxCoeff): function removed. Declarations adapted. All references replaced by `maxNorm()'. git-svn-id: file:///usr/local/Singular/svn/trunk@1223 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 8.0 KB
Line 
1/* emacs edit mode for this file is -*- C++ -*- */
2/* $Id: fac_multivar.cc,v 1.8 1998-03-12 14:30:55 schmidt Exp $ */
3
4#include <config.h>
5
6#include "assert.h"
7#include "debug.h"
8#include "timing.h"
9
10#include "cf_defs.h"
11#include "cf_algorithm.h"
12#include "fac_multivar.h"
13#include "fac_univar.h"
14#include "cf_reval.h"
15#include "cf_map.h"
16#include "fac_util.h"
17#include "cf_binom.h"
18#include "cf_iter.h"
19#include "fac_distrib.h"
20
21
22TIMING_DEFINE_PRINT(fac_content);
23TIMING_DEFINE_PRINT(fac_findeval);
24TIMING_DEFINE_PRINT(fac_distrib);
25TIMING_DEFINE_PRINT(fac_hensel);
26
27static CFArray
28conv_to_factor_array( const CFFList & L )
29{
30    int n;
31    CFFListIterator I = L;
32    bool negate = false;
33
34    if ( ! I.hasItem() )
35        n = 0;
36    else  if ( I.getItem().factor().inBaseDomain() ) {
37        negate = I.getItem().factor().sign() < 0;
38        I++;
39        n = L.length();
40    }
41    else
42        n = L.length() + 1;
43    CFFListIterator J = I;
44    while ( J.hasItem() ) {
45        n += J.getItem().exp() - 1;
46        J++;
47    }
48    CFArray result( 1, n-1 );
49    int i, j, k;
50    i = 1;
51    while ( I.hasItem() ) {
52        k = I.getItem().exp();
53        for ( j = 1; j <= k; j++ ) {
54            result[i] = I.getItem().factor();
55            i++;
56        }
57        I++;
58    }
59    if ( negate )
60        result[1] = -result[1];
61    return result;
62}
63
64static modpk
65coeffBound ( const CanonicalForm & f, int p )
66{
67    int * degs = degrees( f );
68    int M = 0, i, k = f.level();
69    for ( i = 1; i <= k; i++ )
70        M += degs[i];
71    CanonicalForm b = 2 * maxNorm( f ) * power( CanonicalForm( 3 ), M );
72    CanonicalForm B = p;
73    k = 1;
74    while ( B < b ) {
75        B *= p;
76        k++;
77    }
78    return modpk( p, k );
79}
80
81// static bool
82// nonDivisors ( CanonicalForm omega, CanonicalForm delta, const CFArray & F, CFArray & d )
83// {
84//     DEBOUTLN( cerr, "nondivisors omega = " << omega );
85//     DEBOUTLN( cerr, "nondivisors delta = " << delta );
86//     DEBOUTLN( cerr, "nondivisors F = " << F );
87//     CanonicalForm q, r;
88//     int k = F.size();
89//     d = CFArray( 0, k );
90//     d[0] = delta * omega;
91//     for ( int i = 1; i <= k; i++ ) {
92//      q = abs(F[i]);
93//      for ( int j = i-1; j >= 0; j-- ) {
94//          r = d[j];
95//          do {
96//              r = gcd( r, q );
97//              q = q / r;
98//          } while ( r != 1 );
99//          if ( q == 1 )
100//              return false;
101//      }
102//      d[i] = q;
103//     }
104//     return true;
105// }
106
107static void
108findEvaluation ( const CanonicalForm & U, const CanonicalForm & V, const CanonicalForm & omega, const CFFList & F, Evaluation & A, CanonicalForm & U0, CanonicalForm & delta, CFArray & D, int r )
109{
110    DEBINCLEVEL( cerr, "findEvaluation" );
111    CanonicalForm Vn;
112    CFFListIterator I;
113    int j;
114    bool found = false;
115    CFArray FF = CFArray( 1, F.length() );
116    if ( r > 0 )
117        A.nextpoint();
118    while ( ! found ) {
119        Vn = A( V );
120        if ( Vn != 0 ) {
121            U0 = A( U );
122            DEBOUTLN( cerr, "U0 = " << U0 );
123            if ( isSqrFree( U0 ) ) {
124                delta = content( U0 );
125                DEBOUTLN( cerr, "content( U0 ) = " << delta );
126                for ( I = F, j = 1; I.hasItem(); I++, j++ )
127                    FF[j] = A( I.getItem().factor() );
128                found = nonDivisors( omega, delta, FF, D );
129            }
130            else {
131                DEBOUTLN( cerr, "not sqrfree : " << sqrFree( U0 ) );
132            }
133        }
134        if ( ! found )
135            A.nextpoint();
136    }
137    DEBDECLEVEL( cerr, "findEvaluation" );
138}
139
140static CFArray
141ZFactorizeMulti ( const CanonicalForm & arg )
142{
143    DEBINCLEVEL( cerr, "ZFactorizeMulti" );
144    CFMap M;
145    CanonicalForm UU, U = compress( arg, M );
146    CanonicalForm delta, omega, V = LC( U, 1 );
147    int t = U.level();
148    CFFList F = factorize( V );
149    CFFListIterator I, J;
150    CFArray G, lcG, D;
151    int i, j, k, m, r, maxdeg, h;
152    REvaluation A( 2, t, IntRandom( 100 ) );
153    CanonicalForm U0;
154    CanonicalForm ft, ut, gt, d;
155    modpk b;
156    bool negate = false;
157
158    DEBOUTLN( cerr, "-----------------------------------------------------" );
159    DEBOUTLN( cerr, "trying to factorize U = " << U );
160    DEBOUTLN( cerr, "U is a polynomial of level = " << arg.level() );
161    DEBOUTLN( cerr, "U will be factorized with respect to variable " << Variable(1) );
162    DEBOUTLN( cerr, "the leading coefficient of U with respect to that variable is " << V );
163    DEBOUTLN( cerr, "which is factorized as " << F );
164
165    maxdeg = 0;
166    for ( i = 2; i <= t; i++ ) {
167        j = U.degree( Variable( i ) );
168        if ( j > maxdeg ) maxdeg = j;
169    }
170
171    if ( F.getFirst().factor().inCoeffDomain() ) {
172        omega = F.getFirst().factor();
173        F.removeFirst();
174        if ( omega < 0 ) {
175            negate = true;
176            omega = -omega;
177            U = -U;
178        }
179    }
180    else
181        omega = 1;
182
183    bool goodeval = false;
184    r = 0;
185
186//    for ( i = 0; i < 10*t; i++ )
187//      A.nextpoint();
188
189    while ( ! goodeval ) {
190        TIMING_START(fac_findeval);
191        findEvaluation( U, V, omega, F, A, U0, delta, D, r );
192        TIMING_END(fac_findeval);
193        DEBOUTLN( cerr, "the evaluation point to reduce to an univariate problem is " << A );
194        DEBOUTLN( cerr, "corresponding delta = " << delta );
195        DEBOUTLN( cerr, "              omega = " << omega );
196        DEBOUTLN( cerr, "              D     = " << D );
197        DEBOUTLN( cerr, "now factorize the univariate polynomial " << U0 );
198        G = conv_to_factor_array( factorize( U0, false ) );
199        DEBOUTLN( cerr, "which factorizes into " << G );
200        b = coeffBound( U, getZFacModulus().getp() );
201        if ( getZFacModulus().getpk() > b.getpk() )
202            b = getZFacModulus();
203        DEBOUTLN( cerr, "the coefficient bound of the factors of U is " << b.getpk() );
204
205        r = G.size();
206        lcG = CFArray( 1, r );
207        UU = U;
208        DEBOUTLN( cerr, "now trying to distribute the leading coefficients ..." );
209        TIMING_START(fac_distrib);
210        goodeval = distributeLeadingCoeffs( UU, G, lcG, F, D, delta, omega, A, r );
211        TIMING_END(fac_distrib);
212#ifdef DEBUGOUTPUT
213        if ( goodeval ) {
214            DEBOUTLN( cerr, "the univariate factors after distribution are " << G );
215            DEBOUTLN( cerr, "the distributed leading coeffs are " << lcG );
216            DEBOUTLN( cerr, "U may have changed and is now " << UU );
217            DEBOUTLN( cerr, "which has leading coefficient " << LC( UU, Variable(1) ) );
218
219            if ( LC( UU, Variable(1) ) != prod( lcG ) || A(UU) != prod( G ) ) {
220                DEBOUTLN( cerr, "!!! distribution was not correct !!!" );
221                DEBOUTLN( cerr, "product of leading coeffs is " << prod( lcG ) );
222                DEBOUTLN( cerr, "product of univariate factors is " << prod( G ) );
223                DEBOUTLN( cerr, "the new U is evaluated as " << A(UU) );
224            }
225            else
226                DEBOUTLN( cerr, "leading coeffs correct" );
227        }
228        else {
229            DEBOUTLN( cerr, "we have found a bad evaluation point" );
230        }
231#endif
232        if ( goodeval ) {
233            TIMING_START(fac_hensel);
234            goodeval = Hensel( UU, G, lcG, A, b, Variable(1) );
235            TIMING_END(fac_hensel);
236        }
237    }
238    for ( i = 1; i <= r; i++ ) {
239        G[i] /= icontent( G[i] );
240        G[i] = M(G[i]);
241    }
242    // negate noch beachten !
243    if ( negate )
244        G[1] = -G[1];
245    DEBDECLEVEL( cerr, "ZFactorMulti" );
246    return G;
247}
248
249CFFList
250ZFactorizeMultivariate ( const CanonicalForm & f, bool issqrfree )
251{
252    CFFList G, F, R;
253    CFArray GG;
254    CFFListIterator i, j;
255    CFMap M;
256    CanonicalForm g, cont;
257    Variable v1, vm;
258    int k, m, n;
259
260    v1 = Variable(1);
261    if ( issqrfree )
262        F = CFFactor( f, 1 );
263    else
264        F = sqrFree( f );
265
266    for ( i = F; i.hasItem(); i++ ) {
267        if ( i.getItem().factor().inCoeffDomain() ) {
268            if ( ! i.getItem().factor().isOne() )
269                R.append( CFFactor( i.getItem().factor(), i.getItem().exp() ) );
270        }
271        else {
272            TIMING_START(fac_content);
273            g = compress( i.getItem().factor(), M );
274            // now after compress g contains Variable(1)
275            vm = g.mvar();
276            g = swapvar( g, v1, vm );
277            cont = content( g );
278            g = swapvar( g / cont, v1, vm );
279            cont = swapvar( cont, v1, vm );
280            n = i.getItem().exp();
281            TIMING_END(fac_content);
282            DEBOUTLN( cerr, "now after content ..." );
283            if ( g.isUnivariate() ) {
284                G = factorize( g, true );
285                for ( j = G; j.hasItem(); j++ )
286                    if ( ! j.getItem().factor().isOne() )
287                        R.append( CFFactor( M( j.getItem().factor() ), n ) );
288            }
289            else {
290                GG = ZFactorizeMulti( g );
291                m = GG.max();
292                for ( k = GG.min(); k <= m; k++ )
293                    if ( ! GG[k].isOne() )
294                        R.append( CFFactor( M( GG[k] ), n ) );
295            }
296            G = factorize( cont, true );
297            for ( j = G; j.hasItem(); j++ )
298                if ( ! j.getItem().factor().isOne() )
299                    R.append( CFFactor( M( j.getItem().factor() ), n ) );
300        }
301    }
302    return R;
303}
Note: See TracBrowser for help on using the repository browser.