My Project
Loading...
Searching...
No Matches
linearAlgebra.h
Go to the documentation of this file.
1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
4/** @file lineareAlgebra.h
5 *
6 * This file provides basic linear algebra functionality.
7 *
8 * ABSTRACT: This file provides basic algorithms from linear algebra over
9 * any SINGULAR-supported field.
10 * For the time being, the procedures defined in this file expect matrices
11 * containing objects of the SINGULAR type 'number'. This means that, when
12 * 'currentRing' represents some polynomial ring K[x_1, x_2, ..., x_n], then
13 * the entries of the matrices are 'numbers' representing elements of K (and
14 * NOT 'polys' in K[x_1, x_2, ..., x_n]).
15 * This restriction may become obsolete in the future.
16 *
17 * @author Frank Seelisch
18 *
19 *
20 **/
21/*****************************************************************************/
22
23#ifndef LINEAR_ALGEBRA_H
24#define LINEAR_ALGEBRA_H
25
26// include basic SINGULAR structures
27#include "kernel/structs.h"
28
29/**
30 * LU-decomposition of a given (m x n)-matrix.
31 *
32 * Given an (m x n) matrix A, the method computes its LU-decomposition,
33 * that is, it computes matrices P, L, and U such that<br>
34 * - P * A = L * U,<br>
35 * - P is an (m x m) permutation matrix, i.e., its row/columns form the
36 * standard basis of K^m,<br>
37 * - L is an (m x m) matrix in lower triangular form with all diagonal
38 * entries equal to 1,<br>
39 * - U is an (m x n) matrix in upper row echelon form.<br>
40 * From these conditions, it follows immediately that also A = P * L * U,
41 * since P is self-inverse.<br>
42 *
43 * The method will modify all argument matrices except aMat, so that
44 * they will finally contain the matrices P, L, and U as specified
45 * above.
46 **/
47void luDecomp(
48 const matrix aMat, /**< [in] the initial matrix A, */
49 matrix &pMat, /**< [out] the row permutation matrix P */
50 matrix &lMat, /**< [out] the lower triangular matrix L */
51 matrix &uMat, /**< [out] the upper row echelon matrix U */
52 const ring r= currRing /**< [in] current ring */
53 );
54
55/**
56 * Returns a pivot score for any non-zero matrix entry.
57 *
58 * The smaller the score the better will n serve as a pivot element
59 * in subsequent Gauss elimination steps.
60 *
61 * @return the pivot score of n
62 **/
63int pivotScore(
64 number n, /**< [in] a non-zero matrix entry */
65 const ring r= currRing /**< [in] current ring */
66 );
67
68/**
69 * Finds the best pivot element in some part of a given matrix.
70 *
71 * Given any matrix A with valid row indices r1..r2 and valid column
72 * indices c1..c2, this method finds the best pivot element for
73 * subsequent Gauss elimination steps in A[r1..r2, c1..c2]. "Best"
74 * here means best with respect to the implementation of the method
75 * 'pivotScore(number n)'.<br>
76 * In the case that all elements in A[r1..r2, c1..c2] are zero, the
77 * method returns false, otherwise true.
78 *
79 * @return false if all relevant matrix entries are zero, true otherwise
80 * @sa pivotScore(number n)
81 **/
82bool pivot(
83 const matrix aMat, /**< [in] any matrix with number entries */
84 const int r1, /**< [in] lower row index */
85 const int r2, /**< [in] upper row index */
86 const int c1, /**< [in] lower column index */
87 const int c2, /**< [in] upper column index */
88 int* bestR, /**< [out] address of row index of best
89 pivot element */
90 int* bestC, /**< [out] address of column index of
91 best pivot element */
92 const ring r= currRing /**< [in] current ring */
93 );
94
95/**
96 * Computes the inverse of a given (n x n)-matrix.
97 *
98 * This method expects an (n x n)-matrix, that is, it must have as many
99 * rows as columns. Inversion of the first argument is attempted via the
100 * LU-decomposition. There are two cases:<br>
101 * 1) The matrix is not invertible. Then the method returns false, and
102 * &iMat remains unchanged.<br>
103 * 2) The matrix is invertible. Then the method returns true, and the
104 * content of iMat is the inverse of aMat.
105 *
106 * @return true iff aMat is invertible, false otherwise
107 * @sa luInverseFromLUDecomp(const matrix pMat, const matrix lMat,
108 * const matrix uMat, matrix &iMat)
109 **/
110bool luInverse(
111 const matrix aMat, /**< [in] matrix to be inverted */
112 matrix &iMat, /**< [out] inverse of aMat if
113 invertible */
114 const ring r=currRing /**< [in] current ring */
115 );
116
117/**
118 * Computes the inverse of a given (n x n)-matrix in upper right
119 * triangular form.
120 *
121 * This method expects a quadratic matrix, that is, it must have as
122 * many rows as columns. Moreover, uMat[i, j] = 0, at least for all
123 * i > j, that is, u is in upper right triangular form.<br>
124 * If the argument diagonalIsOne is true, then we know additionally,
125 * that uMat[i, i] = 1, for all i. In this case uMat is invertible.
126 * Contrariwise, if diagonalIsOne is false, we do not know anything
127 * about the diagonal entries. (Note that they may still all be
128 * 1.)<br>
129 * In general, there are two cases:<br>
130 * 1) The matrix is not invertible. Then the method returns false,
131 * and &iMat remains unchanged.<br>
132 * 2) The matrix is invertible. Then the method returns true, and
133 * the content of iMat is the inverse of uMat.
134 *
135 * @return true iff uMat is invertible, false otherwise
136 **/
138 const matrix uMat, /**< [in] (n x n)-matrix in upper right
139 triangular form */
140 matrix &iMat, /**< [out] inverse of uMat if invertible */
141 bool diagonalIsOne,/**< [in] if true, then all diagonal
142 entries of uMat are 1 */
143 const ring r= currRing /**< [in] current ring */
144 );
145
146/**
147 * Computes the inverse of a given (n x n)-matrix in lower left
148 * triangular form.
149 *
150 * This method expects an (n x n)-matrix, that is, it must have as
151 * many rows as columns. Moreover, lMat[i,j] = 0, at least for all
152 * j > i, that ism lMat is in lower left triangular form.<br>
153 * If the argument diagonalIsOne is true, then we know additionally,
154 * that lMat[i, i] = 1, for all i. In this case lMat is invertible.
155 * Contrariwise, if diagonalIsOne is false, we do not know anything
156 * about the diagonal entries. (Note that they may still all be
157 * 1.)<br>
158 * In general, there are two cases:<br>
159 * 1) The matrix is not invertible. Then the method returns false,
160 * and &iMat remains unchanged.<br>
161 * 2) The matrix is invertible. Then the method returns true, and
162 * the content of iMat is the inverse of lMat.
163 *
164 * @return true iff lMat is invertible, false otherwise
165 **/
167 const matrix lMat, /**< [in] (n x n)-matrix in lower left
168 triangular form */
169 matrix &iMat, /**< [out] inverse of lMat if invertible */
170 bool diagonalIsOne /**< [in] if true, then all diagonal
171 entries of lMat are 1 */
172 );
173
174/**
175 * Computes the inverse of an (n x n)-matrix which is given by its LU-
176 * decomposition.
177 *
178 * With A denoting the matrix to be inverted, the method expects the
179 * LU-decomposition of A, that is, pMat * A = lMat * uMat, where
180 * the argument matrices have the appropriate properties; see method
181 * 'luDecomp(const matrix aMat, matrix &pMat, matrix &lMat,
182 * matrix &uMat)'.<br>
183 * Furthermore, uMat is expected to be an (n x n)-matrix. Then A^(-1)
184 * is computed according to A^(-1) = uMat^(-1) * lMat^(-1) * pMat,
185 * since pMat is self-inverse. This will work if and only if uMat is
186 * invertible, because lMat and pMat are in any case invertible.<br>
187 * In general, there are two cases:<br>
188 * 1) uMat and hence A is not invertible. Then the method returns
189 * false, and &iMat remains unchanged.<br>
190 * 2) uMat and hence A is invertible. Then the method returns true,
191 * and the content of iMat is the inverse of A.
192 *
193 * @return true if A is invertible, false otherwise
194 * @sa luInverse(const matrix aMat, matrix &iMat)
195 **/
197 const matrix pMat, /**< [in] permutation matrix of an LU-
198 decomposition */
199 const matrix lMat, /**< [in] lower left matrix of an LU-
200 decomposition */
201 const matrix uMat, /**< [in] upper right matrix of an LU-
202 decomposition */
203 matrix &iMat, /**< [out] inverse of A if invertible */
204 const ring r= currRing
205 );
206
207/**
208 * Computes the rank of a given (m x n)-matrix.
209 *
210 * The matrix may already be given in row echelon form, e.g., when
211 * the user has before called luDecomp and passes the upper triangular
212 * matrix U to luRank. In this case, the second argument can be set to
213 * true in order to pass this piece of information to the method.
214 * Otherwise, luDecomp will be called first to compute the matrix U.
215 * The rank is then read off the matrix U.
216 *
217 * @return the rank as an int
218 * @sa luDecomp(const matrix aMat, matrix &pMat, matrix &lMat, matrix &uMat)
219 **/
220int luRank(
221 const matrix aMat, /**< [in] input matrix */
222 const bool isRowEchelon,/**< [in] if true then aMat is assumed to be
223 already in row echelon form */
224 const ring r= currRing
225 );
226
227/**
228 * Solves the linear system A * x = b, where A is an (m x n)-matrix
229 * which is given by its LU-decomposition.
230 *
231 * The method expects the LU-decomposition of A, that is,
232 * pMat * A = lMat * uMat, where the argument matrices have the
233 * appropriate properties; see method
234 * 'luDecomp(const matrix aMat, matrix &pMat, matrix &lMat,
235 * matrix &uMat)'.<br>
236 * Instead of trying to invert A and return x = A^(-1)*b, this
237 * method
238 * 1) computes b' = pMat * b,
239 * 2) solves the simple system L * y = b', and then
240 * 3) solves the simple system U * x = y.
241 * Note that steps 1) and 2) will always work, as L is in any case
242 * invertible. Moreover, y is uniquely determined. Step 3) will only
243 * work if and only if y is in the column span of U. In that case,
244 * there may be infinitely many solutions.
245 * Thus, there are three cases:<br>
246 * 1) There is no solution. Then the method returns false, and &xVec
247 * as well as &H remain unchanged.<br>
248 * 2) There is a unique solution. Then the method returns true and
249 * H is the 1x1 matrix with zero entry.<br>
250 * 3) There are infinitely many solutions. Then the method returns
251 * true and some solution of the given original linear system.
252 * Furthermore, the columns of H span the solution space of the
253 * homogeneous linear system. The dimension of the solution
254 * space is then the number of columns of H.
255 *
256 * @return true if there is at least one solution, false otherwise
257 **/
259 const matrix pMat, /**< [in] permutation matrix of an LU-
260 decomposition */
261 const matrix lMat, /**< [in] lower left matrix of an LU-
262 decomposition */
263 const matrix uMat, /**< [in] upper right matrix of an LU-
264 decomposition */
265 const matrix bVec, /**< [in] right-hand side of the linear
266 system to be solved */
267 matrix &xVec, /**< [out] solution of A*x = b */
268 matrix &H /**< [out] matrix with columns spanning
269 homogeneous solution space */
270 );
271
272/**
273 * Solves the linear system A * x = b, where A is an (m x n)-matrix
274 * which is given by its LDU-decomposition.
275 *
276 * The method expects the LDU-decomposition of A, that is,
277 * pMat * A = lMat * dMat^(-1) * uMat, where the argument matrices have
278 * the appropriate properties; see method
279 * 'lduDecomp(const matrix aMat, matrix &pMat, matrix &lMat,
280 * matrix &dMat, matrix &uMat, poly &l, poly &u, poly &lTimesU)'.<br>
281 * Instead of trying to invert A and return x = A^(-1)*b, this
282 * method
283 * 1) computes b' = l * pMat * b,
284 * 2) solves the simple system L * y = b',
285 * 3) computes y' = u * dMat * y,
286 * 4) solves the simple system U * y'' = y',
287 * 5) computes x = 1/(lTimesU) * y''.
288 * Note that steps 1), 2) and 3) will always work, as L is in any case
289 * invertible. Moreover, y and thus y' are uniquely determined.
290 * Step 4) will only work if and only if y' is in the column span of U.
291 * In that case, there may be infinitely many solutions.
292 * In contrast to luSolveViaLUDecomp, this algorithm guarantees that
293 * all divisions which have to be performed in steps 2) and 4) will
294 * work without remainder. This is due to properties of the given LDU-
295 * decomposition. Only the final step 5) performs a division of a vector
296 * by a member of the ground field. (Suppose the ground field is Q or
297 * some rational function field. Then, if A contains only integers or
298 * polynomials, respectively, then steps 1) - 4) will keep all data
299 * integer or polynomial, respectively. This may speed up computations
300 * considerably.)
301 * For the outcome, there are three cases:<br>
302 * 1) There is no solution. Then the method returns false, and &xVec
303 * as well as &H remain unchanged.<br>
304 * 2) There is a unique solution. Then the method returns true and
305 * H is the 1x1 matrix with zero entry.<br>
306 * 3) There are infinitely many solutions. Then the method returns
307 * true and some solution of the given original linear system.
308 * Furthermore, the columns of H span the solution space of the
309 * homogeneous linear system. The dimension of the solution
310 * space is then the number of columns of H.
311 *
312 * @return true if there is at least one solution, false otherwise
313 **/
315 const matrix pMat, /**< [in] permutation matrix of an LDU-
316 decomposition */
317 const matrix lMat, /**< [in] lower left matrix of an LDU-
318 decomposition */
319 const matrix dMat, /**< [in] diagonal matrix of an LDU-
320 decomposition */
321 const matrix uMat, /**< [in] upper right matrix of an LDU-
322 decomposition */
323 const poly l, /**< [in] pivot product l of an LDU decomp. */
324 const poly u, /**< [in] pivot product u of an LDU decomp. */
325 const poly lTimesU, /**< [in] product of l and u */
326 const matrix bVec, /**< [in] right-hand side of the linear
327 system to be solved */
328 matrix &xVec, /**< [out] solution of A*x = b */
329 matrix &H /**< [out] matrix with columns spanning
330 homogeneous solution space */
331 );
332
333/**
334 * Creates a new matrix which is the (nxn) unit matrix, and returns true
335 * in case of success.
336 *
337 * The method will be successful whenever n >= 1. In this case and only then
338 * true will be returned and the new (nxn) unit matrix will reside inside
339 * the second argument.
340 *
341 * @return true iff the (nxn) unit matrix could be build
342 **/
343bool unitMatrix(
344 const int n, /**< [in] size of the matrix */
345 matrix &unitMat, /**< [out] the new (nxn) unit matrix */
346 const ring r= currRing /** [in] current ring */
347 );
348
349/**
350 * Creates a new matrix which is a submatrix of the first argument, and
351 * returns true in case of success.
352 *
353 * The method will be successful whenever rowIndex1 <= rowIndex2 and
354 * colIndex1 <= colIndex2. In this case and only then true will be
355 * returned and the last argument will afterwards contain a copy of the
356 * respective submatrix of the first argument.
357 * Note that this method may also be used to copy an entire matrix.
358 *
359 * @return true iff the submatrix could be build
360 **/
361bool subMatrix(
362 const matrix aMat, /**< [in] the input matrix */
363 const int rowIndex1, /**< [in] lower row index */
364 const int rowIndex2, /**< [in] higher row index */
365 const int colIndex1, /**< [in] lower column index */
366 const int colIndex2, /**< [in] higher column index */
367 matrix &subMat /**< [out] the new submatrix */
368 );
369
370/**
371 * Swaps two rows of a given matrix and thereby modifies it.
372 *
373 * The method expects two valid, distinct row indices, i.e. in
374 * [1..n], where n is the number of rows in aMat.
375 **/
376void swapRows(
377 int row1, /**< [in] index of first row to swap */
378 int row2, /**< [in] index of second row to swap */
379 matrix& aMat /**< [in/out] matrix subject to swapping */
380 );
381
382/**
383 * Swaps two columns of a given matrix and thereby modifies it.
384 *
385 * The method expects two valid, distinct column indices, i.e. in
386 * [1..n], where n is the number of columns in aMat.
387 **/
388void swapColumns(
389 int column1, /**< [in] index of first column to swap */
390 int column2, /**< [in] index of second column to swap */
391 matrix& aMat /**< [in/out] matrix subject to swapping */
392 );
393
394/**
395 * Creates a new matrix which contains the first argument in the top left
396 * corner, and the second in the bottom right.
397 *
398 * All other entries of the resulting matrix which will be created in the
399 * third argument, are zero.
400 **/
401void matrixBlock(
402 const matrix aMat, /**< [in] the top left input matrix */
403 const matrix bMat, /**< [in] the bottom right input matrix */
404 matrix &block /**< [out] the new block matrix */
405 );
406
407/**
408 * Computes the characteristic polynomial from a quadratic (2x2) matrix
409 * and returns true in case of success.
410 *
411 * The method will be successful whenever the input matrix is a (2x2) matrix.
412 * In this case, the resulting polynomial will be a univariate polynomial in
413 * the first ring variable of degree 2 and it will reside in the second
414 * argument.
415 * The method assumes that the matrix entries are all numbers, i.e., elements
416 * from the ground field/ring.
417 *
418 * @return true iff the input matrix was (2x2)
419 **/
420bool charPoly(
421 const matrix aMat, /**< [in] the input matrix */
422 poly &charPoly /**< [out] the characteristic polynomial */
423 );
424
425/**
426 * Computes the square root of a non-negative real number and returns
427 * it as a new number.
428 *
429 * The method assumes that the current ground field is the complex
430 * numbers, and that the argument has imaginary part zero.
431 * If the real part is negative, then false is returned, otherwise true.
432 * The square root will be computed in the last argument by Heron's
433 * iteration formula with the first argument as the starting value. The
434 * iteration will stop as soon as two successive values (in the resulting
435 * sequence) differ by no more than the given tolerance, which is assumed
436 * to be a positive real number.
437 *
438 * @return the square root of the non-negative argument number
439 **/
440bool realSqrt(
441 const number n, /**< [in] the input number */
442 const number tolerance, /**< [in] accuracy of iteration */
443 number &root /**< [out] the root of n */
444 );
445
446/**
447 * Computes the Hessenberg form of a given square matrix.
448 *
449 * The method assumes that all matrix entries are numbers coming from some
450 * subfield of the reals..
451 * Afterwards, the following conditions will hold:
452 * 1) hessenbergMat = pMat * aMat * pMat^{-1},
453 * 2) hessenbergMat is in Hessenberg form.
454 * During the algorithm, pMat will be constructed as the product of self-
455 * inverse matrices.
456 * The algorithm relies on computing square roots of floating point numbers.
457 * These will be combuted by Heron's iteration formula, with iteration
458 * stopping when two successive approximations of the square root differ by
459 * no more than the given tolerance, which is assumed to be a positive real
460 * number.
461 **/
462void hessenberg(
463 const matrix aMat, /**< [in] the square input matrix */
464 matrix &pMat, /**< [out] the transformation matrix */
465 matrix &hessenbergMat, /**< [out] the Hessenberg form of aMat */
466 const number tolerance, /**< [in] accuracy */
467 const ring r
468 );
469
470/**
471 * Returns all solutions of a quadratic polynomial relation with real
472 * coefficients.
473 *
474 * The method assumes that the polynomial is univariate in the first
475 * ring variable, and that the current ground field is the complex numbers.
476 * The polynomial has degree <= 2. Thus, there may be
477 * a) infinitely many zeros, when p == 0; then -1 is returned;
478 * b) no zero, when p == constant <> 0; then 0 is returned;
479 * c) one zero, when p is linear; then 1 is returned;
480 * d) one double zero; then 2 is returned;
481 * e) two distinct zeros; then 3 is returned;
482 * In the case e), s1 and s2 will contain the two distinct solutions.
483 * In cases c) and d) s1 will contain the single/double solution.
484 *
485 * @return the number of distinct zeros
486 **/
488 const poly p, /**< [in] the polynomial */
489 number &s1, /**< [out] first zero, if any */
490 number &s2, /**< [out] second zero, if any */
491 const number tolerance /**< [in] accuracy */
492 );
493
494
495/**
496 * Computes a factorization of a polynomial h(x, y) in K[[x]][y] up to a
497 * certain degree in x, whenever a factorization of h(0, y) is given.
498 *
499 * The algorithm is based on Hensel's lemma: Let h(x, y) denote a monic
500 * polynomial in y of degree m + n with coefficients in K[[x]]. Suppose there
501 * are two monic factors f_0(y) (of degree n) and g_0(y) of degree (m) such
502 * that h(0, y) = f_0(y) * g_0(y) and <f_0, g_0> = K[y]. Fix an integer d >= 0.
503 * Then there are monic polynomials in y with coefficients in K[[x]], namely
504 * f(x, y) of degree n and g(x, y) of degree m such that
505 * h(x, y) = f(x, y) * g(x, y) modulo <x^(d+1)> (*).
506 *
507 * This implementation expects h, f0, g0, and d as described and computes the
508 * factors f and g. Effectively, h will be given as an element of K[x, y] since
509 * all terms of h with x-degree larger than d can be ignored due to (*).
510 * The method expects the ground ring to contain at least two variables; then
511 * x is the first ring variable, specified by the input index xIndex, and y the
512 * second one, specified by yIndex.
513 *
514 * This code was placed here since the algorithm works by successively solving
515 * d linear equation systems. It is hence an application of other methods
516 * defined in this h-file and its corresponding cc-file.
517 *
518 **/
519void henselFactors(
520 const int xIndex, /**< [in] index of x in {1, ..., nvars(basering)} */
521 const int yIndex, /**< [in] index of y in {1, ..., nvars(basering)} */
522 const poly h, /**< [in] the polynomial h(x, y) as above */
523 const poly f0, /**< [in] the first univariate factor of h(0, y) */
524 const poly g0, /**< [in] the second univariate factor of h(0, y) */
525 const int d, /**< [in] the degree bound, d >= 0 */
526 poly &f, /**< [out] the first factor of h(x, y) */
527 poly &g /**< [out] the second factor of h(x, y) */
528 );
529
530/**
531 * LU-decomposition of a given (m x n)-matrix with performing only those
532 * divisions that yield zero remainders.
533 *
534 * Given an (m x n) matrix A, the method computes its LDU-decomposition,
535 * that is, it computes matrices P, L, D, and U such that<br>
536 * - P * A = L * D^(-1) * U,<br>
537 * - P is an (m x m) permutation matrix, i.e., its row/columns form the
538 * standard basis of K^m,<br>
539 * - L is an (m x m) matrix in lower triangular form of full rank,<br>
540 * - D is an (m x m) diagonal matrix with full rank, and<br>
541 * - U is an (m x n) matrix in upper row echelon form.<br>
542 * From these conditions, it follows immediately that also
543 * A = P * L * D^(-1) * U, since P is self-inverse.<br>
544 *
545 * In contrast to luDecomp, this method only performs those divisions which
546 * yield zero remainders. Hence, when e.g. computing over a rational function
547 * field and starting with polynomial entries only (or over Q and starting
548 * with integer entries only), then any newly computed matrix entry will again
549 * be polynomial (or integer).
550 *
551 * The method will modify all argument matrices except aMat, so that
552 * they will finally contain the matrices P, L, D, and U as specified
553 * above. Moreover, in order to further speed up subsequent calls of
554 * luSolveViaLDUDecomp, the two denominators and their product will also be
555 * returned.
556 **/
557void lduDecomp(
558 const matrix aMat, /**< [in] the initial matrix A, */
559 matrix &pMat, /**< [out] the row permutation matrix P */
560 matrix &lMat, /**< [out] the lower triangular matrix L */
561 matrix &dMat, /**< [out] the diagonal matrix D */
562 matrix &uMat, /**< [out] the upper row echelon matrix U */
563 poly &l, /**< [out] the product of pivots of L */
564 poly &u, /**< [out] the product of pivots of U */
565 poly &lTimesU /**< [out] the product of l and u */
566 );
567
568/* helper for qrDoubleShift */
569bool qrDS( const int n, matrix* queue, int& queueL, number* eigenValues,
570 int& eigenValuesL, const number tol1, const number tol2, const ring R);
571
572/**
573 * Tries to find the number n in the array nn[0..nnLength-1].
574 **/
575int similar(
576 const number* nn, /**< [in] array of numbers to look-up */
577 const int nnLength, /**< [in] length of nn */
578 const number n, /**< [in] number to loop-up in nn */
579 const number tolerance /**< [in] tolerance for comparison */
580 );
581#endif
582/* LINEAR_ALGEBRA_H */
int l
Definition: cfEzgcd.cc:100
int p
Definition: cfModGcd.cc:4078
g
Definition: cfModGcd.cc:4090
FILE * f
Definition: checklibs.c:9
CanonicalForm H
Definition: facAbsFact.cc:60
STATIC_VAR Poly * h
Definition: janet.cc:971
bool luSolveViaLDUDecomp(const matrix pMat, const matrix lMat, const matrix dMat, const matrix uMat, const poly l, const poly u, const poly lTimesU, const matrix bVec, matrix &xVec, matrix &H)
Solves the linear system A * x = b, where A is an (m x n)-matrix which is given by its LDU-decomposit...
bool luInverseFromLUDecomp(const matrix pMat, const matrix lMat, const matrix uMat, matrix &iMat, const ring r=currRing)
Computes the inverse of an (n x n)-matrix which is given by its LU- decomposition.
void luDecomp(const matrix aMat, matrix &pMat, matrix &lMat, matrix &uMat, const ring r=currRing)
LU-decomposition of a given (m x n)-matrix.
bool realSqrt(const number n, const number tolerance, number &root)
Computes the square root of a non-negative real number and returns it as a new number.
void swapColumns(int column1, int column2, matrix &aMat)
Swaps two columns of a given matrix and thereby modifies it.
bool luInverse(const matrix aMat, matrix &iMat, const ring r=currRing)
Computes the inverse of a given (n x n)-matrix.
void hessenberg(const matrix aMat, matrix &pMat, matrix &hessenbergMat, const number tolerance, const ring r)
Computes the Hessenberg form of a given square matrix.
int quadraticSolve(const poly p, number &s1, number &s2, const number tolerance)
Returns all solutions of a quadratic polynomial relation with real coefficients.
bool unitMatrix(const int n, matrix &unitMat, const ring r=currRing)
Creates a new matrix which is the (nxn) unit matrix, and returns true in case of success.
void matrixBlock(const matrix aMat, const matrix bMat, matrix &block)
Creates a new matrix which contains the first argument in the top left corner, and the second in the ...
void henselFactors(const int xIndex, const int yIndex, const poly h, const poly f0, const poly g0, const int d, poly &f, poly &g)
Computes a factorization of a polynomial h(x, y) in K[[x]][y] up to a certain degree in x,...
bool pivot(const matrix aMat, const int r1, const int r2, const int c1, const int c2, int *bestR, int *bestC, const ring r=currRing)
Finds the best pivot element in some part of a given matrix.
int similar(const number *nn, const int nnLength, const number n, const number tolerance)
Tries to find the number n in the array nn[0..nnLength-1].
bool qrDS(const int n, matrix *queue, int &queueL, number *eigenValues, int &eigenValuesL, const number tol1, const number tol2, const ring R)
bool charPoly(const matrix aMat, poly &charPoly)
Computes the characteristic polynomial from a quadratic (2x2) matrix and returns true in case of succ...
void lduDecomp(const matrix aMat, matrix &pMat, matrix &lMat, matrix &dMat, matrix &uMat, poly &l, poly &u, poly &lTimesU)
LU-decomposition of a given (m x n)-matrix with performing only those divisions that yield zero remai...
void swapRows(int row1, int row2, matrix &aMat)
Swaps two rows of a given matrix and thereby modifies it.
int luRank(const matrix aMat, const bool isRowEchelon, const ring r=currRing)
Computes the rank of a given (m x n)-matrix.
bool lowerLeftTriangleInverse(const matrix lMat, matrix &iMat, bool diagonalIsOne)
Computes the inverse of a given (n x n)-matrix in lower left triangular form.
bool subMatrix(const matrix aMat, const int rowIndex1, const int rowIndex2, const int colIndex1, const int colIndex2, matrix &subMat)
Creates a new matrix which is a submatrix of the first argument, and returns true in case of success.
int pivotScore(number n, const ring r=currRing)
Returns a pivot score for any non-zero matrix entry.
bool upperRightTriangleInverse(const matrix uMat, matrix &iMat, bool diagonalIsOne, const ring r=currRing)
Computes the inverse of a given (n x n)-matrix in upper right triangular form.
bool luSolveViaLUDecomp(const matrix pMat, const matrix lMat, const matrix uMat, const matrix bVec, matrix &xVec, matrix &H)
Solves the linear system A * x = b, where A is an (m x n)-matrix which is given by its LU-decompositi...
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:13
#define block
Definition: scanner.cc:646
#define R
Definition: sirandom.c:27