Loading...
Searching...
No Matches
minpoly.h File Reference

Go to the source code of this file.

## Data Structures

class | LinearDependencyMatrix
class | NewVectorMatrix

## Macros

#define | ULONG64 (unsigned long)

## Functions

unsigned long * | computeMinimalPolynomial (unsigned long **matrix, unsigned n, unsigned long p)
unsigned long | modularInverse (long long x, long long p)
void | vectorMatrixMult (unsigned long *vec, unsigned long **mat, unsigned **nonzeroIndices, unsigned *nonzeroCounts, unsigned long *result, unsigned n, unsigned long p)
void | rem (unsigned long *a, unsigned long *q, unsigned long p, int &dega, int degq)
void | quo (unsigned long *a, unsigned long *q, unsigned long p, int &dega, int degq)
void | mult (unsigned long *result, unsigned long *a, unsigned long *b, unsigned long p, int dega, int degb)
int | gcd (unsigned long *g, unsigned long *a, unsigned long *b, unsigned long p, int dega, int degb)
int | lcm (unsigned long *l, unsigned long *a, unsigned long *b, unsigned long p, int dega, int degb)
static unsigned long | multMod (unsigned long a, unsigned long b, unsigned long p)

## Macro Definition Documentation

## ULONG64

#define ULONG64 (unsigned long)

## Function Documentation

## computeMinimalPolynomial()

unsigned long * computeMinimalPolynomial | ( | unsigned long ** | matrix,
|  | unsigned | n,
|  | unsigned long | p )

Definition at line 428 of file minpoly.cc.

{
LinearDependencyMatrix lindepmat (n, p);
NewVectorMatrix newvectormat (n, p);
unsigned long *result = new unsigned long[n + 1];
unsigned long *mpvec = new unsigned long[n + 1];
unsigned long *tmp = new unsigned long[n + 1];
// initialize result = 1
for(unsigned long i = 0; i <= n; i++)
{
result[i] = 0;
}
result[0] = 1;
int degresult = 0;
// Store the indices where the matrix has non-zero entries.
// This has a huge impact on spares matrices.
unsigned* nonzeroCounts = new unsigned[n];
unsigned** nonzeroIndices = new unsigned*[n];
for (unsigned long i = 0; i < n; i++)
{
nonzeroIndices[i] = new unsigned[n];
nonzeroCounts[i] = 0;
for (unsigned long j = 0; j < n; j++)
{
if (matrix[j][i] != 0)
{
nonzeroIndices[i][nonzeroCounts[i]] = j;
nonzeroCounts[i]++;
}
}
}
int i = n-1;
unsigned long *vec = new unsigned long[n];
unsigned long *vecnew = new unsigned long[n];
unsigned loopsEven = true;
while(i != -1)
{
for(unsigned long j = 0; j < n; j++)
{
vec[j] = 0;
}
vec[i] = 1;
lindepmat.resetMatrix ();
while(true)
{
bool ld = lindepmat.findLinearDependency (vec, mpvec);
if(ld)
{
break;
}
vectorMatrixMult (vec, matrix, nonzeroIndices, nonzeroCounts, vecnew, n, p);
unsigned long *swap = vec;
vec = vecnew;
vecnew = swap;
}
unsigned degmpvec = n;
while(mpvec[degmpvec] == 0)
{
degmpvec--;
}
// if the dimension is already maximal, i.e., the minimal polynomial of vec has the highest possible degree,
// then we are done.
if(degmpvec == n)
{
unsigned long *swap = result;
result = mpvec;
mpvec = swap;
i = -1;
}
else
{
// otherwise, we have to compute the lcm of mpvec and prior result
// tmp = zeropol
for(unsigned long j = 0; j <= n; j++)
{
tmp[j] = 0;
}
degresult = lcm (tmp, result, mpvec, p, degresult, degmpvec);
unsigned long *swap = result;
result = tmp;
tmp = swap;
if(degresult == n)
{
// minimal polynomial has highest possible degree, so we are done.
i = -1;
}
else
{
newvectormat.insertMatrix (lindepmat);
// choose new unit vector from the front or the end, alternating
// for each round. If the matrix is the companion matrix for x^n,
// then taking vectors from the end is best. If it is the transpose,
// taking vectors from the front is best.
// This tries to take the middle way
if (loopsEven)
{
i = newvectormat.findSmallestNonpivot ();
}
else
{
i = newvectormat.findLargestNonpivot ();
}
}
}
loopsEven = !loopsEven;
}
for (unsigned long i = 0; i < n; i++)
{
delete[] nonzeroIndices[i];
}
delete[] nonzeroIndices;
delete[] nonzeroCounts;
delete[]vecnew;
delete[]vec;
delete[]tmp;
delete[]mpvec;
return result;
}
swap
#define swap(_i, _j)
i
int i
Definition cfEzgcd.cc:132
p
int p
Definition cfModGcd.cc:4086
LinearDependencyMatrix
Definition minpoly.h:68
NewVectorMatrix
Definition minpoly.h:105
result
return result
Definition facAbsBiFact.cc:76
vec
fq_nmod_poly_t * vec
Definition facHensel.cc:108
j
int j
Definition facHensel.cc:110
matrix
ip_smatrix * matrix
Definition matpol.h:43
vectorMatrixMult
void vectorMatrixMult(unsigned long *vec, unsigned long **mat, unsigned **nonzeroIndices, unsigned *nonzeroCounts, unsigned long *result, unsigned n, unsigned long p)
Definition minpoly.cc:393
lcm
int lcm(unsigned long *l, unsigned long *a, unsigned long *b, unsigned long p, int dega, int degb)
Definition minpoly.cc:709

## gcd()

int gcd | ( | unsigned long * | g,
|  | unsigned long * | a,
|  | unsigned long * | b,
|  | unsigned long | p,
|  | int | dega,
|  | int | degb )

Definition at line 666 of file minpoly.cc.

{
unsigned long *tmp1 = new unsigned long[dega + 1];
unsigned long *tmp2 = new unsigned long[degb + 1];
for(int i = 0; i <= dega; i++)
{
tmp1[i] = a[i];
}
for(int i = 0; i <= degb; i++)
{
tmp2[i] = b[i];
}
int degtmp1 = dega;
int degtmp2 = degb;
unsigned long *swappol;
int swapdeg;
while(degtmp2 >= 0)
{
rem (tmp1, tmp2, p, degtmp1, degtmp2);
swappol = tmp1;
tmp1 = tmp2;
tmp2 = swappol;
swapdeg = degtmp1;
degtmp1 = degtmp2;
degtmp2 = swapdeg;
}
for(int i = 0; i <= degtmp1; i++)
{
g[i] = tmp1[i];
}
delete[]tmp1;
delete[]tmp2;
return degtmp1;
}
g
g
Definition cfModGcd.cc:4098
b
CanonicalForm b
Definition cfModGcd.cc:4111
tmp1
CFList tmp1
Definition facFqBivar.cc:75
tmp2
CFList tmp2
Definition facFqBivar.cc:75
rem
void rem(unsigned long *a, unsigned long *q, unsigned long p, int &dega, int degq)
Definition minpoly.cc:572

## lcm()

int lcm | ( | unsigned long * | l,
|  | unsigned long * | a,
|  | unsigned long * | b,
|  | unsigned long | p,
|  | int | dega,
|  | int | degb )

Definition at line 709 of file minpoly.cc.

{
unsigned long *g = new unsigned long[dega + 1];
// initialize entries of g to zero
for(int i = 0; i <= dega; i++)
{
g[i] = 0;
}
int degg = gcd (g, a, b, p, dega, degb);
if(degg > 0)
{
// non-trivial gcd, so compute a = (a/g)
quo (a, g, p, dega, degg);
}
mult (l, a, b, p, dega, degb);
// normalize
if(l[dega + degb + 1] != 1)
{
unsigned long inv = modularInverse (l[dega + degb], p);
for(int i = 0; i <= dega + degb; i++)
{
l[i] = multMod (l[i], inv, p);
}
}
return dega + degb;
}
l
int l
Definition cfEzgcd.cc:100
degg
int degg
Definition facAlgExt.cc:64
mult
void mult(unsigned long *result, unsigned long *a, unsigned long *b, unsigned long p, int dega, int degb)
Definition minpoly.cc:647
quo
void quo(unsigned long *a, unsigned long *q, unsigned long p, int &dega, int degq)
Definition minpoly.cc:597
modularInverse
unsigned long modularInverse(long long x, long long p)
Definition minpoly.cc:744
gcd
int gcd(unsigned long *g, unsigned long *a, unsigned long *b, unsigned long p, int dega, int degb)
Definition minpoly.cc:666
multMod
static unsigned long multMod(unsigned long a, unsigned long b, unsigned long p)
Definition minpoly.h:202

## modularInverse()

unsigned long modularInverse | ( | long long | x,
|  | long long | p )

Definition at line 744 of file minpoly.cc.

{
long long u1 = 1;
long long u2 = 0;
long long u3 = x;
long long v1 = 0;
long long v2 = 1;
long long v3 = p;
long long q, t1, t2, t3;
while(v3 != 0)
{
q = u3 / v3;
t1 = u1 - q * v1;
t2 = u2 - q * v2;
t3 = u3 - q * v3;
u1 = v1;
u2 = v2;
u3 = v3;
v1 = t1;
v2 = t2;
v3 = t3;
}
if(u1 < 0)
{
u1 += p;
}
return (unsigned long) u1;
}
x
Variable x
Definition cfModGcd.cc:4090

## mult()

void mult | ( | unsigned long * | result,
|  | unsigned long * | a,
|  | unsigned long * | b,
|  | unsigned long | p,
|  | int | dega,
|  | int | degb )

Definition at line 647 of file minpoly.cc.

{
// NOTE: we assume that every entry in result is preinitialized to zero!
for(int i = 0; i <= dega; i++)
{
for(int j = 0; j <= degb; j++)
{
result[i + j] += multMod (a[i], b[j], p);
if (result[i + j] >= p)
{
result[i + j] -= p;
}
}
}
}

## multMod()

unsigned long multMod | ( | unsigned long | a,
|  | unsigned long | b,
|  | unsigned long | p )

inline static

Definition at line 202 of file minpoly.h.

{
#if SIZEOF_LONG == 4
#define ULONG64 (unsigned long long)
#else
#define ULONG64 (unsigned long)
#endif
return (unsigned long)((ULONG64 a)*(ULONG64 b) % (ULONG64 p));
}
ULONG64
#define ULONG64

## quo()

void quo | ( | unsigned long * | a,
|  | unsigned long * | q,
|  | unsigned long | p,
|  | int & | dega,
|  | int | degq )

Definition at line 597 of file minpoly.cc.

{
unsigned degres = dega - degq;
unsigned long *result = new unsigned long[degres + 1];
// initialize to zero
for (unsigned i = 0; i <= degres; i++)
{
result[i] = 0;
}
while(degq <= dega)
{
unsigned d = dega - degq;
unsigned long inv = modularInverse (q[degq], p);
result[d] = multMod (a[dega], inv, p);
for(int i = degq; i >= 0; i--)
{
unsigned long tmp = p - multMod (result[d], q[i], p);
a[d + i] += tmp;
if (a[d + i] >= p)
{
a[d + i] -= p;
}
}
while(dega >= 0 && a[dega] == 0)
{
dega--;
}
}
// copy result into a
for(unsigned i = 0; i <= degres; i++)
{
a[i] = result[i];
}
// set all following entries of a to zero
for(unsigned i = degres + 1; i <= degq + degres; i++)
{
a[i] = 0;
}
dega = degres;
delete[]result;
}

## rem()

void rem | ( | unsigned long * | a,
|  | unsigned long * | q,
|  | unsigned long | p,
|  | int & | dega,
|  | int | degq )

Definition at line 572 of file minpoly.cc.

{
while(degq <= dega)
{
unsigned d = dega - degq;
long factor = multMod (a[dega], modularInverse (q[degq], p), p);
for(int i = degq; i >= 0; i--)
{
long tmp = p - multMod (factor, q[i], p);
a[d + i] += tmp;
if (a[d + i] >= p)
{
a[d + i] -= p;
}
}
while(dega >= 0 && a[dega] == 0)
{
dega--;
}
}
}
factor
CanonicalForm factor
Definition facAbsFact.cc:97

## vectorMatrixMult()

void vectorMatrixMult | ( | unsigned long * | vec,
|  | unsigned long ** | mat,
|  | unsigned ** | nonzeroIndices,
|  | unsigned * | nonzeroCounts,
|  | unsigned long * | result,
|  | unsigned | n,
|  | unsigned long | p )

Definition at line 393 of file minpoly.cc.

{
unsigned long tmp;
for(unsigned long i = 0; i < n; i++)
{
result[i] = 0;
for(unsigned j = 0; j < nonzeroCounts[i]; j++)
{
tmp = multMod (vec[nonzeroIndices[i][j]], mat[nonzeroIndices[i][j]][i], p);
result[i] += tmp;
if (result[i] >= p)
{
result[i] -= p;
}
}
}
}