1
#ifndef LINALG_H_INCLUDED
2
#define LINALG_H_INCLUDED
6
#include "field_rationals.h"
11
Unfortunately it seems that these classes cannot share code with the
12
Vektor and Matrix templates since the classes in this file carry
13
typeinfo around. In particular it is not possible to use the Matrix
14
template on the FieldElement type since a FieldElement has no
15
default initialization from a double -- its initialization really
16
depends on which Field we are using.
22
@brief A vector containing FieldElement s.
24
When a FieldVector is constructed the Field of which the entries of
25
the vector will be members must be specified either explicitly or
31
The field of which the entries of the vector are members.
35
The length of the vector.
39
The entries of the vector.
41
vector<FieldElement> v;
44
Constructs a FieldVector with entries initialized to 0.
46
@param _theField the Field of which the entries of the vector will be members.
47
@param _n the number of entries in the vector.
49
FieldVector(Field const &_theField, int _n);
51
* Produces a vector with all entries equal to 1.
53
static FieldVector allOnes(Field const &_theField, int _n);
57
inline Field const&getField()const{return theField;}
59
Provides access to i th entry of the vector.
61
FieldElement& operator[](int i);
63
Provides access to i th entry of the vector.
65
FieldElement const& operator[](int i)const;
67
Returns true if all entries of the vector equal the 0-element of the Field, false otherwise.
71
Computes the first integer vector on the half ray spanned by the
72
vector. Assumes that the Field is Q, asserts if this is not the
73
case. Asserts if the computed vector has entries that do not fit
76
IntegerVector primitive()const;
78
Returns the number of entries in the vector.
82
Returns the subvector (v_begin,\dots,v_end-1).
84
FieldVector subvector(int begin, int end)const;
86
Returns the subvector .
88
FieldVector subvector(list<int> const &l)const;
90
Returns the concatenation of the two vectors.
92
friend FieldVector concatenation(FieldVector const &a, FieldVector const &b);
94
Computes the dot product of the vectors a and b.
96
friend FieldElement dot(FieldVector const &a, FieldVector const &b);
98
Computes the vector q scaled by a factor s.
100
friend FieldVector operator*(FieldElement s, const FieldVector& q);
102
Computes the coordinate-wise subtraction a-b.
104
friend FieldVector operator-(FieldVector const &a, FieldVector const &b);
106
Computes the coordinate-wise division of a by b.
108
friend FieldVector operator/(FieldVector const &a, FieldVector const &b);
110
Adds q to the current vector coordinatewise. Asserts if the
111
current vector and q do not have the same number of elements.
113
FieldVector& operator+=(const FieldVector& q);
115
Adds s*q to the current vector coordinatewise. Asserts if the
116
current vector and q do not have the same number of elements.
118
FieldVector& madd(const FieldElement &s, const FieldVector& q);
120
Prints the vector using the Printer P.
122
void print(Printer &P)const;
124
friend IntegerVector fieldVectorToIntegerVector(FieldVector const &v);//MUST CONTAIN INTEGER ENTRIES ONLY // should this be a member?
125
bool operator==(FieldVector const &b)const;
130
@brief A matrix containing FieldElement s.
132
When a FieldMatrix is constructed the Field of which the entries of
133
the matrix will be members must be specified either explicitly or
134
implicitly. The FieldMatrix class contains routines for performing
141
vector<FieldVector> rows;
144
Constructs a FieldMatrix with entries initialized to 0.
146
@param _theField the Field of which the entries of the matrix will be members.
147
@param _height the number of rows of the matrix.
148
@param _width the number of columns of the matrix.
150
FieldMatrix(Field const &_theField, int _height, int _width);
152
Returns the nxn identity matrix.
154
static FieldMatrix identity(Field const &_theField, int _n);
158
inline Field const&getField()const{return theField;}
160
Provides access to i th row of the matrix which is a FieldVector.
162
FieldVector& operator[](int i);
164
Provides access to i th row of the matrix which is a FieldVector.
166
FieldVector const& operator[](int i)const;
168
Returns the number of rows of the matrix.
170
int getHeight()const;
172
Returns the number of columns of the matrix.
176
Swaps the i th and the j th row.
178
void swapRows(int i, int j);
180
Returns the submatrix consisting of rows index by rowIndices.
182
FieldMatrix submatrixRows(list<int> const &rowIndices)const;
184
Adds a times the i th row to the j th row.
186
void madd(int i, FieldElement a, int j);
188
Returns the index of a row whose index is at least currentRow and
189
which has a non-zero element on the column th entry. If no such
190
row exists then -1 is returned. This routine is used in the Gauss
191
reduction. To make the reduction more efficient the routine
192
chooses its row with as few non-zero entries as possible.
194
int findRowIndex(int column, int currentRow)const;
196
Prints the matrix using the Printer P.
198
void printMatrix(Printer &P)const;
200
This method is used for iterating through the pivots in a matrix
201
in row echelon form. To find the first pivot put i=-1 and
202
j=-1 and call this routine. The (i,j) th entry of the matrix is a
203
pivot. Call the routine again to get the next pivot. When no more
204
pivots are found the routine returns false.
206
bool nextPivot(int &i, int &j)const;
208
Performs a Gauss reduction and returns the number of row swaps
209
done. The result is a matrix in row echelon form. The pivots may
210
not be all 1. In terms of Groebner bases, what is computed is a
211
minimal (not necessarily reduced) Groebner basis of the linear
212
ideal generated by the rows. The number of swaps is need if one
213
wants to compute the determinant afterwards. In this case it is
214
also a good idea to set the flag returnIfZeroDeterminant which
215
make the routine terminate before completion if it discovers that
216
the determinant is zero.
218
int reduce(bool returnIfZeroDeterminant=false, bool hermite=false); //bool reducedRowEcholonForm,
220
Calls reduce() and returns the determinant of the matrix. If
221
the matrix is not square the routine asserts.
223
FieldElement reduceAndComputeDeterminant();
225
Calls reduce() and returns the rank of the matrix.
227
int reduceAndComputeRank();
229
Scales every row of the matrix so that the first non-zero entry of
232
void scaleFirstNonZeroEntryToOne();
234
Assumes that the matrix has a kernel of dimension 1.
235
Calls reduce() and returns a non-zero vector in the kernel.
236
If the matrix is an (n-1)x(n) matrix then the returned vector has
237
the property that if appended as a row to the original matrix
238
then the determinant of that matrix would be positive. Of course
239
this property, as described here, only makes sense for ordered fields.
241
FieldVector reduceAndComputeVectorInKernel();
243
Calls reduce() and constructs matrix whose rows forms a basis of
244
the kernel of the linear map defined by the original matrix. The
245
return value is the new matrix.
247
FieldMatrix reduceAndComputeKernel();
249
Returns the transposed of the the matrix.
251
FieldMatrix transposed()const;
253
Returns the matrix with rows in reversed order.
255
FieldMatrix flipped()const;
257
Takes two matrices with the same number of columns and construct
258
a new matrix which has the rows of the matrix top on the top and
259
the rows of the matrix bottom at the bottom. The return value is
260
the constructed matrix.
262
friend FieldMatrix combineOnTop(FieldMatrix const &top, FieldMatrix const &bottom);
264
Takes two matrices with the same number of rows and construct
265
a new matrix which has the columns of the matrix left on the left and
266
the columns of the matrix right on the right. The return value is
267
the constructed matrix.
269
friend FieldMatrix combineLeftRight(FieldMatrix const &left, FieldMatrix const &right);
271
Computes a reduced row echelon form from a row echelon form. In
272
Groebner basis terms this is the same as tranforming a minimal
273
Groebner basis to a reduced one except that we do not force
274
pivots to be 1 unless the scalePivotsToOne parameter is set.
276
int REformToRREform(bool scalePivotsToOne=false);
278
This function may be called if the FieldMatrix is in Row Echelon
279
Form. The input is a FieldVector which is rewritten modulo the
280
rows of the matrix. The result is unique and is the same as the
281
normal form of the input vector modulo the Groebner basis of the
282
linear ideal generated by the rows of the matrix.
284
FieldVector canonicalize(FieldVector v)const;
286
This routine removes the zero rows of the matrix.
288
void removeZeroRows();
290
This routine produces a matrix in row Echelon form solving a
291
system of the kind Ax=b, where A equals *this. The way this is
292
done is by considering [A^t -I] and computing its RE form. To
293
solve the system afterwards simply call canonicalize of the new
294
matrix on (b_1,\dots,b_m,0,\dots,0) to get
295
(0,\dots,0,x_1,\dots,x_n). If no solution exists then the first
296
m entries of the returned vector are not all 0.
298
FieldMatrix solver()const;
300
Computes the vector m scaled by a factor s.
302
friend FieldMatrix operator*(FieldElement s, const FieldMatrix& m);
304
Returns the matrix product a*b.
306
friend FieldMatrix operator*(FieldMatrix const &a, const FieldMatrix &b);
308
Returns the vector a*b, where a is considered as a row vector.
310
friend FieldVector operator*(FieldVector const &a, const FieldMatrix &b);
312
Compares two matrices.
314
bool operator==(FieldMatrix const &b)const;
316
Returns the rank of the matrix. This function is const and slow -
317
will copy the matrix and perform Gauss reduction.
321
Returns the indices of the columns containing a pivot.
322
The returned list is sorted.
323
The matrix must be in row echelon form.
325
list<int> pivotColumns()const;
327
Returns the indices of the columns not containing a pivot.
328
The returned list is sorted.
329
The matrix must be in row echelon form.
331
list<int> nonPivotColumns()const;
333
* Returns the specified submatrix with size (endRow-startRow)x(endColumn-startColumn).
335
FieldMatrix submatrix(int startRow, int startColumn, int endRow, int endColumn)const;
337
* Thinking of an upper triangular matrix as a Groebner basis this routine returns the normal form
338
* of the linear polynomial represented by v. If coeff is non-zero, then that vector is assigned
339
* the coefficients used in the reduction. That is, if the normal form is zero, then multiplying coeff
340
* and *this gives v. The matrix must be in row-echelon form. coeff does not have to have the right length before the call.
342
FieldVector normalForm(FieldVector const &v, FieldVector *coeff)const;
346
FieldMatrix integerMatrixToFieldMatrix(IntegerMatrix const &m, Field f);
347
IntegerMatrix fieldMatrixToIntegerMatrix(FieldMatrix const &m);//MUST CONTAIN INTEGER ENTRIES ONLY
348
IntegerMatrix fieldMatrixToIntegerMatrixPrimitive(FieldMatrix const &m);//Scales to primitive vector
349
FieldVector integerVectorToFieldVector(IntegerVector const &v, Field f);
352
Computes a non-zero vector in the kernel of the given matrix. The
353
dimension of the kernel is assumed to be 1. The returned vector is primitive.
355
IntegerVector vectorInKernel(IntegerMatrix const &m);
358
* Using exact arithmetics this routine chooses a subset of the vectors in m, which generate the span of m.
360
IntegerVectorList subsetBasis(IntegerVectorList const &m);