Context Navigation

source: git/factory/facAbsBiFact.cc @ 2234726

fieker-DuValspielwiese

Last change on this file since 2234726 was 0341f1, checked in by Hans Schoenemann <hannes@…>, 3 years ago
fix: absFactorize with Flint, w/o NTL
Property mode set to `100644`
File size: 21.7 KB

Line
1	/*****************************************************************************\
2	* Computer Algebra System SINGULAR
3	\*****************************************************************************/
4	/** @file facAbsBiFact.cc
5	*
6	* @author Martin Lee
7	*
8	**/
9	/*****************************************************************************/
10
11
12	#include "config.h"
13
14
15	#include "timing.h"
16	#include "debug.h"
17
18	#include "facAbsBiFact.h"
19	#include "facBivar.h"
20	#include "facFqBivar.h"
21	#include "cf_reval.h"
22	#include "cf_primes.h"
23	#include "cf_algorithm.h"
24	#ifdef HAVE_FLINT
25	#include "FLINTconvert.h"
26	#include <flint/fmpz_poly_factor.h>
27	#endif
28	#ifdef HAVE_NTL
29	#include "NTLconvert.h"
30	#include <NTL/LLL.h>
31	#endif
32
33	#if defined(HAVE_NTL) \|\| defined(HAVE_FLINT)
34	TIMING_DEFINE_PRINT(fac_Qa_factorize)
35	TIMING_DEFINE_PRINT(fac_evalpoint)
36
37	CFAFList uniAbsFactorize (const CanonicalForm& F, bool full)
38	{
39	CFAFList result;
40	if (degree (F) == 1)
41	{
42	bool isRat= isOn (SW_RATIONAL);
43	On (SW_RATIONAL);
44	result= CFAFList (CFAFactor (F/Lc(F), 1, 1));
45	result.insert (CFAFactor (Lc (F), 1, 1));
46	if (!isRat)
47	Off (SW_RATIONAL);
48	return result;
49	}
50	CanonicalForm LcF= 1;
51	Variable alpha;
52	CFFList QaFactors;
53	CFFListIterator iter;
54	alpha= rootOf (F);
55	QaFactors= factorize (F, alpha);
56	iter= QaFactors;
57	if (iter.getItem().factor().inCoeffDomain())
58	{
59	LcF = iter.getItem().factor();
60	iter++;
61	}
62	for (;iter.hasItem(); iter++)
63	{
64	if (full)
65	result.append (CFAFactor (iter.getItem().factor(), getMipo (alpha),
66	iter.getItem().exp()));
67	if (!full && degree (iter.getItem().factor()) == 1)
68	{
69	result.append (CFAFactor (iter.getItem().factor(), getMipo (alpha),
70	iter.getItem().exp()));
71	break;
72	}
73	}
74	result.insert (CFAFactor (LcF, 1, 1));
75	return result;
76	}
77
78	//TODO optimize choice of p -> choose p as large as possible (better than small p since factorization mod p does not require field extension, also less lifting)
79	int
80	choosePoint (const CanonicalForm& F, int tdegF, CFArray& eval, bool rec,
81	int absValue)
82	{
83	REvaluation E1 (1, 1, IntRandom (absValue));
84	REvaluation E2 (2, 2, IntRandom (absValue));
85	if (rec)
86	{
87	E1.nextpoint();
88	E2.nextpoint();
89	}
90
91	CanonicalForm f, f1, f2, Fp;
92	int i, p;
93	CFFList f1Factors, f2Factors;
94	CFFListIterator iter;
95	int count= 0;
96	while (1)
97	{
98	count++;
99	f1= E1 (F);
100	if (!f1.isZero() && degree (f1) == degree (F,2))
101	{
102	f1Factors= factorize (f1);
103	if (f1Factors.getFirst().factor().inCoeffDomain())
104	f1Factors.removeFirst();
105	if (f1Factors.length() == 1 && f1Factors.getFirst().exp() == 1)
106	{
107	f= E2(f1);
108	f2= E2 (F);
109	f2Factors= factorize (f2);
110	Off (SW_RATIONAL);
111	if (f2Factors.getFirst().factor().inCoeffDomain())
112	f2Factors.removeFirst();
113	if (f2Factors.length() == 1 && f2Factors.getFirst().exp() == 1)
114	{
115	#ifdef HAVE_FLINT
116	// init
117	fmpz_t FLINTD1,FLINTD2;
118	fmpz_init(FLINTD1);
119	fmpz_init(FLINTD2);
120	fmpz_poly_t FLINTf1;
121	fmpz_poly_t FLINTf2;
122	// conversion
123	convertFacCF2Fmpz_poly_t(FLINTf1,f1);
124	convertFacCF2Fmpz_poly_t(FLINTf2,f2);
125	// discriminant
126	fmpz_poly_discriminant(FLINTD1,FLINTf1);
127	fmpz_poly_discriminant(FLINTD2,FLINTf2);
128	// conversion
129	CanonicalForm D1= convertFmpz2CF(FLINTD1);
130	CanonicalForm D2= convertFmpz2CF(FLINTD2);
131	// clean up
132	fmpz_poly_clear(FLINTf1);
133	fmpz_poly_clear(FLINTf2);
134	fmpz_clear(FLINTD1);
135	fmpz_clear(FLINTD2);
136	#elif defined(HAVE_NTL)
137	ZZX NTLf1= convertFacCF2NTLZZX (f1);
138	ZZX NTLf2= convertFacCF2NTLZZX (f2);
139	ZZ NTLD1= discriminant (NTLf1);
140	ZZ NTLD2= discriminant (NTLf2);
141	CanonicalForm D1= convertZZ2CF (NTLD1);
142	CanonicalForm D2= convertZZ2CF (NTLD2);
143	#else
144	factoryError("NTL/FLINT missing: choosePoint");
145	#endif
146	if ((!f.isZero()) &&
147	(abs(f)>cf_getSmallPrime (cf_getNumSmallPrimes()-1)))
148	{
149	for (i= cf_getNumPrimes()-1; i >= 0; i--)
150	{
151	if (f % CanonicalForm (cf_getPrime (i)) == 0)
152	{
153	p= cf_getPrime(i);
154	Fp= mod (F,p);
155	if (totaldegree (Fp) == tdegF &&
156	degree (mod (f2,p), 1) == degree (F,1) &&
157	degree (mod (f1, p),2) == degree (F,2))
158	{
159	if (mod (D1, p) != 0 && mod (D2, p) != 0)
160	{
161	eval[0]= E1[1];
162	eval[1]= E2[2];
163	return p;
164	}
165	}
166	}
167	}
168	}
169	else if (!f.isZero())
170	{
171	for (i= cf_getNumSmallPrimes()-1; i >= 0; i--)
172	{
173	if (f % CanonicalForm (cf_getSmallPrime (i)) == 0)
174	{
175	p= cf_getSmallPrime (i);
176	Fp= mod (F,p);
177	if (totaldegree (Fp) == tdegF &&
178	degree (mod (f2, p),1) == degree (F,1) &&
179	degree (mod (f1,p),2) == degree (F,2))
180	{
181	if (mod (D1, p) != 0 && mod (D2, p) != 0)
182	{
183	eval[0]= E1[1];
184	eval[1]= E2[2];
185	return p;
186	}
187	}
188	}
189	}
190	}
191	}
192	E2.nextpoint();
193	On (SW_RATIONAL);
194	}
195	}
196	E1.nextpoint();
197	if (count == 2)
198	{
199	count= 0;
200	absValue++;
201	E1=REvaluation (1, 1, IntRandom (absValue));
202	E2=REvaluation (2, 2, IntRandom (absValue));
203	E1.nextpoint();
204	E2.nextpoint();
205	}
206	}
207	return 0;
208	}
209
210	//G is assumed to be bivariate, irreducible over Q, primitive wrt x and y, compressed
211	CFAFList absBiFactorizeMain (const CanonicalForm& G, bool full)
212	{
213	CanonicalForm F= bCommonDen (G)*G;
214	bool isRat= isOn (SW_RATIONAL);
215	Off (SW_RATIONAL);
216	F /= icontent (F);
217	On (SW_RATIONAL);
218
219	mpz_t * M=new mpz_t [4];
220	mpz_init (M[0]);
221	mpz_init (M[1]);
222	mpz_init (M[2]);
223	mpz_init (M[3]);
224
225	mpz_t * S=new mpz_t [2];
226	mpz_init (S[0]);
227	mpz_init (S[1]);
228
229	F= compress (F, M, S);
230
231	if (F.isUnivariate())
232	{
233	if (degree (F) == 1)
234	{
235	mpz_clear (M[0]);
236	mpz_clear (M[1]);
237	mpz_clear (M[2]);
238	mpz_clear (M[3]);
239	delete [] M;
240
241	mpz_clear (S[0]);
242	mpz_clear (S[1]);
243	delete [] S;
244	if (!isRat)
245	Off (SW_RATIONAL);
246	return CFAFList (CFAFactor (G, 1, 1));
247	}
248	CFAFList result= uniAbsFactorize (F, full);
249	if (result.getFirst().factor().inCoeffDomain())
250	result.removeFirst();
251	for (CFAFListIterator iter=result; iter.hasItem(); iter++)
252	iter.getItem()= CFAFactor (decompress (iter.getItem().factor(), M, S),
253	iter.getItem().minpoly(),iter.getItem().exp());
254	mpz_clear (M[0]);
255	mpz_clear (M[1]);
256	mpz_clear (M[2]);
257	mpz_clear (M[3]);
258	delete [] M;
259
260	mpz_clear (S[0]);
261	mpz_clear (S[1]);
262	delete [] S;
263	if (!isRat)
264	Off (SW_RATIONAL);
265	return result;
266	}
267
268	if (degree (F, 1) == 1 \|\| degree (F, 2) == 1)
269	{
270	mpz_clear (M[0]);
271	mpz_clear (M[1]);
272	mpz_clear (M[2]);
273	mpz_clear (M[3]);
274	delete [] M;
275
276	mpz_clear (S[0]);
277	mpz_clear (S[1]);
278	delete [] S;
279	if (!isRat)
280	Off (SW_RATIONAL);
281	return CFAFList (CFAFactor (G, 1, 1));
282	}
283	int minTdeg, tdegF= totaldegree (F);
284	CanonicalForm Fp, smallestFactor;
285	int p;
286	CFFList factors;
287	Variable alpha;
288	bool rec= false;
289	Variable x= Variable (1);
290	Variable y= Variable (2);
291	CanonicalForm bufF= F;
292	CFFListIterator iter;
293	CFArray eval= CFArray (2);
294	int absValue= 1;
295	differentevalpoint:
296	while (1)
297	{
298	TIMING_START (fac_evalpoint);
299	p= choosePoint (F, tdegF, eval, rec, absValue);
300	TIMING_END_AND_PRINT (fac_evalpoint, "time to find eval point: ");
301
302	//after here isOn (SW_RATIONAL)==false
303	setCharacteristic (p);
304	Fp=F.mapinto();
305	factors= factorize (Fp);
306
307	if (factors.getFirst().factor().inCoeffDomain())
308	factors.removeFirst();
309
310	if (factors.length() == 1 && factors.getFirst().exp() == 1)
311	{
312	if (absIrredTest (Fp))
313	{
314	if (isRat)
315	On (SW_RATIONAL);
316	setCharacteristic(0);
317	mpz_clear (M[0]);
318	mpz_clear (M[1]);
319	mpz_clear (M[2]);
320	mpz_clear (M[3]);
321	delete [] M;
322
323	mpz_clear (S[0]);
324	mpz_clear (S[1]);
325	delete [] S;
326	return CFAFList (CFAFactor (G, 1, 1));
327	}
328	else
329	{
330	setCharacteristic (0);
331	if (modularIrredTestWithShift (F))
332	{
333	if (isRat)
334	On (SW_RATIONAL);
335	mpz_clear (M[0]);
336	mpz_clear (M[1]);
337	mpz_clear (M[2]);
338	mpz_clear (M[3]);
339	delete [] M;
340
341	mpz_clear (S[0]);
342	mpz_clear (S[1]);
343	delete [] S;
344	return CFAFList (CFAFactor (G, 1, 1));
345	}
346	rec= true;
347	continue;
348	}
349	}
350	iter= factors;
351	smallestFactor= iter.getItem().factor();
352	while (smallestFactor.isUnivariate() && iter.hasItem())
353	{
354	iter++;
355	if (!iter.hasItem())
356	break;
357	smallestFactor= iter.getItem().factor();
358	}
359
360	minTdeg= totaldegree (smallestFactor);
361	if (iter.hasItem())
362	iter++;
363	for (; iter.hasItem(); iter++)
364	{
365	if (!iter.getItem().factor().isUnivariate() &&
366	(totaldegree (iter.getItem().factor()) < minTdeg))
367	{
368	smallestFactor= iter.getItem().factor();
369	minTdeg= totaldegree (smallestFactor);
370	}
371	}
372	if (tdegF % minTdeg == 0)
373	break;
374	setCharacteristic(0);
375	rec=true;
376	}
377	CanonicalForm Gp= Fp/smallestFactor;
378	Gp= Gp /Lc (Gp);
379
380	CanonicalForm Gpy= Gp (eval[0].mapinto(), 1);
381	CanonicalForm smallestFactorEvaly= smallestFactor (eval[0].mapinto(),1);
382	CanonicalForm Gpx= Gp (eval[1].mapinto(), 2);
383	CanonicalForm smallestFactorEvalx= smallestFactor (eval[1].mapinto(),2);
384
385	bool xValid= !(Gpx.inCoeffDomain() \|\| smallestFactorEvalx.inCoeffDomain() \|\|
386	!gcd (Gpx, smallestFactorEvalx).inCoeffDomain());
387	bool yValid= !(Gpy.inCoeffDomain() \|\| smallestFactorEvaly.inCoeffDomain() \|\|
388	!gcd (Gpy, smallestFactorEvaly).inCoeffDomain());
389	if (!xValid \|\| !yValid)
390	{
391	rec= true;
392	setCharacteristic (0);
393	goto differentevalpoint;
394	}
395
396	setCharacteristic (0);
397
398	CanonicalForm mipo;
399
400	CFArray mipos= CFArray (2);
401	CFFList mipoFactors;
402	for (int i= 1; i < 3; i++)
403	{
404	CanonicalForm Fi= F(eval[i-1],i);
405
406	int s= tdegF/minTdeg + 1;
407	CanonicalForm bound= power (maxNorm (Fi), 2*(s-1));
408	bound *= power (CanonicalForm (s),s-1);
409	bound = power (CanonicalForm (2), ((s-1)(s-1))/2); //possible int overflow
410
411	CanonicalForm B = p;
412	long k = 1L;
413	while ( B < bound ) {
414	B *= p;
415	k++;
416	}
417
418	//TODO take floor (log2(k))
419	k= k+1;
420	#ifdef HAVE_FLINT
421	fmpz_poly_t FLINTFi;
422	convertFacCF2Fmpz_poly_t (FLINTFi, Fi);
423	setCharacteristic (p);
424	nmod_poly_t FLINTFpi, FLINTGpi;
425	if (i == 2)
426	{
427	convertFacCF2nmod_poly_t (FLINTFpi,
428	smallestFactorEvalx/lc (smallestFactorEvalx));
429	convertFacCF2nmod_poly_t (FLINTGpi, Gpx/lc (Gpx));
430	}
431	else
432	{
433	convertFacCF2nmod_poly_t (FLINTFpi,
434	smallestFactorEvaly/lc (smallestFactorEvaly));
435	convertFacCF2nmod_poly_t (FLINTGpi, Gpy/lc (Gpy));
436	}
437	nmod_poly_factor_t nmodFactors;
438	nmod_poly_factor_init (nmodFactors);
439	nmod_poly_factor_insert (nmodFactors, FLINTFpi, 1L);
440	nmod_poly_factor_insert (nmodFactors, FLINTGpi, 1L);
441
442	// the following fix is due to interface changes from FLINT 2.3 -> FLINT 2.4
443	# ifndef slong
444	# define slong long
445	# endif
446
447	slong * link= new slong [2];
448	fmpz_poly_t *v= new fmpz_poly_t[2];
449	fmpz_poly_t *w= new fmpz_poly_t[2];
450	fmpz_poly_init(v[0]);
451	fmpz_poly_init(v[1]);
452	fmpz_poly_init(w[0]);
453	fmpz_poly_init(w[1]);
454
455	fmpz_poly_factor_t liftedFactors;
456	fmpz_poly_factor_init (liftedFactors);
457	_fmpz_poly_hensel_start_lift (liftedFactors, link, v, w, FLINTFi,
458	nmodFactors, k);
459
460	nmod_poly_factor_clear (nmodFactors);
461	nmod_poly_clear (FLINTFpi);
462	nmod_poly_clear (FLINTGpi);
463
464	setCharacteristic(0);
465	CanonicalForm liftedSmallestFactor=
466	convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[0],x);
467
468	CanonicalForm otherFactor=
469	convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[1],x);
470	modpk pk= modpk (p, k);
471	#elif defined(HAVE_NTL)
472	modpk pk= modpk (p, k);
473	ZZX NTLFi=convertFacCF2NTLZZX (pk (Fi*pk.inverse (lc(Fi))));
474	setCharacteristic (p);
475
476	if (fac_NTL_char != p)
477	{
478	fac_NTL_char= p;
479	zz_p::init (p);
480	}
481	zz_pX NTLFpi, NTLGpi;
482	if (i == 2)
483	{
484	NTLFpi=convertFacCF2NTLzzpX (smallestFactorEvalx/lc(smallestFactorEvalx));
485	NTLGpi=convertFacCF2NTLzzpX (Gpx/lc (Gpx));
486	}
487	else
488	{
489	NTLFpi=convertFacCF2NTLzzpX (smallestFactorEvaly/lc(smallestFactorEvaly));
490	NTLGpi=convertFacCF2NTLzzpX (Gpy/lc (Gpy));
491	}
492	vec_zz_pX modFactors;
493	modFactors.SetLength(2);
494	modFactors[0]= NTLFpi;
495	modFactors[1]= NTLGpi;
496	vec_ZZX liftedFactors;
497	MultiLift (liftedFactors, modFactors, NTLFi, (long) k);
498	setCharacteristic(0);
499	CanonicalForm liftedSmallestFactor=
500	convertNTLZZX2CF (liftedFactors[0], x);
501
502	CanonicalForm otherFactor=
503	convertNTLZZX2CF (liftedFactors[1], x);
504	#else
505	factoryError("absBiFactorizeMain: NTL/FLINT missing");
506	#endif
507
508	Off (SW_SYMMETRIC_FF);
509	liftedSmallestFactor= pk (liftedSmallestFactor);
510	if (i == 2)
511	liftedSmallestFactor= liftedSmallestFactor (eval[0],1);
512	else
513	liftedSmallestFactor= liftedSmallestFactor (eval[1],1);
514
515	On (SW_SYMMETRIC_FF);
516	CFMatrix *M= new CFMatrix (s, s);
517	(*M)(s,s)= power (CanonicalForm (p), k);
518	for (int j= 1; j < s; j++)
519	{
520	(*M) (j,j)= 1;
521	(*M) (j+1,j)= -liftedSmallestFactor;
522	}
523
524	#ifdef HAVE_FLINT
525	fmpz_mat_t FLINTM;
526	convertFacCFMatrix2Fmpz_mat_t(FLINTM,*M);
527	fmpq_t delta,eta;
528	fmpq_init(delta); fmpq_set_si(delta,1,1);
529	fmpq_init(eta); fmpq_set_si(eta,3,4);
530	fmpz_mat_transpose(FLINTM,FLINTM);
531	fmpz_mat_lll_storjohann(FLINTM,delta,eta);
532	fmpz_mat_transpose(FLINTM,FLINTM);
533	delete M;
534	M=convertFmpz_mat_t2FacCFMatrix(FLINTM);
535	fmpz_mat_clear(FLINTM);
536	#elif defined(HAVE_NTL)
537	mat_ZZ * NTLM= convertFacCFMatrix2NTLmat_ZZ (*M);
538
539	ZZ det;
540
541	transpose (NTLM, NTLM);
542	(void) LLL (det, *NTLM, 1L, 1L); //use floating point LLL ?
543	transpose (NTLM, NTLM);
544	delete M;
545	M= convertNTLmat_ZZ2FacCFMatrix (*NTLM);
546	delete NTLM;
547	#else
548	factoryError("NTL/FLINT missing: absBiFactorizeMain");
549	#endif
550
551	mipo= 0;
552	for (int j= 1; j <= s; j++)
553	mipo += (M) (j,1)power (x,s-j);
554
555	delete M;
556	mipoFactors= factorize (mipo);
557	if (mipoFactors.getFirst().factor().inCoeffDomain())
558	mipoFactors.removeFirst();
559
560	#ifdef HAVE_FLINT
561	fmpz_poly_clear (v[0]);
562	fmpz_poly_clear (v[1]);
563	fmpz_poly_clear (w[0]);
564	fmpz_poly_clear (w[1]);
565	delete [] v;
566	delete [] w;
567	delete [] link;
568	fmpz_poly_factor_clear (liftedFactors);
569	#endif
570
571	if (mipoFactors.length() > 1 \|\|
572	(mipoFactors.length() == 1 && mipoFactors.getFirst().exp() > 1) \|\|
573	mipo.inCoeffDomain())
574	{
575	rec=true;
576	goto differentevalpoint;
577	}
578	else
579	mipos[i-1]= mipo;
580	}
581
582	if (degree (mipos[0]) != degree (mipos[1]))
583	{
584	rec=true;
585	goto differentevalpoint;
586	}
587
588	On (SW_RATIONAL);
589	if (maxNorm (mipos[0]) < maxNorm (mipos[1]))
590	alpha= rootOf (mipos[0]);
591	else
592	alpha= rootOf (mipos[1]);
593
594	int wrongMipo= 0;
595
596	Variable beta;
597	if (maxNorm (mipos[0]) < maxNorm (mipos[1]))
598	{
599	mipoFactors= factorize (mipos[1], alpha);
600	if (mipoFactors.getFirst().factor().inCoeffDomain())
601	mipoFactors.removeFirst();
602	for (iter= mipoFactors; iter.hasItem(); iter++)
603	{
604	if (degree (iter.getItem().factor()) > 1)
605	wrongMipo++;
606	}
607	if (wrongMipo == mipoFactors.length())
608	{
609	rec=true;
610	goto differentevalpoint;
611	}
612	wrongMipo= 0;
613	beta= rootOf (mipos[1]);
614	mipoFactors= factorize (mipos[0], beta);
615	if (mipoFactors.getFirst().factor().inCoeffDomain())
616	mipoFactors.removeFirst();
617	for (iter= mipoFactors; iter.hasItem(); iter++)
618	{
619	if (degree (iter.getItem().factor()) > 1)
620	wrongMipo++;
621	}
622	if (wrongMipo == mipoFactors.length())
623	{
624	rec=true;
625	goto differentevalpoint;
626	}
627	}
628	else
629	{
630	mipoFactors= factorize (mipos[0], alpha);
631	if (mipoFactors.getFirst().factor().inCoeffDomain())
632	mipoFactors.removeFirst();
633	for (iter= mipoFactors; iter.hasItem(); iter++)
634	{
635	if (degree (iter.getItem().factor()) > 1)
636	wrongMipo++;
637	}
638	if (wrongMipo == mipoFactors.length())
639	{
640	rec=true;
641	goto differentevalpoint;
642	}
643	wrongMipo= 0;
644	beta= rootOf (mipos[0]);
645	mipoFactors= factorize (mipos[1], beta);
646	if (mipoFactors.getFirst().factor().inCoeffDomain())
647	mipoFactors.removeFirst();
648	for (iter= mipoFactors; iter.hasItem(); iter++)
649	{
650	if (degree (iter.getItem().factor()) > 1)
651	wrongMipo++;
652	}
653	if (wrongMipo == mipoFactors.length())
654	{
655	rec=true;
656	goto differentevalpoint;
657	}
658	}
659
660
661	CanonicalForm F1;
662	if (degree (F,1) > minTdeg)
663	F1= F (eval[1], 2);
664	else
665	F1= F (eval[0], 1);
666
667	CFFList QaF1Factors;
668	bool swap= false;
669	if (F1.level() == 2)
670	{
671	swap= true;
672	F1=swapvar (F1, x, y);
673	F= swapvar (F, x, y);
674	}
675
676	wrongMipo= 0;
677	QaF1Factors= factorize (F1, alpha);
678	if (QaF1Factors.getFirst().factor().inCoeffDomain())
679	QaF1Factors.removeFirst();
680	for (iter= QaF1Factors; iter.hasItem(); iter++)
681	{
682	if (degree (iter.getItem().factor()) > minTdeg)
683	wrongMipo++;
684	}
685
686	if (wrongMipo == QaF1Factors.length())
687	{
688	rec= true;
689	F= bufF;
690	goto differentevalpoint;
691	}
692
693	CanonicalForm evaluation;
694	if (swap)
695	evaluation= eval[0];
696	else
697	evaluation= eval[1];
698
699	F *= bCommonDen (F);
700	F= F (y + evaluation, y);
701
702	int liftBound= degree (F,y) + 1;
703
704	modpk b= modpk();
705
706	CanonicalForm den= 1;
707
708	mipo= getMipo (alpha);
709
710	CFList uniFactors;
711	for (iter=QaF1Factors; iter.hasItem(); iter++)
712	uniFactors.append (iter.getItem().factor());
713
714	F /= Lc (F1);
715	#ifdef HAVE_FLINT
716	// init
717	fmpz_t FLINTf,FLINTD;
718	fmpz_init(FLINTf);
719	fmpz_init(FLINTD);
720	fmpz_poly_t FLINTmipo;
721	fmpz_poly_t FLINTLcf;
722	//conversion
723	convertFacCF2Fmpz_poly_t(FLINTmipo,mipo);
724	convertFacCF2Fmpz_poly_t(FLINTLcf,Lc (F*bCommonDen (F)));
725	// resultant, discriminant
726	fmpz_poly_resultant(FLINTf,FLINTmipo,FLINTLcf);
727	fmpz_poly_discriminant(FLINTD,FLINTmipo);
728	fmpz_mul(FLINTf,FLINTD,FLINTf);
729	den= abs (convertFmpz2CF(FLINTf));
730	// clean up
731	fmpz_clear(FLINTf);
732	// FLINTD is used below
733	fmpz_poly_clear(FLINTLcf);
734	fmpz_poly_clear(FLINTmipo);
735	#elif defined(HAVE_NTL)
736	ZZX NTLmipo= convertFacCF2NTLZZX (mipo);
737	ZZX NTLLcf= convertFacCF2NTLZZX (Lc (F*bCommonDen (F)));
738	ZZ NTLf= resultant (NTLmipo, NTLLcf);
739	ZZ NTLD= discriminant (NTLmipo);
740	den= abs (convertZZ2CF (NTLD*NTLf));
741	#else
742	factoryError("NTL/FLINT missing: absBiFactorizeMain");
743	#endif
744
745	// make factors elements of Z(a)[x] disable for modularDiophant
746	CanonicalForm multiplier= 1;
747	for (CFListIterator i= uniFactors; i.hasItem(); i++)
748	{
749	multiplier *= bCommonDen (i.getItem());
750	i.getItem()= i.getItem()*bCommonDen(i.getItem());
751	}
752	F *= multiplier;
753	F *= bCommonDen (F);
754
755	Off (SW_RATIONAL);
756	int ii= 0;
757	#ifdef HAVE_FLINT
758	CanonicalForm discMipo= convertFmpz2CF(FLINTD);
759	fmpz_clear(FLINTD);
760	#elif defined(HAVE_NTL)
761	CanonicalForm discMipo= convertZZ2CF (NTLD);
762	#endif
763	findGoodPrime (bufF*discMipo,ii);
764	findGoodPrime (F1*discMipo,ii);
765	findGoodPrime (F*discMipo,ii);
766
767	int pp=cf_getBigPrime(ii);
768	b = coeffBound( F, pp, mipo );
769	modpk bb= coeffBound (F1, pp, mipo);
770	if (bb.getk() > b.getk() ) b=bb;
771	bb= coeffBound (F, pp, mipo);
772	if (bb.getk() > b.getk() ) b=bb;
773
774	ExtensionInfo dummy= ExtensionInfo (alpha, false);
775	DegreePattern degs= DegreePattern (uniFactors);
776
777	bool earlySuccess= false;
778	CFList earlyFactors;
779	uniFactors= henselLiftAndEarly
780	(F, earlySuccess, earlyFactors, degs, liftBound,
781	uniFactors, dummy, evaluation, b, den);
782
783	DEBOUTLN (cerr, "lifted factors= " << uniFactors);
784
785	CanonicalForm MODl= power (y, liftBound);
786
787	On (SW_RATIONAL);
788	F *= bCommonDen (F);
789	Off (SW_RATIONAL);
790
791	CFList biFactors;
792
793	biFactors= factorRecombination (uniFactors, F, MODl, degs, evaluation, 1,
794	uniFactors.length()/2, b, den);
795
796	On (SW_RATIONAL);
797
798	if (earlySuccess)
799	biFactors= Union (earlyFactors, biFactors);
800	else if (!earlySuccess && degs.getLength() == 1)
801	biFactors= earlyFactors;
802
803	bool swap2= false;
804	appendSwapDecompress (biFactors, CFList(), CFList(), swap, swap2, CFMap());
805	if (isOn (SW_RATIONAL))
806	normalize (biFactors);
807
808	CFAFList result;
809	bool found= false;
810
811	for (CFListIterator i= biFactors; i.hasItem(); i++)
812	{
813	if (full)
814	result.append (CFAFactor (decompress (i.getItem(), M, S),
815	getMipo (alpha), 1));
816
817	if (totaldegree (i.getItem()) == minTdeg)
818	{
819	if (!full)
820	result= CFAFList (CFAFactor (decompress (i.getItem(), M, S),
821	getMipo (alpha), 1));
822	found= true;
823
824	if (!full)
825	break;
826	}
827	}
828
829	if (!found)
830	{
831	rec= true;
832	F= bufF;
833	goto differentevalpoint;
834	}
835
836	if (isRat)
837	On (SW_RATIONAL);
838	else
839	Off (SW_RATIONAL);
840
841	mpz_clear (M[0]);
842	mpz_clear (M[1]);
843	mpz_clear (M[2]);
844	mpz_clear (M[3]);
845	delete [] M;
846
847	mpz_clear (S[0]);
848	mpz_clear (S[1]);
849	delete [] S;
850
851	return result;
852	}
853	#endif

Note: See TracBrowser for help on using the repository browser.

Download in other formats: