Loading...
Searching...
No Matches
facSparseHensel.cc
Go to the documentation of this file.
/*****************************************************************************\
* Computer Algebra System SINGULAR
\*****************************************************************************/
/** @file facSparseHensel.cc
*
* This file implements functions for sparse heuristic Hensel lifting
*
* ABSTRACT: "A fast implementation of polynomial factorization" by M. Lucks and
* "Effective polynomial computation" by R. Zippel
*
* @author Martin Lee
*
**/
/*****************************************************************************/
#include "config.h"
#include "cf_assert.h"
#include "facSparseHensel.h"
#include "cf_algorithm.h"
#include "cfModGcd.h"
#include "facFqFactorize.h"
int
LucksWangSparseHeuristic (const CanonicalForm& F, const CFList& factors,
int level, const CFList& leadingCoeffs, CFList& result)
{
int threshold= 450;
CFArray termsF= getBiTerms (F, threshold);
int si=termsF.size();
int fl=factors.length();
//printf("size:%d, length=%d, char=%d\n",si,fl,getCharacteristic());
if ((si > threshold)
|| (si <= fl)
)
return 0;
sort (termsF, level);
CFArray* monoms= new CFArray [factors.length()];
int i= 0;
int num= 0;
for (CFListIterator iter= factors; iter.hasItem(); iter++, i++)
{
monoms[i]= getTerms (iter.getItem());
num += monoms [i].size();
sort (monoms [i]);
}
i= 0;
CFArray* monomsLead= new CFArray [leadingCoeffs.length()];
for (CFListIterator iter= leadingCoeffs; iter.hasItem(); iter++, i++)
{
monomsLead[i]= getTerms (iter.getItem());
sort (monomsLead [i]);
groupTogether (monomsLead [i], level);
strip (monomsLead [i], level);
}
CFArray solution= CFArray (num);
int k, d, count, j= F.level() + 1;
num= 0;
i= 0;
for (CFListIterator iter= factors; iter.hasItem(); i++, iter++)
{
d= degree (iter.getItem(), 1);
count= 0;
for (k= 0; k < monoms[i].size(); k++, j++, num++)
{
monoms [i][k] *= Variable (j);
if (degree (monoms[i][k], 1) == d)
{
solution[num]= monomsLead [i] [count];
count++;
}
}
}
delete [] monomsLead;
CFList tmp;
CFArray* stripped2= new CFArray [factors.length()];
for (i= factors.length() - 1; i > -1; i--)
{
tmp.insert (buildPolyFromArray (monoms [i]));
strip (monoms[i], stripped2 [i], level);
}
delete [] monoms;
CanonicalForm H= prod (tmp);
CFArray monomsH= getMonoms (H);
sort (monomsH,F.level());
groupTogether (monomsH, F.level());
if (monomsH.size() != termsF.size())
{
delete [] stripped2;
return 0;
}
CFArray strippedH;
strip (monomsH, strippedH, level);
CFArray strippedF;
strip (termsF, strippedF, level);
if (!isEqual (strippedH, strippedF))
{
delete [] stripped2;
return 0;
}
CFArray A= getEquations (monomsH, termsF);
CFArray startingSolution= solution;
CFArray newSolution= CFArray (solution.size());
result= CFList();
do
{
evaluate (A, solution, F.level() + 1);
if (isZero (A))
break;
if (!simplify (A, newSolution, F.level() + 1))
{
delete [] stripped2;
return 0;
}
if (isZero (newSolution))
break;
if (!merge (solution, newSolution))
break;
} while (1);
if (isEqual (startingSolution, solution))
{
delete [] stripped2;
return 0;
}
CanonicalForm factor;
num= 0;
for (i= 0; i < factors.length(); i++)
{
k= stripped2[i].size();
factor= 0;
for (j= 0; j < k; j++, num++)
{
if (solution [num].isZero())
continue;
factor += solution [num]*stripped2[i][j];
}
result.append (factor);
}
delete [] stripped2;
if (result.length() > 0)
return 1;
return 0;
}
CFList
sparseHeuristic (const CanonicalForm& A, const CFList& biFactors,
CFList*& moreBiFactors, const CFList& evaluation,
int minFactorsLength)
{
int j= A.level() - 1;
int i;
//initialize storage
CFArray *** storeFactors= new CFArray** [j];
for (i= 0; i < j; i++)
storeFactors [i]= new CFArray* [2];
CFArray eval= CFArray (j);
i= j - 1;
for (CFListIterator iter= evaluation; iter.hasItem(); iter++,i--)
eval[i]= iter.getItem();
storeFactors [0] [0]= new CFArray [minFactorsLength];
storeFactors [0] [1]= new CFArray [minFactorsLength];
for (i= 1; i < j; i++)
{
storeFactors[i] [0]= new CFArray [minFactorsLength];
storeFactors[i] [1]= new CFArray [minFactorsLength];
}
//
CFList * normalizingFactors= new CFList [j];
CFList uniFactors;
normalizingFactors [0]= findNormalizingFactor1 (biFactors,
evaluation.getLast(), uniFactors);
for (i= j - 1; i > 0; i--)
{
if (moreBiFactors[i-1].length() != minFactorsLength)
{
moreBiFactors[i-1]=
recombination (moreBiFactors [i-1], uniFactors, 1,
moreBiFactors[i-1].length()-uniFactors.length()+1,
eval[i], Variable (i + 2)
);
}
normalizingFactors [i]= findNormalizingFactor2 (moreBiFactors [i - 1],
eval[i], uniFactors);
}
CFList tmp;
tmp= normalize (biFactors, normalizingFactors[0]);
getTerms2 (tmp, storeFactors [0] [0]);
storeFactors [0] [1]= evaluate (storeFactors [0] [0], minFactorsLength,
evaluation.getLast(), Variable (2));
for (i= j - 1; i > 0; i--)
{
tmp= normalize (moreBiFactors [i-1], normalizingFactors [i]);
getTerms2 (tmp, storeFactors [i] [0]);
storeFactors [i] [1]= evaluate (storeFactors [i] [0], minFactorsLength,
eval[i], Variable (i + 2));
}
int k, l, m, mm, count, sizeOfUniFactors= 0;
int*** seperator= new int** [j];
Variable x= Variable (1);
for (i= 0; i < j; i++)
seperator [i]= new int* [minFactorsLength];
for (k= 0; k < minFactorsLength; k++)
{
for (i= 0; i < j; i++)
{
count= 0;
for (l= 0; l < storeFactors [i][0][k].size() - 1; l++)
{
if (degree (storeFactors[i][0][k][l], x) <
degree (storeFactors[i][0][k][l+1], x))
count++;
}
if (i == 0)
sizeOfUniFactors= count;
else if (sizeOfUniFactors != count)
{
for (m= 0; m < j; m++)
{
delete [] storeFactors [m] [0];
delete [] storeFactors [m] [1];
delete [] storeFactors [m];
for (mm= 0; mm < k; mm++)
delete [] seperator [m][mm];
delete [] seperator [m];
}
delete [] storeFactors;
delete [] seperator;
delete [] normalizingFactors;
return CFList();
}
seperator [i][k]= new int [count + 3];
seperator [i][k][0]= count + 1;
seperator [i][k][1]= 0;
count= 2;
for (l= 0; l < storeFactors [i][0][k].size() - 1; l++)
{
if (degree (storeFactors[i][0][k][l], x) <
degree (storeFactors[i][0][k][l+1], x))
{
seperator[i][k][count]=l + 1;
count++;
}
}
seperator [i][k][count]= storeFactors[i][0][k].size();
}
}
CanonicalForm tmp1, factor, quot;
CanonicalForm B= A;
CFList result;
int maxTerms, n, index1, index2, mmm, found, columns, oneCount;
int ** mat;
for (k= 0; k < minFactorsLength; k++)
{
factor= 0;
sizeOfUniFactors= seperator [0][k][0];
for (n= 1; n <= sizeOfUniFactors; n++)
{
columns= 0;
maxTerms= 1;
index1= j - 1;
for (i= j - 1; i >= 0; i--)
{
if (maxTerms < seperator[i][k][n+1]-seperator[i][k][n])
{
maxTerms= seperator[i][k][n + 1]-seperator[i][k][n];
index1= i;
}
}
for (i= j - 1; i >= 0; i--)
{
if (i == index1)
continue;
columns += seperator [i][k][n+1]-seperator[i][k][n];
}
mat= new int *[maxTerms];
mm= 0;
for (m= seperator[index1][k][n]; m < seperator[index1][k][n+1]; m++, mm++)
{
tmp1= storeFactors [index1][1][k][m];
mat[mm]= new int [columns];
for (i= 0; i < columns; i++)
mat[mm][i]= 0;
index2= 0;
for (i= j - 1; i >= 0; i--)
{
if (i == index1)
continue;
found= -1;
if ((found= search (storeFactors[i][1][k], tmp1,
seperator[i][k][n], seperator[i][k][n+1])) >= 0)
mat[mm][index2 + found - seperator [i][k][n]]= 1;
index2 += seperator [i][k][n+1]-seperator[i][k][n];
}
}
index2= 0;
for (i= j - 1; i >= 0; i--)
{
if (i == index1)
continue;
oneCount= 0;
for (mm= 0; mm < seperator [i][k][n + 1] - seperator [i][k][n]; mm++)
{
for (m= 0; m < maxTerms; m++)
{
if (mat[m][mm+index2] == 1)
oneCount++;
}
}
if (oneCount == seperator [i][k][n+1]-seperator[i][k][n] - 1)
{
for (mm= 0; mm < seperator [i][k][n+1]-seperator[i][k][n]; mm++)
{
oneCount= 0;
for (m= 0; m < maxTerms; m++)
if (mat[m][mm+index2] == 1)
oneCount++;
if (oneCount > 0)
continue;
for (m= 0; m < maxTerms; m++)
{
oneCount= 0;
for (mmm= 0; mmm < seperator[i][k][n+1]-seperator[i][k][n]; mmm++)
{
if (mat[m][mmm+index2] == 1)
oneCount++;
}
if (oneCount > 0)
continue;
mat[m][mm+index2]= 1;
}
}
}
index2 += seperator [i][k][n+1] - seperator [i][k][n];
}
//read off solution
mm= 0;
for (m= seperator[index1][k][n]; m < seperator[index1][k][n+1]; m++, mm++)
{
tmp1= storeFactors [index1][0][k][m];
index2= 0;
for (i= j - 1; i > -1; i--)
{
if (i == index1)
continue;
for (mmm= 0; mmm < seperator [i][k][n+1]-seperator[i][k][n]; mmm++)
if (mat[mm][mmm+index2] == 1)
tmp1= patch (tmp1, storeFactors[i][0][k][seperator[i][k][n]+mmm],
eval[i]);
index2 += seperator [i][k][n+1]-seperator[i][k][n];
}
factor += tmp1;
}
for (m= 0; m < maxTerms; m++)
delete [] mat [m];
delete [] mat;
}
if (fdivides (factor, B, quot))
{
result.append (factor);
B= quot;
if (result.length() == biFactors.length() - 1)
{
result.append (quot);
break;
}
}
}
//delete
for (i= 0; i < j; i++)
{
delete [] storeFactors [i] [0];
delete [] storeFactors [i] [1];
delete [] storeFactors [i];
for (k= 0; k < minFactorsLength; k++)
delete [] seperator [i][k];
delete [] seperator [i];
}
delete [] seperator;
delete [] storeFactors;
delete [] normalizingFactors;
//
return result;
}
degree
int degree(const CanonicalForm &f)
Definition canonicalform.h:316
CFArray
Array< CanonicalForm > CFArray
Definition canonicalform.h:397
CFListIterator
ListIterator< CanonicalForm > CFListIterator
Definition canonicalform.h:396
num
CanonicalForm num(const CanonicalForm &f)
Definition canonicalform.h:337
level
int level(const CanonicalForm &f)
Definition canonicalform.h:331
CFList
List< CanonicalForm > CFList
Definition canonicalform.h:395
l
int l
Definition cfEzgcd.cc:100
m
int m
Definition cfEzgcd.cc:128
i
int i
Definition cfEzgcd.cc:132
k
int k
Definition cfEzgcd.cc:99
isEqual
bool isEqual(int *a, int *b, int lower, int upper)
Definition cfGcdAlgExt.cc:947
getMonoms
CFArray getMonoms(const CanonicalForm &F)
extract monomials of F, parts in algebraic variable are considered coefficients
Definition cfModGcd.cc:1965
x
Variable x
Definition cfModGcd.cc:4090
evaluate
CFArray evaluate(const CFArray &A, const CFList &evalPoints)
Definition cfModGcd.cc:2039
cfModGcd.h
modular and sparse modular GCD algorithms over finite fields and Z.
sort
static void sort(int **points, int sizePoints)
Definition cfNewtonPolygon.cc:100
merge
int ** merge(int **points1, int sizePoints1, int **points2, int sizePoints2, int &sizeResult)
Definition cfNewtonPolygon.cc:230
cf_algorithm.h
declarations of higher level algorithms.
getTerms
void getTerms(const CanonicalForm &f, const CanonicalForm &t, CFList &result)
get_Terms: Split the polynomial in the containing terms.
Definition cf_factor.cc:283
fdivides
bool fdivides(const CanonicalForm &f, const CanonicalForm &g)
cf_assert.h
assertions for Factory
Array::size
int size() const
Definition ftmpl_array.cc:92
CanonicalForm
factory's main class
Definition canonicalform.h:86
CanonicalForm::level
int level() const
level() returns the level of CO.
Definition canonicalform.cc:576
List::length
int length() const
Definition ftmpl_list.cc:273
List::insert
void insert(const T &)
Definition ftmpl_list.cc:193
Variable
factory's class for variables
Definition factory.h:127
iter
CFFListIterator iter
Definition facAbsBiFact.cc:54
result
return result
Definition facAbsBiFact.cc:76
evaluation
const CanonicalForm int const CFList & evaluation
Definition facAbsFact.cc:52
factor
CanonicalForm factor
Definition facAbsFact.cc:97
H
CanonicalForm H
Definition facAbsFact.cc:60
B
b *CanonicalForm B
Definition facBivar.cc:52
found
bool found
Definition facFactorize.cc:55
eval
CFList & eval
Definition facFactorize.cc:47
minFactorsLength
CFList int & minFactorsLength
[in,out] minimal length of bivariate factors
Definition facFactorize.h:33
tmp1
CFList tmp1
Definition facFqBivar.cc:75
recombination
CFList recombination(const CFList &factors1, const CFList &factors2, int s, int thres, const CanonicalForm &evalPoint, const Variable &x)
recombination of bivariate factors factors1 s. t. the result evaluated at evalPoint coincides with fa...
Definition facFqFactorize.cc:2114
facFqFactorize.h
This file provides functions for factorizing a multivariate polynomial over , or GF.
j
int j
Definition facHensel.cc:110
prod
fq_nmod_poly_t prod
Definition facHensel.cc:100
LucksWangSparseHeuristic
int LucksWangSparseHeuristic(const CanonicalForm &F, const CFList &factors, int level, const CFList &leadingCoeffs, CFList &result)
sparse heuristic lifting by Wang and Lucks
Definition facSparseHensel.cc:26
sparseHeuristic
CFList sparseHeuristic(const CanonicalForm &A, const CFList &biFactors, CFList *&moreBiFactors, const CFList &evaluation, int minFactorsLength)
sparse heuristic which patches together bivariate factors of A wrt. different second variables by the...
Definition facSparseHensel.cc:160
facSparseHensel.h
This file provides functions for sparse heuristic Hensel lifting.
buildPolyFromArray
CanonicalForm buildPolyFromArray(const CFArray &A)
build a poly from entries in A
Definition facSparseHensel.h:267
groupTogether
void groupTogether(CFArray &A, int level)
group together elements in A, where entries in A are put together if they coincide up to level level
Definition facSparseHensel.h:278
getTerms2
CFArray getTerms2(const CanonicalForm &F)
get terms of F wrt. Variable (1)
Definition facSparseHensel.h:492
getBiTerms
CFArray getBiTerms(const CanonicalForm &F, int threshold)
get terms of F where F is considered a bivariate poly in Variable(1), Variable (2)
Definition facSparseHensel.h:236
simplify
CanonicalForm simplify(const CanonicalForm &A, int level)
simplify A if possible, i.e. A consists of 2 terms and contains only one variable of level greater or...
Definition facSparseHensel.h:390
isZero
bool isZero(const CFArray &A)
checks if entries of A are zero
Definition facSparseHensel.h:468
findNormalizingFactor1
CFList findNormalizingFactor1(const CFList &biFactors, const CanonicalForm &evalPoint, CFList &uniFactors)
find normalizing factors for biFactors and build monic univariate factors from biFactors
Definition facSparseHensel.h:123
findNormalizingFactor2
CFList findNormalizingFactor2(CFList &biFactors, const CanonicalForm &evalPoint, const CFList &uniFactors)
find normalizing factors for biFactors and sort biFactors s.t. the returned biFactors evaluated at ev...
Definition facSparseHensel.h:140
search
int search(const CFArray &A, const CanonicalForm &F, int i, int j)
search for F in A between index i and j
Definition facSparseHensel.h:566
patch
CanonicalForm patch(const CanonicalForm &F1, const CanonicalForm &F2, const CanonicalForm &eval)
patch together F1 and F2 and normalize by a power of eval F1 and F2 are assumed to be bivariate with ...
Definition facSparseHensel.h:577
strip
void strip(CFArray &F, CFArray &G, int level)
strip off those parts of entries in F whose level is less than or equal than level and store the stri...
Definition facSparseHensel.h:310
getEquations
CFArray getEquations(const CFArray &A, const CFArray &B)
get equations for LucksWangSparseHeuristic
Definition facSparseHensel.h:350
length
static BOOLEAN length(leftv result, leftv arg)
Definition interval.cc:257
count
int status int void size_t count
Definition si_signals.h:69
A
#define A
Definition sirandom.c:24
normalize
static poly normalize(poly next_p, ideal add_generators, syStrategy syzstr, int *g_l, int *p_l, int crit_comp)
Definition syz3.cc:1027