comparison doc/interpreter/linalg.txi @ 3372:f16c2ce14886

[project @ 1999-11-23 19:07:09 by jwe]
author jwe
date Tue, 23 Nov 1999 19:07:18 +0000
parents bfe1573bd2ae
children 36405da8e173
comparison
equal deleted inserted replaced
3371:86873384cd10 3372:f16c2ce14886
17 @end menu 17 @end menu
18 18
19 @node Basic Matrix Functions, Matrix Factorizations, Linear Algebra, Linear Algebra 19 @node Basic Matrix Functions, Matrix Factorizations, Linear Algebra, Linear Algebra
20 @section Basic Matrix Functions 20 @section Basic Matrix Functions
21 21
22 @deftypefn {Loadable Function} {@var{aa} =} balance (@var{a}, @var{opt}) 22 DOCSTRING(balance)
23 @deftypefnx {Loadable Function} {[@var{dd}, @var{aa}] =} balance (@var{a}, @var{opt})
24 @deftypefnx {Loadable Function} {[@var{cc}, @var{dd}, @var{aa}, @var{bb]} =} balance (@var{a}, @var{b}, @var{opt})
25 23
26 @code{[dd, aa] = balance (a)} returns @code{aa = dd \ a * dd}. 24 DOCSTRING(cond)
27 @code{aa} is a matrix whose row and column norms are roughly equal in
28 magnitude, and @code{dd} = @code{p * d}, where @code{p} is a permutation
29 matrix and @code{d} is a diagonal matrix of powers of two. This allows
30 the equilibration to be computed without roundoff. Results of
31 eigenvalue calculation are typically improved by balancing first.
32 25
33 @code{[cc, dd, aa, bb] = balance (a, b)} returns @code{aa = cc*a*dd} and 26 DOCSTRING(det)
34 @code{bb = cc*b*dd)}, where @code{aa} and @code{bb} have non-zero
35 elements of approximately the same magnitude and @code{cc} and @code{dd}
36 are permuted diagonal matrices as in @code{dd} for the algebraic
37 eigenvalue problem.
38 27
39 The eigenvalue balancing option @code{opt} is selected as follows: 28 DOCSTRING(eig)
40 29
41 @table @asis 30 DOCSTRING(givens)
42 @item @code{"N"}, @code{"n"}
43 No balancing; arguments copied, transformation(s) set to identity.
44 31
45 @item @code{"P"}, @code{"p"} 32 DOCSTRING(inv)
46 Permute argument(s) to isolate eigenvalues where possible.
47 33
48 @item @code{"S"}, @code{"s"} 34 DOCSTRING(norm)
49 Scale to improve accuracy of computed eigenvalues.
50 35
51 @item @code{"B"}, @code{"b"} 36 DOCSTRING(null)
52 Permute and scale, in that order. Rows/columns of a (and b)
53 that are isolated by permutation are not scaled. This is the default
54 behavior.
55 @end table
56 37
57 Algebraic eigenvalue balancing uses standard @sc{Lapack} routines. 38 DOCSTRING(orth)
58 39
59 Generalized eigenvalue problem balancing uses Ward's algorithm 40 DOCSTRING(pinv)
60 (SIAM Journal on Scientific and Statistical Computing, 1981).
61 @end deftypefn
62 41
63 @deftypefn {} {} cond (@var{a}) 42 DOCSTRING(rank)
64 Compute the (two-norm) condition number of a matrix. @code{cond (a)} is
65 defined as @code{norm (a) * norm (inv (a))}, and is computed via a
66 singular value decomposition.
67 @end deftypefn
68 43
69 @deftypefn {Loadable Function} {} det (@var{a}) 44 DOCSTRING(trace)
70 Compute the determinant of @var{a} using @sc{Linpack}.
71 @end deftypefn
72
73 @deftypefn {Loadable Function} {@var{lambda} =} eig (@var{a})
74 @deftypefnx {Loadable Function} {[@var{v}, @var{lambda}] =} eig (@var{a})
75 The eigenvalues (and eigenvectors) of a matrix are computed in a several
76 step process which begins with a Hessenberg decomposition, followed by a
77 Schur decomposition, from which the eigenvalues are apparent. The
78 eigenvectors, when desired, are computed by further manipulations of the
79 Schur decomposition.
80 @end deftypefn
81
82 @deftypefn {Loadable Function} {@var{G} =} givens (@var{x}, @var{y})
83 @deftypefnx {Loadable Function} {[@var{c}, @var{s}] =} givens (@var{x}, @var{y})
84 @iftex
85 @tex
86 Return a $2\times 2$ orthogonal matrix
87 $$
88 G = \left[\matrix{c & s\cr -s'& c\cr}\right]
89 $$
90 such that
91 $$
92 G \left[\matrix{x\cr y}\right] = \left[\matrix{\ast\cr 0}\right]
93 $$
94 with $x$ and $y$ scalars.
95 @end tex
96 @end iftex
97 @ifinfo
98 Return a 2 by 2 orthogonal matrix
99 @code{@var{G} = [@var{c} @var{s}; -@var{s}' @var{c}]} such that
100 @code{@var{G} [@var{x}; @var{y}] = [*; 0]} with @var{x} and @var{y} scalars.
101 @end ifinfo
102
103 For example,
104
105 @example
106 @group
107 givens (1, 1)
108 @result{} 0.70711 0.70711
109 -0.70711 0.70711
110 @end group
111 @end example
112 @end deftypefn
113
114 @deftypefn {Loadable Function} {} inv (@var{a})
115 @deftypefnx {Loadable Function} {} inverse (@var{a})
116 Compute the inverse of the square matrix @var{a}.
117 @end deftypefn
118
119 @deftypefn {Function File} {} norm (@var{a}, @var{p})
120 Compute the p-norm of the matrix @var{a}. If the second argument is
121 missing, @code{p = 2} is assumed.
122
123 If @var{a} is a matrix:
124
125 @table @asis
126 @item @var{p} = @code{1}
127 1-norm, the largest column sum of @var{a}.
128
129 @item @var{p} = @code{2}
130 Largest singular value of @var{a}.
131
132 @item @var{p} = @code{Inf}
133 @cindex infinity norm
134 Infinity norm, the largest row sum of @var{a}.
135
136 @item @var{p} = @code{"fro"}
137 @cindex Frobenius norm
138 Frobenius norm of @var{a}, @code{sqrt (sum (diag (@var{a}' * @var{a})))}.
139 @end table
140
141 If @var{a} is a vector or a scalar:
142
143 @table @asis
144 @item @var{p} = @code{Inf}
145 @code{max (abs (@var{a}))}.
146
147 @item @var{p} = @code{-Inf}
148 @code{min (abs (@var{a}))}.
149
150 @item other
151 p-norm of @var{a}, @code{(sum (abs (@var{a}) .^ @var{p})) ^ (1/@var{p})}.
152 @end table
153 @end deftypefn
154
155 @deftypefn {Function File} {} null (@var{a}, @var{tol})
156 Return an orthonormal basis of the null space of @var{a}.
157
158 The dimension of the null space is taken as the number of singular
159 values of @var{a} not greater than @var{tol}. If the argument @var{tol}
160 is missing, it is computed as
161
162 @example
163 max (size (@var{a})) * max (svd (@var{a})) * eps
164 @end example
165 @end deftypefn
166
167 @deftypefn {Function File} {} orth (@var{a}, @var{tol})
168 Return an orthonormal basis of the range space of @var{a}.
169
170 The dimension of the range space is taken as the number of singular
171 values of @var{a} greater than @var{tol}. If the argument @var{tol} is
172 missing, it is computed as
173
174 @example
175 max (size (@var{a})) * max (svd (@var{a})) * eps
176 @end example
177 @end deftypefn
178
179 @deftypefn {Function File} {} pinv (@var{x}, @var{tol})
180 Return the pseudoinverse of @var{x}. Singular values less than
181 @var{tol} are ignored.
182
183 If the second argument is omitted, it is assumed that
184
185 @example
186 tol = max (size (@var{x})) * sigma_max (@var{x}) * eps,
187 @end example
188
189 @noindent
190 where @code{sigma_max (@var{x})} is the maximal singular value of @var{x}.
191 @end deftypefn
192
193 @deftypefn {Function File} {} rank (@var{a}, @var{tol})
194 Compute the rank of @var{a}, using the singular value decomposition.
195 The rank is taken to be the number of singular values of @var{a} that
196 are greater than the specified tolerance @var{tol}. If the second
197 argument is omitted, it is taken to be
198
199 @example
200 tol = max (size (@var{a})) * sigma (1) * eps;
201 @end example
202
203 @noindent
204 where @code{eps} is machine precision and @code{sigma} is the largest
205 singular value of @var{a}.
206 @end deftypefn
207
208 @deftypefn {Function File} {} trace (@var{a})
209 Compute the trace of @var{a}, @code{sum (diag (@var{a}))}.
210 @end deftypefn
211 45
212 @node Matrix Factorizations, Functions of a Matrix, Basic Matrix Functions, Linear Algebra 46 @node Matrix Factorizations, Functions of a Matrix, Basic Matrix Functions, Linear Algebra
213 @section Matrix Factorizations 47 @section Matrix Factorizations
214 48
215 @deftypefn {Loadable Function} {} chol (@var{a}) 49 DOCSTRING(chol)
216 @cindex Cholesky factorization
217 Compute the Cholesky factor, @var{r}, of the symmetric positive definite
218 matrix @var{a}, where
219 @iftex
220 @tex
221 $ R^T R = A $.
222 @end tex
223 @end iftex
224 @ifinfo
225 50
226 @example 51 DOCSTRING(hess)
227 r' * r = a.
228 @end example
229 @end ifinfo
230 @end deftypefn
231 52
232 @deftypefn {Loadable Function} {@var{h} =} hess (@var{a}) 53 DOCSTRING(lu)
233 @deftypefnx {Loadable Function} {[@var{p}, @var{h}] =} hess (@var{a})
234 @cindex Hessenberg decomposition
235 Compute the Hessenberg decomposition of the matrix @var{a}.
236 54
237 The Hessenberg decomposition is usually used as the first step in an 55 DOCSTRING(qr)
238 eigenvalue computation, but has other applications as well (see Golub,
239 Nash, and Van Loan, IEEE Transactions on Automatic Control, 1979. The
240 Hessenberg decomposition is
241 @iftex
242 @tex
243 $$
244 A = PHP^T
245 $$
246 where $P$ is a square unitary matrix ($P^HP = I$), and $H$
247 is upper Hessenberg ($H_{i,j} = 0, \forall i \ge j+1$).
248 @end tex
249 @end iftex
250 @ifinfo
251 @code{p * h * p' = a} where @code{p} is a square unitary matrix
252 (@code{p' * p = I}, using complex-conjugate transposition) and @code{h}
253 is upper Hessenberg (@code{i >= j+1 => h (i, j) = 0}).
254 @end ifinfo
255 @end deftypefn
256 56
257 @deftypefn {Loadable Function} {[@var{l}, @var{u}, @var{p}] =} lu (@var{a}) 57 DOCSTRING(qz)
258 @cindex LU decomposition
259 Compute the LU decomposition of @var{a}, using subroutines from
260 @sc{Lapack}. The result is returned in a permuted form, according to
261 the optional return value @var{p}. For example, given the matrix
262 @code{a = [1, 2; 3, 4]},
263 58
264 @example 59 DOCSTRING(qzhess)
265 [l, u, p] = lu (a)
266 @end example
267 60
268 @noindent 61 DOCSTRING(schur)
269 returns
270 62
271 @example 63 DOCSTRING(svd)
272 l =
273
274 1.00000 0.00000
275 0.33333 1.00000
276
277 u =
278
279 3.00000 4.00000
280 0.00000 0.66667
281
282 p =
283
284 0 1
285 1 0
286 @end example
287 @end deftypefn
288
289 @deftypefn {Loadable Function} {[@var{q}, @var{r}, @var{p}] =} qr (@var{a})
290 @cindex QR factorization
291 Compute the QR factorization of @var{a}, using standard @sc{Lapack}
292 subroutines. For example, given the matrix @code{a = [1, 2; 3, 4]},
293
294 @example
295 [q, r] = qr (a)
296 @end example
297
298 @noindent
299 returns
300
301 @example
302 q =
303
304 -0.31623 -0.94868
305 -0.94868 0.31623
306
307 r =
308
309 -3.16228 -4.42719
310 0.00000 -0.63246
311 @end example
312
313 The @code{qr} factorization has applications in the solution of least
314 squares problems
315 @iftex
316 @tex
317 $$
318 \min_x \left\Vert A x - b \right\Vert_2
319 $$
320 @end tex
321 @end iftex
322 @ifinfo
323
324 @example
325 @code{min norm(A x - b)}
326 @end example
327
328 @end ifinfo
329 for overdetermined systems of equations (i.e.,
330 @iftex
331 @tex
332 $A$
333 @end tex
334 @end iftex
335 @ifinfo
336 @code{a}
337 @end ifinfo
338 is a tall, thin matrix). The QR factorization is
339 @iftex
340 @tex
341 $QR = A$ where $Q$ is an orthogonal matrix and $R$ is upper triangular.
342 @end tex
343 @end iftex
344 @ifinfo
345 @code{q * r = a} where @code{q} is an orthogonal matrix and @code{r} is
346 upper triangular.
347 @end ifinfo
348
349 The permuted QR factorization @code{[@var{q}, @var{r}, @var{p}] =
350 qr (@var{a})} forms the QR factorization such that the diagonal
351 entries of @code{r} are decreasing in magnitude order. For example,
352 given the matrix @code{a = [1, 2; 3, 4]},
353
354 @example
355 [q, r, pi] = qr(a)
356 @end example
357
358 @noindent
359 returns
360
361 @example
362 q =
363
364 -0.44721 -0.89443
365 -0.89443 0.44721
366
367 r =
368
369 -4.47214 -3.13050
370 0.00000 0.44721
371
372 p =
373
374 0 1
375 1 0
376 @end example
377
378 The permuted @code{qr} factorization @code{[q, r, p] = qr (a)}
379 factorization allows the construction of an orthogonal basis of
380 @code{span (a)}.
381 @end deftypefn
382
383
384 @deftypefn {Function File} {@var{lambda} =} qz (@var{a}, @var{b})
385 Generalized eigenvalue problem @math{A x = s B x},
386 @var{QZ} decomposition. Three ways to call:
387 @enumerate
388 @item @code{lambda = qz(A,B)}
389
390 Computes the generalized eigenvalues @var{lambda} of @math{(A - sB)}.
391
392 @item @code{[AA, BB, Q, Z @{, V, W, lambda@}] = qz (A, B)}
393
394 Computes qz decomposition, generalized eigenvectors, and
395 generalized eigenvalues of @math{(A - sB)}
396 @example
397 @group
398 A V = B V diag(lambda)
399 W' A = diag(lambda) W' B
400 AA = Q'*A*Z, BB = Q'*B*Z with Q, Z orthogonal (unitary)= I
401 @end group
402 @end example
403
404 @item @code{[AA,BB,Z@{,lambda@}] = qz(A,B,opt)}
405
406 As in form [2], but allows ordering of generalized eigenpairs
407 for (e.g.) solution of discrete time algebraic Riccati equations.
408 Form 3 is not available for complex matrices and does not compute
409 the generalized eigenvectors V, W, nor the orthogonal matrix Q.
410 @table @var
411 @item opt
412 for ordering eigenvalues of the GEP pencil. The leading block
413 of the revised pencil contains all eigenvalues that satisfy:
414 @table @code
415 @item "N"
416 = unordered (default)
417
418 @item "S"
419 = small: leading block has all |lambda| <=1
420
421 @item "B"
422 = big: leading block has all |lambda >= 1
423
424 @item "-"
425 = negative real part: leading block has all eigenvalues
426 in the open left half-plant
427
428 @item "+"
429 = nonnegative real part: leading block has all eigenvalues
430 in the closed right half-plane
431 @end table
432 @end table
433 @end enumerate
434
435 Note: qz performs permutation balancing, but not scaling (see balance).
436 Order of output arguments was selected for compatibility with MATLAB
437
438 See also: balance, dare, eig, schur
439 @end deftypefn
440
441 @deftypefn {Function File} {[@var{aa}, @var{bb}, @var{q}, @var{z}] =} qzhess (@var{a}, @var{b})
442 Compute the Hessenberg-triangular decomposition of the matrix pencil
443 @code{(@var{a}, @var{b})}, returning
444 @code{@var{aa} = @var{q} * @var{a} * @var{z}},
445 @code{@var{bb} = @var{q} * @var{b} * @var{z}}, with @var{q} and @var{z}
446 orthogonal. For example,
447
448 @example
449 @group
450 [aa, bb, q, z] = qzhess ([1, 2; 3, 4], [5, 6; 7, 8])
451 @result{} aa = [ -3.02244, -4.41741; 0.92998, 0.69749 ]
452 @result{} bb = [ -8.60233, -9.99730; 0.00000, -0.23250 ]
453 @result{} q = [ -0.58124, -0.81373; -0.81373, 0.58124 ]
454 @result{} z = [ 1, 0; 0, 1 ]
455 @end group
456 @end example
457
458 The Hessenberg-triangular decomposition is the first step in
459 Moler and Stewart's QZ decomposition algorithm.
460
461 Algorithm taken from Golub and Van Loan, @cite{Matrix Computations, 2nd
462 edition}.
463 @end deftypefn
464
465 @deftypefn {Loadable Function} {} qzval (@var{a}, @var{b})
466 Compute generalized eigenvalues of the matrix pencil
467 @iftex
468 @tex
469 $a - \lambda b$.
470 @end tex
471 @end iftex
472 @ifinfo
473 @code{@var{a} - lambda @var{b}}.
474 @end ifinfo
475
476 The arguments @var{a} and @var{b} must be real matrices.
477 @end deftypefn
478
479 @deftypefn {Loadable Function} {@var{s} =} schur (@var{a})
480 @deftypefnx {Loadable Function} {[@var{u}, @var{s}] =} schur (@var{a}, @var{opt})
481 @cindex Schur decomposition
482 The Schur decomposition is used to compute eigenvalues of a
483 square matrix, and has applications in the solution of algebraic
484 Riccati equations in control (see @code{are} and @code{dare}).
485 @code{schur} always returns
486 @iftex
487 @tex
488 $S = U^T A U$
489 @end tex
490 @end iftex
491 @ifinfo
492 @code{s = u' * a * u}
493 @end ifinfo
494 where
495 @iftex
496 @tex
497 $U$
498 @end tex
499 @end iftex
500 @ifinfo
501 @code{u}
502 @end ifinfo
503 is a unitary matrix
504 @iftex
505 @tex
506 ($U^T U$ is identity)
507 @end tex
508 @end iftex
509 @ifinfo
510 (@code{u'* u} is identity)
511 @end ifinfo
512 and
513 @iftex
514 @tex
515 $S$
516 @end tex
517 @end iftex
518 @ifinfo
519 @code{s}
520 @end ifinfo
521 is upper triangular. The eigenvalues of
522 @iftex
523 @tex
524 $A$ (and $S$)
525 @end tex
526 @end iftex
527 @ifinfo
528 @code{a} (and @code{s})
529 @end ifinfo
530 are the diagonal elements of
531 @iftex
532 @tex
533 $S$
534 @end tex
535 @end iftex
536 @ifinfo
537 @code{s}
538 @end ifinfo
539 If the matrix
540 @iftex
541 @tex
542 $A$
543 @end tex
544 @end iftex
545 @ifinfo
546 @code{a}
547 @end ifinfo
548 is real, then the real Schur decomposition is computed, in which the
549 matrix
550 @iftex
551 @tex
552 $U$
553 @end tex
554 @end iftex
555 @ifinfo
556 @code{u}
557 @end ifinfo
558 is orthogonal and
559 @iftex
560 @tex
561 $S$
562 @end tex
563 @end iftex
564 @ifinfo
565 @code{s}
566 @end ifinfo
567 is block upper triangular
568 with blocks of size at most
569 @iftex
570 @tex
571 $2\times 2$
572 @end tex
573 @end iftex
574 @ifinfo
575 @code{2 x 2}
576 @end ifinfo
577 blocks along the diagonal. The diagonal elements of
578 @iftex
579 @tex
580 $S$
581 @end tex
582 @end iftex
583 @ifinfo
584 @code{s}
585 @end ifinfo
586 (or the eigenvalues of the
587 @iftex
588 @tex
589 $2\times 2$
590 @end tex
591 @end iftex
592 @ifinfo
593 @code{2 x 2}
594 @end ifinfo
595 blocks, when
596 appropriate) are the eigenvalues of
597 @iftex
598 @tex
599 $A$
600 @end tex
601 @end iftex
602 @ifinfo
603 @code{a}
604 @end ifinfo
605 and
606 @iftex
607 @tex
608 $S$.
609 @end tex
610 @end iftex
611 @ifinfo
612 @code{s}.
613 @end ifinfo
614
615 The eigenvalues are optionally ordered along the diagonal according to
616 the value of @code{opt}. @code{opt = "a"} indicates that all
617 eigenvalues with negative real parts should be moved to the leading
618 block of
619 @iftex
620 @tex
621 $S$
622 @end tex
623 @end iftex
624 @ifinfo
625 @code{s}
626 @end ifinfo
627 (used in @code{are}), @code{opt = "d"} indicates that all eigenvalues
628 with magnitude less than one should be moved to the leading block of
629 @iftex
630 @tex
631 $S$
632 @end tex
633 @end iftex
634 @ifinfo
635 @code{s}
636 @end ifinfo
637 (used in @code{dare}), and @code{opt = "u"}, the default, indicates that
638 no ordering of eigenvalues should occur. The leading
639 @iftex
640 @tex
641 $k$
642 @end tex
643 @end iftex
644 @ifinfo
645 @code{k}
646 @end ifinfo
647 columns of
648 @iftex
649 @tex
650 $U$
651 @end tex
652 @end iftex
653 @ifinfo
654 @code{u}
655 @end ifinfo
656 always span the
657 @iftex
658 @tex
659 $A$-invariant
660 @end tex
661 @end iftex
662 @ifinfo
663 @code{a}-invariant
664 @end ifinfo
665 subspace corresponding to the
666 @iftex
667 @tex
668 $k$
669 @end tex
670 @end iftex
671 @ifinfo
672 @code{k}
673 @end ifinfo
674 leading eigenvalues of
675 @iftex
676 @tex
677 $S$.
678 @end tex
679 @end iftex
680 @ifinfo
681 @code{s}.
682 @end ifinfo
683 @end deftypefn
684
685 @deftypefn {Loadable Function} {@var{s} =} svd (@var{a})
686 @deftypefnx {Loadable Function} {[@var{u}, @var{s}, @var{v}] =} svd (@var{a})
687 @cindex singular value decomposition
688 Compute the singular value decomposition of @var{a}
689 @iftex
690 @tex
691 $$
692 A = U\Sigma V^H
693 $$
694 @end tex
695 @end iftex
696 @ifinfo
697
698 @example
699 a = u * sigma * v'
700 @end example
701 @end ifinfo
702
703 The function @code{svd} normally returns the vector of singular values.
704 If asked for three return values, it computes
705 @iftex
706 @tex
707 $U$, $S$, and $V$.
708 @end tex
709 @end iftex
710 @ifinfo
711 U, S, and V.
712 @end ifinfo
713 For example,
714
715 @example
716 svd (hilb (3))
717 @end example
718
719 @noindent
720 returns
721
722 @example
723 ans =
724
725 1.4083189
726 0.1223271
727 0.0026873
728 @end example
729
730 @noindent
731 and
732
733 @example
734 [u, s, v] = svd (hilb (3))
735 @end example
736
737 @noindent
738 returns
739
740 @example
741 u =
742
743 -0.82704 0.54745 0.12766
744 -0.45986 -0.52829 -0.71375
745 -0.32330 -0.64901 0.68867
746
747 s =
748
749 1.40832 0.00000 0.00000
750 0.00000 0.12233 0.00000
751 0.00000 0.00000 0.00269
752
753 v =
754
755 -0.82704 0.54745 0.12766
756 -0.45986 -0.52829 -0.71375
757 -0.32330 -0.64901 0.68867
758 @end example
759
760 If given a second argument, @code{svd} returns an economy-sized
761 decomposition, eliminating the unnecessary rows or columns of @var{u} or
762 @var{v}.
763 @end deftypefn
764 64
765 @node Functions of a Matrix, , Matrix Factorizations, Linear Algebra 65 @node Functions of a Matrix, , Matrix Factorizations, Linear Algebra
766 @section Functions of a Matrix 66 @section Functions of a Matrix
767 67
768 @deftypefn {Loadable Function} {} expm (@var{a}) 68 DOCSTRING(expm)
769 Return the exponential of a matrix, defined as the infinite Taylor
770 series
771 @iftex
772 @tex
773 $$
774 \exp (A) = I + A + {A^2 \over 2!} + {A^3 \over 3!} + \cdots
775 $$
776 @end tex
777 @end iftex
778 @ifinfo
779 69
780 @example 70 DOCSTRING(logm)
781 expm(a) = I + a + a^2/2! + a^3/3! + ...
782 @end example
783 71
784 @end ifinfo 72 DOCSTRING(sqrtm)
785 The Taylor series is @emph{not} the way to compute the matrix
786 exponential; see Moler and Van Loan, @cite{Nineteen Dubious Ways to
787 Compute the Exponential of a Matrix}, SIAM Review, 1978. This routine
788 uses Ward's diagonal
789 @iftex
790 @tex
791 Pad\'e
792 @end tex
793 @end iftex
794 @ifinfo
795 Pade'
796 @end ifinfo
797 approximation method with three step preconditioning (SIAM Journal on
798 Numerical Analysis, 1977). Diagonal
799 @iftex
800 @tex
801 Pad\'e
802 @end tex
803 @end iftex
804 @ifinfo
805 Pade'
806 @end ifinfo
807 approximations are rational polynomials of matrices
808 @iftex
809 @tex
810 $D_q(a)^{-1}N_q(a)$
811 @end tex
812 @end iftex
813 @ifinfo
814 73
815 @example 74 DOCSTRING(kron)
816 -1
817 D (a) N (a)
818 @end example
819 75
820 @end ifinfo 76 DOCSTRING(syl)
821 whose Taylor series matches the first
822 @iftex
823 @tex
824 $2 q + 1 $
825 @end tex
826 @end iftex
827 @ifinfo
828 @code{2q+1}
829 @end ifinfo
830 terms of the Taylor series above; direct evaluation of the Taylor series
831 (with the same preconditioning steps) may be desirable in lieu of the
832 @iftex
833 @tex
834 Pad\'e
835 @end tex
836 @end iftex
837 @ifinfo
838 Pade'
839 @end ifinfo
840 approximation when
841 @iftex
842 @tex
843 $D_q(a)$
844 @end tex
845 @end iftex
846 @ifinfo
847 @code{Dq(a)}
848 @end ifinfo
849 is ill-conditioned.
850 @end deftypefn
851
852 @deftypefn {Loadable Function} {} logm (@var{a})
853 Compute the matrix logarithm of the square matrix @var{a}. Note that
854 this is currently implemented in terms of an eigenvalue expansion and
855 needs to be improved to be more robust.
856 @end deftypefn
857
858 @deftypefn {Loadable Function} {} sqrtm (@var{a})
859 Compute the matrix square root of the square matrix @var{a}. Note that
860 this is currently implemented in terms of an eigenvalue expansion and
861 needs to be improved to be more robust.
862 @end deftypefn
863
864 @deftypefn {Function File} {} kron (@var{a}, @var{b})
865 Form the kronecker product of two matrices, defined block by block as
866
867 @example
868 x = [a(i, j) b]
869 @end example
870
871 For example,
872
873 @example
874 @group
875 kron (1:4, ones (3, 1))
876 @result{} 1 2 3 4
877 1 2 3 4
878 1 2 3 4
879 @end group
880 @end example
881 @end deftypefn
882
883 @deftypefn {Loadable Function} {@var{x} =} syl (@var{a}, @var{b}, @var{c})
884 Solve the Sylvester equation
885 @iftex
886 @tex
887 $$
888 A X + X B + C = 0
889 $$
890 @end tex
891 @end iftex
892 @ifinfo
893
894 @example
895 A X + X B + C = 0
896 @end example
897 @end ifinfo
898 using standard @sc{Lapack} subroutines. For example,
899
900 @example
901 @group
902 syl ([1, 2; 3, 4], [5, 6; 7, 8], [9, 10; 11, 12])
903 @result{} [ -0.50000, -0.66667; -0.66667, -0.50000 ]
904 @end group
905 @end example
906 @end deftypefn