| blob | fdf9185e961a0dcff662280e6df3f12e22d15445 |
1 @menu
2 * Introduction to linearalgebra::
3 * Functions and Variables for linearalgebra::
4 @end menu
6 @c -----------------------------------------------------------------------------
7 @node Introduction to linearalgebra, Functions and Variables for linearalgebra
8 @section Introduction to linearalgebra
9 @c -----------------------------------------------------------------------------
11 @code{linearalgebra} is a collection of functions for linear algebra.
13 Example:
15 @c ===beg===
16 @c M : matrix ([1, 2], [1, 2]);
17 @c nullspace (M);
18 @c columnspace (M);
19 @c ptriangularize (M - z*ident(2), z);
20 @c M : matrix ([1, 2, 3], [4, 5, 6], [7, 8, 9]) - z*ident(3);
21 @c MM : ptriangularize (M, z);
22 @c algebraic : true;
23 @c tellrat (MM [3, 3]);
24 @c MM : ratsimp (MM);
25 @c nullspace (MM);
26 @c M : matrix ([1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12],
27 @c [13, 14, 15, 16]);
28 @c columnspace (M);
29 @c apply ('orthogonal_complement, args (nullspace (transpose (M))));
30 @c ===end===
31 @example
32 (%i1) M : matrix ([1, 2], [1, 2]);
33 [ 1 2 ]
34 (%o1) [ ]
35 [ 1 2 ]
36 (%i2) nullspace (M);
37 [ 1 ]
38 [ ]
39 (%o2) span([ 1 ])
40 [ - - ]
41 [ 2 ]
42 (%i3) columnspace (M);
43 [ 1 ]
44 (%o3) span([ ])
45 [ 1 ]
46 (%i4) ptriangularize (M - z*ident(2), z);
47 [ 1 2 - z ]
48 (%o4) [ ]
49 [ 2 ]
50 [ 0 3 z - z ]
51 (%i5) M : matrix ([1, 2, 3], [4, 5, 6], [7, 8, 9]) - z*ident(3);
52 [ 1 - z 2 3 ]
53 [ ]
54 (%o5) [ 4 5 - z 6 ]
55 [ ]
56 [ 7 8 9 - z ]
57 (%i6) MM : ptriangularize (M, z);
58 [ 4 5 - z 6 ]
59 [ ]
60 [ 2 ]
61 [ 66 z 102 z 132 ]
62 [ 0 -- - -- + ----- + --- ]
63 (%o6) [ 49 7 49 49 ]
64 [ ]
65 [ 3 2 ]
66 [ 49 z 245 z 147 z ]
67 [ 0 0 ----- - ------ - ----- ]
68 [ 264 88 44 ]
69 (%i7) algebraic : true;
70 (%o7) true
71 (%i8) tellrat (MM [3, 3]);
72 3 2
73 (%o8) [z - 15 z - 18 z]
74 (%i9) MM : ratsimp (MM);
75 [ 4 5 - z 6 ]
76 [ ]
77 [ 2 ]
78 (%o9) [ 66 7 z - 102 z - 132 ]
79 [ 0 -- - ------------------ ]
80 [ 49 49 ]
81 [ ]
82 [ 0 0 0 ]
83 (%i10) nullspace (MM);
84 [ 1 ]
85 [ ]
86 [ 2 ]
87 [ z - 14 z - 16 ]
88 [ -------------- ]
89 (%o10) span([ 8 ])
90 [ ]
91 [ 2 ]
92 [ z - 18 z - 12 ]
93 [ - -------------- ]
94 [ 12 ]
95 (%i11) M : matrix ([1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12],
96 [13, 14, 15, 16]);
97 [ 1 2 3 4 ]
98 [ ]
99 [ 5 6 7 8 ]
100 (%o11) [ ]
101 [ 9 10 11 12 ]
102 [ ]
103 [ 13 14 15 16 ]
104 (%i12) columnspace (M);
105 [ 1 ] [ 2 ]
106 [ ] [ ]
107 [ 5 ] [ 6 ]
108 (%o12) span([ ], [ ])
109 [ 9 ] [ 10 ]
110 [ ] [ ]
111 [ 13 ] [ 14 ]
112 (%i13) apply ('orthogonal_complement, args (nullspace (transpose (M))));
113 [ 0 ] [ 1 ]
114 [ ] [ ]
115 [ 1 ] [ 0 ]
116 (%o13) span([ ], [ ])
117 [ 2 ] [ - 1 ]
118 [ ] [ ]
119 [ 3 ] [ - 2 ]
120 @end example
122 @opencatbox{Categories:}
123 @category{Linear algebra}
124 @category{Share packages}
125 @category{Package linearalgebra}
126 @closecatbox
128 @c -----------------------------------------------------------------------------
129 @need 800
130 @node Functions and Variables for linearalgebra, , Introduction to linearalgebra
131 @section Functions and Variables for linearalgebra
132 @c -----------------------------------------------------------------------------
134 @c -----------------------------------------------------------------------------
135 @anchor{addmatrices}
136 @deffn {Function} addmatrices (@var{f}, @var{M_1}, @dots{}, @var{M_n})
138 @c REWORD -- THE RESULT IS NOT GENERALLY THE SUM OF M_1, ..., M_N
139 Using the function @var{f} as the addition function, return the sum of the
140 matrices @var{M_1}, @dots{}, @var{M_n}. The function @var{f} must accept any
141 number of arguments (a Maxima nary function).
143 Examples:
145 @c ===beg===
146 @c m1 : matrix([1,2],[3,4])$
147 @c m2 : matrix([7,8],[9,10])$
148 @c addmatrices('max,m1,m2);
149 @c addmatrices('max,m1,m2,5*m1);
150 @c ===end===
151 @example
152 (%i1) m1 : matrix([1,2],[3,4])$
153 (%i2) m2 : matrix([7,8],[9,10])$
154 (%i3) addmatrices('max,m1,m2);
155 (%o3) matrix([7,8],[9,10])
156 (%i4) addmatrices('max,m1,m2,5*m1);
157 (%o4) matrix([7,10],[15,20])
158 @end example
160 @opencatbox{Categories:}
161 @category{Package linearalgebra}
162 @closecatbox
163 @end deffn
165 @c -----------------------------------------------------------------------------
166 @anchor{blockmatrixp}
167 @deffn {Function} blockmatrixp (@var{M})
169 Return true if and only if @var{M} is a matrix and every entry of
170 @var{M} is a matrix.
172 @opencatbox{Categories:}
173 @category{Package linearalgebra}
174 @category{Predicate functions}
175 @closecatbox
176 @end deffn
178 @c -----------------------------------------------------------------------------
179 @anchor{columnop}
180 @deffn {Function} columnop (@var{M}, @var{i}, @var{j}, @var{theta})
182 If @var{M} is a matrix, return the matrix that results from doing the column
183 operation @code{C_i <- C_i - @var{theta} * C_j}. If @var{M} doesn't have a row
184 @var{i} or @var{j}, signal an error.
186 @opencatbox{Categories:}
187 @category{Package linearalgebra}
188 @closecatbox
189 @end deffn
191 @c -----------------------------------------------------------------------------
192 @anchor{columnswap}
193 @deffn {Function} columnswap (@var{M}, @var{i}, @var{j})
195 If @var{M} is a matrix, swap columns @var{i} and @var{j}. If @var{M} doesn't
196 have a column @var{i} or @var{j}, signal an error.
198 @opencatbox{Categories:}
199 @category{Package linearalgebra}
200 @closecatbox
201 @end deffn
203 @c -----------------------------------------------------------------------------
204 @anchor{columnspace}
205 @deffn {Function} columnspace (@var{M})
207 If @var{M} is a matrix, return @code{span (v_1, ..., v_n)}, where the set
208 @code{@{v_1, ..., v_n@}} is a basis for the column space of @var{M}. The span
209 of the empty set is @code{@{0@}}. Thus, when the column space has only
210 one member, return @code{span ()}.
212 @opencatbox{Categories:}
213 @category{Package linearalgebra}
214 @closecatbox
215 @end deffn
217 @c -----------------------------------------------------------------------------
218 @anchor{cholesky}
219 @deffn {Function} cholesky @
220 @fname{cholesky} (@var{M}) @
221 @fname{cholesky} (@var{M}, @var{field})
223 Return the Cholesky factorization of the matrix selfadjoint (or hermitian)
224 matrix @var{M}. The second argument defaults to 'generalring.' For a
225 description of the possible values for @var{field}, see @code{lu_factor}.
227 @opencatbox{Categories:}
228 @category{Matrix decompositions}
229 @category{Package linearalgebra}
230 @closecatbox
231 @end deffn
233 @c -----------------------------------------------------------------------------
234 @anchor{ctranspose}
235 @deffn {Function} ctranspose (@var{M})
237 Return the complex conjugate transpose of the matrix @var{M}. The function
238 @code{ctranspose} uses @code{matrix_element_transpose} to transpose each matrix
239 element.
241 @opencatbox{Categories:}
242 @category{Package linearalgebra}
243 @closecatbox
244 @end deffn
246 @c -----------------------------------------------------------------------------
247 @anchor{diag_matrix}
248 @deffn {Function} diag_matrix (@var{d_1}, @var{d_2}, @dots{}, @var{d_n})
250 Return a diagonal matrix with diagonal entries @var{d_1}, @var{d_2}, @dots{},
251 @var{d_n}. When the diagonal entries are matrices, the zero entries of the
252 returned matrix are zero matrices of the appropriate size; for example:
254 @c ===beg===
255 @c diag_matrix(diag_matrix(1,2),diag_matrix(3,4));
256 @c diag_matrix(p,q);
257 @c ===end===
258 @example
259 (%i1) diag_matrix(diag_matrix(1,2),diag_matrix(3,4));
261 [ [ 1 0 ] [ 0 0 ] ]
262 [ [ ] [ ] ]
263 [ [ 0 2 ] [ 0 0 ] ]
264 (%o1) [ ]
265 [ [ 0 0 ] [ 3 0 ] ]
266 [ [ ] [ ] ]
267 [ [ 0 0 ] [ 0 4 ] ]
268 (%i2) diag_matrix(p,q);
270 [ p 0 ]
271 (%o2) [ ]
272 [ 0 q ]
273 @end example
275 @opencatbox{Categories:}
276 @category{Package linearalgebra}
277 @closecatbox
278 @end deffn
280 @c -----------------------------------------------------------------------------
281 @anchor{dotproduct}
282 @deffn {Function} dotproduct (@var{u}, @var{v})
284 Return the dotproduct of vectors @var{u} and @var{v}. This is the same as
285 @code{conjugate (transpose (@var{u})) . @var{v}}. The arguments @var{u} and
286 @var{v} must be column vectors.
288 @opencatbox{Categories:}
289 @category{Package linearalgebra}
290 @closecatbox
291 @end deffn
293 @c -----------------------------------------------------------------------------
294 @anchor{eigens_by_jacobi}
295 @deffn {Function} eigens_by_jacobi @
296 @fname{eigens_by_jacobi} (@var{A}) @
297 @fname{eigens_by_jacobi} (@var{A}, @var{field_type})
299 Computes the eigenvalues and eigenvectors of @var{A} by the method of Jacobi
300 rotations. @var{A} must be a symmetric matrix (but it need not be positive
301 definite nor positive semidefinite). @var{field_type} indicates the
302 computational field, either @code{floatfield} or @code{bigfloatfield}.
303 If @var{field_type} is not specified, it defaults to @code{floatfield}.
305 The elements of @var{A} must be numbers or expressions which evaluate to numbers
306 via @code{float} or @code{bfloat} (depending on @var{field_type}).
308 Examples:
310 @c ===beg===
311 @c S : matrix ([1/sqrt(2), 1/sqrt(2)], [- 1/sqrt(2), 1/sqrt(2)]);
312 @c L : matrix ([sqrt(3), 0], [0, sqrt(5)]);
313 @c M : S . L . transpose (S);
314 @c eigens_by_jacobi (M);
315 @c float ([[sqrt(3), sqrt(5)], S]);
316 @c eigens_by_jacobi (M, bigfloatfield);
317 @c ===end===
318 @example
319 (%i1) S: matrix([1/sqrt(2), 1/sqrt(2)],[-1/sqrt(2), 1/sqrt(2)]);
320 [ 1 1 ]
321 [ ------- ------- ]
322 [ sqrt(2) sqrt(2) ]
323 (%o1) [ ]
324 [ 1 1 ]
325 [ - ------- ------- ]
326 [ sqrt(2) sqrt(2) ]
327 (%i2) L : matrix ([sqrt(3), 0], [0, sqrt(5)]);
328 [ sqrt(3) 0 ]
329 (%o2) [ ]
330 [ 0 sqrt(5) ]
331 (%i3) M : S . L . transpose (S);
332 [ sqrt(5) sqrt(3) sqrt(5) sqrt(3) ]
333 [ ------- + ------- ------- - ------- ]
334 [ 2 2 2 2 ]
335 (%o3) [ ]
336 [ sqrt(5) sqrt(3) sqrt(5) sqrt(3) ]
337 [ ------- - ------- ------- + ------- ]
338 [ 2 2 2 2 ]
339 (%i4) eigens_by_jacobi (M);
340 The largest percent change was 0.1454972243679
341 The largest percent change was 0.0
342 number of sweeps: 2
343 number of rotations: 1
344 (%o4) [[1.732050807568877, 2.23606797749979],
345 [ 0.70710678118655 0.70710678118655 ]
346 [ ]]
347 [ - 0.70710678118655 0.70710678118655 ]
348 (%i5) float ([[sqrt(3), sqrt(5)], S]);
349 (%o5) [[1.732050807568877, 2.23606797749979],
350 [ 0.70710678118655 0.70710678118655 ]
351 [ ]]
352 [ - 0.70710678118655 0.70710678118655 ]
353 (%i6) eigens_by_jacobi (M, bigfloatfield);
354 The largest percent change was 1.454972243679028b-1
355 The largest percent change was 0.0b0
356 number of sweeps: 2
357 number of rotations: 1
358 (%o6) [[1.732050807568877b0, 2.23606797749979b0],
359 [ 7.071067811865475b-1 7.071067811865475b-1 ]
360 [ ]]
361 [ - 7.071067811865475b-1 7.071067811865475b-1 ]
362 @end example
364 @opencatbox{Categories:}
365 @category{Matrix decompositions}
366 @category{Package linearalgebra}
367 @closecatbox
368 @end deffn
370 @c -----------------------------------------------------------------------------
371 @anchor{get_lu_factors}
372 @deffn {Function} get_lu_factors (@var{x})
374 When @code{@var{x} = lu_factor (@var{A})}, then @code{get_lu_factors} returns a
375 list of the form @code{[P, L, U]}, where @var{P} is a permutation matrix,
376 @var{L} is lower triangular with ones on the diagonal, and @var{U} is upper
377 triangular, and @code{@var{A} = @var{P} @var{L} @var{U}}.
379 @opencatbox{Categories:}
380 @category{Package linearalgebra}
381 @closecatbox
382 @end deffn
384 @c -----------------------------------------------------------------------------
385 @anchor{hankel}
386 @deffn {Function} hankel @
387 @fname{hankel} (@var{col}) @
388 @fname{hankel} (@var{col}, @var{row})
390 Return a Hankel matrix @var{H}. The first column of @var{H} is @var{col};
391 except for the first entry, the last row of @var{H} is @var{row}. The
392 default for @var{row} is the zero vector with the same length as @var{col}.
394 @opencatbox{Categories:}
395 @category{Package linearalgebra}
396 @closecatbox
397 @end deffn
399 @c -----------------------------------------------------------------------------
400 @anchor{hessian}
401 @deffn {Function} hessian (@var{f}, @var{x})
403 Returns the Hessian matrix of @var{f} with respect to the list of variables
404 @var{x}. The @code{(i, j)}-th element of the Hessian matrix is
405 @code{diff(@var{f}, @var{x}[i], 1, @var{x}[j], 1)}.
407 Examples:
409 @c ===beg===
410 @c hessian (x * sin (y), [x, y]);
411 @c depends (F, [a, b]);
412 @c hessian (F, [a, b]);
413 @c ===end===
414 @example
415 (%i1) hessian (x * sin (y), [x, y]);
416 [ 0 cos(y) ]
417 (%o1) [ ]
418 [ cos(y) - x sin(y) ]
419 (%i2) depends (F, [a, b]);
420 (%o2) [F(a, b)]
421 (%i3) hessian (F, [a, b]);
422 [ 2 2 ]
423 [ d F d F ]
424 [ --- ----- ]
425 [ 2 da db ]
426 [ da ]
427 (%o3) [ ]
428 [ 2 2 ]
429 [ d F d F ]
430 [ ----- --- ]
431 [ da db 2 ]
432 [ db ]
433 @end example
435 @opencatbox{Categories:}
436 @category{Differential calculus}
437 @category{Package linearalgebra}
438 @closecatbox
439 @end deffn
441 @c -----------------------------------------------------------------------------
442 @anchor{hilbert_matrix}
443 @deffn {Function} hilbert_matrix (@var{n})
445 Return the @var{n} by @var{n} Hilbert matrix. When @var{n} isn't a positive
446 integer, signal an error.
448 @opencatbox{Categories:}
449 @category{Package linearalgebra}
450 @closecatbox
451 @end deffn
453 @c -----------------------------------------------------------------------------
454 @anchor{identfor}
455 @deffn {Function} identfor @
456 @fname{identfor} (@var{M}) @
457 @fname{identfor} (@var{M}, @var{fld})
459 Return an identity matrix that has the same shape as the matrix
460 @var{M}. The diagonal entries of the identity matrix are the
461 multiplicative identity of the field @var{fld}; the default for
462 @var{fld} is @var{generalring}.
464 The first argument @var{M} should be a square matrix or a non-matrix. When
465 @var{M} is a matrix, each entry of @var{M} can be a square matrix -- thus
466 @var{M} can be a blocked Maxima matrix. The matrix can be blocked to any
467 (finite) depth.
469 See also @mref{zerofor}
471 @opencatbox{Categories:}
472 @category{Package linearalgebra}
473 @closecatbox
474 @end deffn
476 @c -----------------------------------------------------------------------------
477 @anchor{invert_by_lu}
478 @deffn {Function} invert_by_lu (@var{M}, @var{(rng generalring)})
480 Invert a matrix @var{M} by using the LU factorization. The LU factorization
481 is done using the ring @var{rng}.
483 @opencatbox{Categories:}
484 @category{Package linearalgebra}
485 @closecatbox
486 @end deffn
488 @c -----------------------------------------------------------------------------
489 @anchor{jacobian}
490 @deffn {Function} jacobian (@var{f}, @var{x})
492 Returns the Jacobian matrix of the list of functions @var{f} with respect to
493 the list of variables @var{x}. The @code{(i, j)}-th element of the Jacobian
494 matrix is @code{diff(@var{f}[i], @var{x}[j])}.
496 Examples:
498 @c ===beg===
499 @c jacobian ([sin (u - v), sin (u * v)], [u, v]);
500 @c depends ([F, G], [y, z]);
501 @c jacobian ([F, G], [y, z]);
502 @c ===end===
503 @example
504 (%i1) jacobian ([sin (u - v), sin (u * v)], [u, v]);
505 [ cos(v - u) - cos(v - u) ]
506 (%o1) [ ]
507 [ v cos(u v) u cos(u v) ]
508 (%i2) depends ([F, G], [y, z]);
509 (%o2) [F(y, z), G(y, z)]
510 (%i3) jacobian ([F, G], [y, z]);
511 [ dF dF ]
512 [ -- -- ]
513 [ dy dz ]
514 (%o3) [ ]
515 [ dG dG ]
516 [ -- -- ]
517 [ dy dz ]
518 @end example
520 @opencatbox{Categories:}
521 @category{Differential calculus}
522 @category{Package linearalgebra}
523 @closecatbox
524 @end deffn
526 @c -----------------------------------------------------------------------------
527 @anchor{kronecker_product}
528 @deffn {Function} kronecker_product (@var{A}, @var{B})
530 Return the Kronecker product of the matrices @var{A} and @var{B}.
532 @opencatbox{Categories:}
533 @category{Package linearalgebra}
534 @closecatbox
535 @end deffn
537 @c -----------------------------------------------------------------------------
538 @anchor{listp}
539 @deffn {Function} listp @
540 @fname{listp} (@var{e}, @var{p}) @
541 @fname{listp} (@var{e})
543 Given an optional argument @var{p}, return @code{true} if @var{e} is a Maxima
544 list and @var{p} evaluates to @code{true} for every list element. When
545 @code{listp} is not given the optional argument, return @code{true} if @var{e}
546 is a Maxima list. In all other cases, return @code{false}.
548 @opencatbox{Categories:}
549 @category{Package linearalgebra}
550 @category{Predicate functions}
551 @closecatbox
552 @end deffn
554 @c -----------------------------------------------------------------------------
555 @anchor{locate_matrix_entry}
556 @deffn {Function} locate_matrix_entry (@var{M}, @var{r_1}, @var{c_1}, @var{r_2}, @var{c_2}, @var{f}, @var{rel})
558 The first argument must be a matrix; the arguments
559 @var{r_1} through @var{c_2} determine a sub-matrix of @var{M} that consists of
560 rows @var{r_1} through @var{r_2} and columns @var{c_1} through @var{c_2}.
562 Find an entry in the sub-matrix @var{M} that satisfies some property.
563 Three cases:
565 (1) @code{@var{rel} = 'bool} and @var{f} a predicate:
567 Scan the sub-matrix from left to right then top to bottom,
568 and return the index of the first entry that satisfies the
569 predicate @var{f}. If no matrix entry satisfies @var{f}, return @code{false}.
571 (2) @code{@var{rel} = 'max} and @var{f} real-valued:
573 Scan the sub-matrix looking for an entry that maximizes @var{f}.
574 Return the index of a maximizing entry.
576 (3) @code{@var{rel} = 'min} and @var{f} real-valued:
578 Scan the sub-matrix looking for an entry that minimizes @var{f}.
579 Return the index of a minimizing entry.
581 @opencatbox{Categories:}
582 @category{Package linearalgebra}
583 @closecatbox
584 @end deffn
586 @c -----------------------------------------------------------------------------
587 @anchor{lu_backsub}
588 @deffn {Function} lu_backsub (@var{M}, @var{b})
590 When @code{@var{M} = lu_factor (@var{A}, @var{field})},
591 then @code{lu_backsub (@var{M}, @var{b})} solves the linear
592 system @code{@var{A} @var{x} = @var{b}}.
594 The @var{n} by @var{m} matrix @code{@var{b}}, with @var{n} the number of
595 rows of the matrix @code{@var{A}}, contains one right hand side per column. If
596 there is only one right hand side then @code{@var{b}} must be a @var{n} by 1
597 matrix.
599 Each column of the matrix @code{@var{x}=lu_backsub (@var{M}, @var{b})} is the
600 solution corresponding to the respective column of @code{@var{b}}.
602 Examples:
604 @c ===beg===
605 @c A : matrix ([1 - z, 3], [3, 8 - z]);
606 @c M : lu_factor (A,generalring);
607 @c b : matrix([a],[c]);
608 @c x : lu_backsub(M,b);
609 @c ratsimp(A . x - b);
610 @c B : matrix([a,d],[c,f]);
611 @c x : lu_backsub(M,B);
612 @c ratsimp(A . x - B);
613 @c ===end===
614 @example
615 (%i1) A : matrix ([1 - z, 3], [3, 8 - z]);
616 [ 1 - z 3 ]
617 (%o1) [ ]
618 [ 3 8 - z ]
619 (%i2) M : lu_factor (A,generalring);
620 [ 1 - z 3 ]
621 [ ]
622 (%o2) [[ 3 9 ], [1, 2], generalring]
623 [ ----- (- z) - ----- + 8 ]
624 [ 1 - z 1 - z ]
625 (%i3) b : matrix([a],[c]);
626 [ a ]
627 (%o3) [ ]
628 [ c ]
629 (%i4) x : lu_backsub(M,b);
630 [ 3 a ]
631 [ 3 (c - -----) ]
632 [ 1 - z ]
633 [ a - ----------------- ]
634 [ 9 ]
635 [ (- z) - ----- + 8 ]
636 [ 1 - z ]
637 [ --------------------- ]
638 (%o4) [ 1 - z ]
639 [ ]
640 [ 3 a ]
641 [ c - ----- ]
642 [ 1 - z ]
643 [ ----------------- ]
644 [ 9 ]
645 [ (- z) - ----- + 8 ]
646 [ 1 - z ]
647 (%i5) ratsimp(A . x - b);
648 [ 0 ]
649 (%o5) [ ]
650 [ 0 ]
651 (%i6) B : matrix([a,d],[c,f]);
652 [ a d ]
653 (%o6) [ ]
654 [ c f ]
655 (%i7) x : lu_backsub(M,B);
656 [ 3 a 3 d ]
657 [ 3 (c - -----) 3 (f - -----) ]
658 [ 1 - z 1 - z ]
659 [ a - ----------------- d - ----------------- ]
660 [ 9 9 ]
661 [ (- z) - ----- + 8 (- z) - ----- + 8 ]
662 [ 1 - z 1 - z ]
663 [ --------------------- --------------------- ]
664 (%o7) [ 1 - z 1 - z ]
665 [ ]
666 [ 3 a 3 d ]
667 [ c - ----- f - ----- ]
668 [ 1 - z 1 - z ]
669 [ ----------------- ----------------- ]
670 [ 9 9 ]
671 [ (- z) - ----- + 8 (- z) - ----- + 8 ]
672 [ 1 - z 1 - z ]
673 (%i8) ratsimp(A . x - B);
674 [ 0 0 ]
675 (%o8) [ ]
676 [ 0 0 ]
677 @end example
679 @opencatbox{Categories:}
680 @category{Package linearalgebra}
681 @closecatbox
682 @end deffn
684 @c -----------------------------------------------------------------------------
685 @anchor{lu_factor}
686 @deffn {Function} lu_factor (@var{M}, @var{field})
688 Return a list of the form @code{[@var{LU}, @var{perm}, @var{fld}]}, or
689 @code{[@var{LU}, @var{perm}, @var{fld}, @var{lower-cnd} @var{upper-cnd}]}, where
691 (1) The matrix @var{LU} contains the factorization of @var{M} in a packed form.
692 Packed form means three things: First, the rows of @var{LU} are permuted
693 according to the list @var{perm}. If, for example, @var{perm} is the list
694 @code{[3,2,1]}, the actual first row of the @var{LU} factorization is the
695 third row of the matrix @var{LU}. Second, the lower triangular factor of
696 m is the lower triangular part of @var{LU} with the diagonal entries
697 replaced by all ones. Third, the upper triangular factor of @var{M} is the
698 upper triangular part of @var{LU}.
700 (2) When the field is either @code{floatfield} or @code{complexfield}, the
701 numbers @var{lower-cnd} and @var{upper-cnd} are lower and upper bounds for
702 the infinity norm condition number of @var{M}. For all fields, the
703 condition number might not be estimated; for such fields, @code{lu_factor}
704 returns a two item list. Both the lower and upper bounds can differ from
705 their true values by arbitrarily large factors. (See also @mref{mat_cond}.)
707 The argument @var{M} must be a square matrix.
709 The optional argument @var{fld} must be a symbol that determines a ring or
710 field. The pre-defined fields and rings are:
712 (a) @code{generalring} -- the ring of Maxima expressions,
714 (b) @code{floatfield} -- the field of floating point numbers of the
715 type double,
717 (c) @code{complexfield} -- the field of complex floating point numbers of
718 the type double,
720 (d) @code{crering} -- the ring of Maxima CRE expressions,
722 (e) @code{rationalfield} -- the field of rational numbers,
724 (f) @code{runningerror} -- track the all floating point rounding errors,
726 (g) @code{noncommutingring} -- the ring of Maxima expressions where
727 multiplication is the non-commutative dot
728 operator.
730 When the field is @code{floatfield}, @code{complexfield}, or
731 @code{runningerror}, the algorithm uses partial pivoting; for all
732 other fields, rows are switched only when needed to avoid a zero
733 pivot.
735 Floating point addition arithmetic isn't associative, so the meaning
736 of 'field' differs from the mathematical definition.
738 A member of the field @code{runningerror} is a two member Maxima list
739 of the form @code{[x,n]},where @var{x} is a floating point number and
740 @code{n} is an integer. The relative difference between the 'true'
741 value of @code{x} and @code{x} is approximately bounded by the machine
742 epsilon times @code{n}. The running error bound drops some terms that
743 of the order the square of the machine epsilon.
745 There is no user-interface for defining a new field. A user that is
746 familiar with Common Lisp should be able to define a new field. To do
747 this, a user must define functions for the arithmetic operations and
748 functions for converting from the field representation to Maxima and
749 back. Additionally, for ordered fields (where partial pivoting will be
750 used), a user must define functions for the magnitude and for
751 comparing field members. After that all that remains is to define a
752 Common Lisp structure @code{mring}. The file @code{mring} has many
753 examples.
755 To compute the factorization, the first task is to convert each matrix
756 entry to a member of the indicated field. When conversion isn't
757 possible, the factorization halts with an error message. Members of
758 the field needn't be Maxima expressions. Members of the
759 @code{complexfield}, for example, are Common Lisp complex numbers. Thus
760 after computing the factorization, the matrix entries must be
761 converted to Maxima expressions.
763 See also @mref{get_lu_factors}.
765 Examples:
767 @c ===beg===
768 @c w[i,j] := random (1.0) + %i * random (1.0);
769 @c showtime : true$
770 @c M : genmatrix (w, 100, 100)$
771 @c lu_factor (M, complexfield)$
772 @c lu_factor (M, generalring)$
773 @c showtime : false$
774 @c M : matrix ([1 - z, 3], [3, 8 - z]);
775 @c lu_factor (M, generalring);
776 @c get_lu_factors (%);
777 @c %[1] . %[2] . %[3];
778 @c ===end===
779 @example
780 (%i1) w[i,j] := random (1.0) + %i * random (1.0);
781 (%o1) w := random(1.) + %i random(1.)
782 i, j
783 (%i2) showtime : true$
784 Evaluation took 0.00 seconds (0.00 elapsed)
785 (%i3) M : genmatrix (w, 100, 100)$
786 Evaluation took 7.40 seconds (8.23 elapsed)
787 (%i4) lu_factor (M, complexfield)$
788 Evaluation took 28.71 seconds (35.00 elapsed)
789 (%i5) lu_factor (M, generalring)$
790 Evaluation took 109.24 seconds (152.10 elapsed)
791 (%i6) showtime : false$
793 (%i7) M : matrix ([1 - z, 3], [3, 8 - z]);
794 [ 1 - z 3 ]
795 (%o7) [ ]
796 [ 3 8 - z ]
797 (%i8) lu_factor (M, generalring);
798 [ 1 - z 3 ]
799 [ ]
800 (%o8) [[ 3 9 ], [1, 2], generalring]
801 [ ----- - z - ----- + 8 ]
802 [ 1 - z 1 - z ]
803 (%i9) get_lu_factors (%);
804 [ 1 0 ] [ 1 - z 3 ]
805 [ 1 0 ] [ ] [ ]
806 (%o9) [[ ], [ 3 ], [ 9 ]]
807 [ 0 1 ] [ ----- 1 ] [ 0 - z - ----- + 8 ]
808 [ 1 - z ] [ 1 - z ]
809 (%i10) %[1] . %[2] . %[3];
810 [ 1 - z 3 ]
811 (%o10) [ ]
812 [ 3 8 - z ]
813 @end example
815 @opencatbox{Categories:}
816 @category{Matrix decompositions}
817 @category{Package linearalgebra}
818 @closecatbox
819 @end deffn
821 @c -----------------------------------------------------------------------------
822 @anchor{mat_cond}
823 @deffn {Function} mat_cond @
824 @fname{mat_cond} (@var{M}, 1) @
825 @fname{mat_cond} (@var{M}, inf)
827 Return the @var{p}-norm matrix condition number of the matrix
828 @var{m}. The allowed values for @var{p} are 1 and @var{inf}. This
829 function uses the LU factorization to invert the matrix @var{m}. Thus
830 the running time for @code{mat_cond} is proportional to the cube of
831 the matrix size; @code{lu_factor} determines lower and upper bounds
832 for the infinity norm condition number in time proportional to the
833 square of the matrix size.
835 @opencatbox{Categories:}
836 @category{Package linearalgebra}
837 @closecatbox
838 @end deffn
840 @c -----------------------------------------------------------------------------
841 @anchor{mat_norm}
842 @deffn {Function} mat_norm @
843 @fname{mat_norm} (@var{M}, 1) @
844 @fname{mat_norm} (@var{M}, inf) @
845 @fname{mat_norm} (@var{M}, frobenius)
847 Return the matrix @var{p}-norm of the matrix @var{M}. The allowed values for
848 @var{p} are 1, @code{inf}, and @code{frobenius} (the Frobenius matrix norm).
849 The matrix @var{M} should be an unblocked matrix.
851 @opencatbox{Categories:}
852 @category{Package linearalgebra}
853 @closecatbox
854 @end deffn
856 @c -----------------------------------------------------------------------------
857 @anchor{linearalgebra_matrixp}
858 @deffn {Function} matrixp @
859 @fname{matrixp} (@var{e}, @var{p}) @
860 @fname{matrixp} (@var{e})
862 Given an optional argument @var{p}, return @code{true} if @var{e} is
863 a matrix and @var{p} evaluates to @code{true} for every matrix element.
864 When @code{matrixp} is not given an optional argument, return @code{true}
865 if @code{e} is a matrix. In all other cases, return @code{false}.
867 See also @mref{blockmatrixp}
869 @opencatbox{Categories:}
870 @category{Package linearalgebra}
871 @category{Predicate functions}
872 @closecatbox
873 @end deffn
875 @c -----------------------------------------------------------------------------
876 @anchor{matrix_size}
877 @deffn {Function} matrix_size (@var{M})
879 Return a two member list that gives the number of rows and columns, respectively
880 of the matrix @var{M}.
882 @opencatbox{Categories:}
883 @category{Package linearalgebra}
884 @closecatbox
885 @end deffn
887 @c -----------------------------------------------------------------------------
888 @anchor{mat_fullunblocker}
889 @deffn {Function} mat_fullunblocker (@var{M})
891 If @var{M} is a block matrix, unblock the matrix to all levels. If @var{M} is
892 a matrix, return @var{M}; otherwise, signal an error.
894 @opencatbox{Categories:}
895 @category{Package linearalgebra}
896 @closecatbox
897 @end deffn
899 @c -----------------------------------------------------------------------------
900 @anchor{mat_trace}
901 @deffn {Function} mat_trace (@var{M})
903 Return the trace of the matrix @var{M}. If @var{M} isn't a matrix, return a
904 noun form. When @var{M} is a block matrix, @code{mat_trace(M)} returns
905 the same value as does @code{mat_trace(mat_unblocker(m))}.
907 @opencatbox{Categories:}
908 @category{Package linearalgebra}
909 @closecatbox
910 @end deffn
912 @c -----------------------------------------------------------------------------
913 @anchor{mat_unblocker}
914 @deffn {Function} mat_unblocker (@var{M})
916 If @var{M} is a block matrix, unblock @var{M} one level. If @var{M} is a
917 matrix, @code{mat_unblocker (M)} returns @var{M}; otherwise, signal an error.
919 Thus if each entry of @var{M} is matrix, @code{mat_unblocker (M)} returns an
920 unblocked matrix, but if each entry of @var{M} is a block matrix,
921 @code{mat_unblocker (M)} returns a block matrix with one less level of blocking.
923 If you use block matrices, most likely you'll want to set
924 @code{matrix_element_mult} to @code{"."} and @code{matrix_element_transpose} to
925 @code{'transpose}. See also @mref{mat_fullunblocker}.
927 Example:
929 @c ===beg===
930 @c A : matrix ([1, 2], [3, 4]);
931 @c B : matrix ([7, 8], [9, 10]);
932 @c matrix ([A, B]);
933 @c mat_unblocker (%);
934 @c ===end===
935 @example
936 (%i1) A : matrix ([1, 2], [3, 4]);
937 [ 1 2 ]
938 (%o1) [ ]
939 [ 3 4 ]
940 (%i2) B : matrix ([7, 8], [9, 10]);
941 [ 7 8 ]
942 (%o2) [ ]
943 [ 9 10 ]
944 (%i3) matrix ([A, B]);
945 @group
946 [ [ 1 2 ] [ 7 8 ] ]
947 (%o3) [ [ ] [ ] ]
948 [ [ 3 4 ] [ 9 10 ] ]
949 @end group
950 (%i4) mat_unblocker (%);
951 [ 1 2 7 8 ]
952 (%o4) [ ]
953 [ 3 4 9 10 ]
954 @end example
956 @opencatbox{Categories:}
957 @category{Package linearalgebra}
958 @closecatbox
959 @end deffn
961 @c -----------------------------------------------------------------------------
962 @anchor{nullspace}
963 @deffn {Function} nullspace (@var{M})
965 If @var{M} is a matrix, return @code{span (v_1, ..., v_n)}, where the set
966 @code{@{v_1, ..., v_n@}} is a basis for the nullspace of @var{M}. The span of
967 the empty set is @code{@{0@}}. Thus, when the nullspace has only one member,
968 return @code{span ()}.
970 @opencatbox{Categories:}
971 @category{Package linearalgebra}
972 @closecatbox
973 @end deffn
975 @c -----------------------------------------------------------------------------
976 @anchor{nullity}
977 @deffn {Function} nullity (@var{M})
979 If @var{M} is a matrix, return the dimension of the nullspace of @var{M}.
981 @opencatbox{Categories:}
982 @category{Package linearalgebra}
983 @closecatbox
984 @end deffn
986 @c -----------------------------------------------------------------------------
987 @anchor{orthogonal_complement}
988 @deffn {Function} orthogonal_complement (@var{v_1}, @dots{}, @var{v_n})
990 Return @code{span (u_1, ..., u_m)}, where the set @code{@{u_1, ..., u_m@}} is a
991 basis for the orthogonal complement of the set @code{(v_1, ..., v_n)}.
993 Each vector @var{v_1} through @var{v_n} must be a column vector.
995 @opencatbox{Categories:}
996 @category{Package linearalgebra}
997 @closecatbox
998 @end deffn
1000 @c -----------------------------------------------------------------------------
1001 @anchor{polytocompanion}
1002 @deffn {Function} polytocompanion (@var{p}, @var{x})
1004 If @var{p} is a polynomial in @var{x}, return the companion matrix of @var{p}.
1005 For a monic polynomial @var{p} of degree @var{n}, we have
1006 @code{@var{p} = (-1)^@var{n} charpoly (polytocompanion (@var{p}, @var{x}))}.
1008 When @var{p} isn't a polynomial in @var{x}, signal an error.
1010 @opencatbox{Categories:}
1011 @category{Package linearalgebra}
1012 @closecatbox
1013 @end deffn
1015 @c -----------------------------------------------------------------------------
1016 @anchor{ptringularize}
1017 @deffn {Function} ptriangularize (@var{M}, @var{v})
1019 If @var{M} is a matrix with each entry a polynomial in @var{v}, return
1020 a matrix @var{M2} such that
1022 (1) @var{M2} is upper triangular,
1024 (2) @code{@var{M2} = @var{E_n} ... @var{E_1} @var{M}},
1025 where @var{E_1} through @var{E_n} are elementary matrices
1026 whose entries are polynomials in @var{v},
1028 (3) @code{|det (@var{M})| = |det (@var{M2})|},
1030 Note: This function doesn't check that every entry is a polynomial in @var{v}.
1032 @opencatbox{Categories:}
1033 @category{Package linearalgebra}
1034 @closecatbox
1035 @end deffn
1037 @c -----------------------------------------------------------------------------
1038 @anchor{rowp}
1039 @deffn {Function} rowop (@var{M}, @var{i}, @var{j}, @var{theta})
1041 If @var{M} is a matrix, return the matrix that results from doing the
1042 row operation @code{R_i <- R_i - theta * R_j}. If @var{M} doesn't have a row
1043 @var{i} or @var{j}, signal an error.
1045 @opencatbox{Categories:}
1046 @category{Package linearalgebra}
1047 @closecatbox
1049 @end deffn
1051 @anchor{linalg_rank}
1052 @deffn {Function} linalg_rank (@var{M})
1054 Return the rank of the matrix @var{M}. This function is equivalent to
1055 function @mref{rank}, but it uses a different algorithm: it finds the
1056 @mref{columnspace} of the matrix and counts its elements, since the rank
1057 of a matrix is the dimension of its column space.
1059 @c ===beg===
1060 @c linalg_rank(matrix([1,2],[2,4]));
1061 @c linalg_rank(matrix([1,b],[c,d]));
1062 @c ===end===
1063 @example
1064 @group
1065 (%i1) linalg_rank(matrix([1,2],[2,4]));
1066 (%o1) 1
1067 @end group
1068 @group
1069 (%i2) linalg_rank(matrix([1,b],[c,d]));
1070 (%o2) 2
1071 @end group
1072 @end example
1074 @opencatbox{Categories:}
1075 @category{Package linearalgebra}
1076 @closecatbox
1077 @end deffn
1079 @c -----------------------------------------------------------------------------
1080 @anchor{rowswap}
1081 @deffn {Function} rowswap (@var{M}, @var{i}, @var{j})
1083 If @var{M} is a matrix, swap rows @var{i} and @var{j}. If @var{M} doesn't
1084 have a row @var{i} or @var{j}, signal an error.
1086 @opencatbox{Categories:}
1087 @category{Package linearalgebra}
1088 @closecatbox
1089 @end deffn
1091 @c -----------------------------------------------------------------------------
1092 @anchor{toeplitz}
1093 @deffn {Function} toeplitz @
1094 @fname{toeplitz} (@var{col}) @
1095 @fname{toeplitz} (@var{col}, @var{row})
1097 Return a Toeplitz matrix @var{T}. The first first column of @var{T} is
1098 @var{col}; except for the first entry, the first row of @var{T} is @var{row}.
1099 The default for @var{row} is complex conjugate of @var{col}. Example:
1101 @c ===beg===
1102 @c toeplitz([1,2,3],[x,y,z]);
1103 @c toeplitz([1,1+%i]);
1104 @c ===end===
1105 @example
1106 (%i1) toeplitz([1,2,3],[x,y,z]);
1107 @group
1108 [ 1 y z ]
1109 [ ]
1110 (%o1) [ 2 1 y ]
1111 [ ]
1112 [ 3 2 1 ]
1113 @end group
1114 (%i2) toeplitz([1,1+%i]);
1116 [ 1 1 - %I ]
1117 (%o2) [ ]
1118 [ %I + 1 1 ]
1119 @end example
1121 @opencatbox{Categories:}
1122 @category{Package linearalgebra}
1123 @closecatbox
1124 @end deffn
1126 @c -----------------------------------------------------------------------------
1127 @anchor{vandermonde_matrix}
1128 @deffn {Function} vandermonde_matrix ([@var{x_1}, ..., @var{x_n}])
1130 Return a @var{n} by @var{n} matrix whose @var{i}-th row is
1131 @code{[1, @var{x_i}, @var{x_i}^2, ... @var{x_i}^(@var{n}-1)]}.
1133 @opencatbox{Categories:}
1134 @category{Package linearalgebra}
1135 @closecatbox
1136 @end deffn
1138 @c -----------------------------------------------------------------------------
1139 @anchor{zerofor}
1140 @deffn {Function} zerofor @
1141 @fname{zerofor} (@var{M}) @
1142 @fname{zerofor} (@var{M}, @var{fld})
1144 Return a zero matrix that has the same shape as the matrix
1145 @var{M}. Every entry of the zero matrix is the
1146 additive identity of the field @var{fld}; the default for
1147 @var{fld} is @var{generalring}.
1149 The first argument @var{M} should be a square matrix or a
1150 non-matrix. When @var{M} is a matrix, each entry of @var{M} can be a
1151 square matrix -- thus @var{M} can be a blocked Maxima matrix. The
1152 matrix can be blocked to any (finite) depth.
1154 See also @mref{identfor}
1156 @opencatbox{Categories:}
1157 @category{Package linearalgebra}
1158 @closecatbox
1159 @end deffn
1161 @c -----------------------------------------------------------------------------
1162 @anchor{zeromatrixp}
1163 @deffn {Function} zeromatrixp (@var{M})
1165 If @var{M} is not a block matrix, return @code{true} if
1166 @code{is (equal (@var{e}, 0))} is true for each element @var{e} of the matrix
1167 @var{M}. If @var{M} is a block matrix, return @code{true} if @code{zeromatrixp}
1168 evaluates to @code{true} for each element of @var{e}.
1170 @opencatbox{Categories:}
1171 @category{Package linearalgebra}
1172 @category{Predicate functions}
1173 @closecatbox
1174 @end deffn