functions in linalg.i -
_sv_get_dcmp
|
local sigma, u, vt; _sv_get_dcmp(a); -or- local sigma, u, vt; _sv_get_dcmp(a, full) Private routine used to extract Singular Value Decomposition of A into external symbols: SIGMA, U and VT. | |
SEE ALSO: | sv_dcmp, sv_solve_trunc, sv_solve_wiener |
cholesky
|
cholesky(a) -or- cholesky(a, 0/1) Given a symmetric positive definite matrix A, returns a lower triangular matrix C which is the Cholesky decomposition of A such that: A = C(+,)*C(+,); If optional second argument is true (non-nil and non-zero), the scratch values in the upper triangular part of C are left unchanged; otherwise (the default behavior), the upper triangular part of C is filled with zeros. |
diag
|
diag(a) Returns the diagonal of matrix A (if A is a 2D array) or a square diagonal matrix with diagonal elements equal to those of vector A (if A is a 1D array). | |
SEE ALSO: | trace, unit |
euclidean_norm
|
euclidean_norm(x) Returns the Euclidian norm of X: sqrt(sum(X*X)), taking care of overflows. |
gram_schmidt_orthonormalization
|
gram_schmidt_orthonormalization, b; -or- gram_schmidt_orthonormalization(b) Performs Gram-Schmidt orthonormalization of basis functions B. If B is an array of pointers, then the input basis vectors are *B(i) for i=1,..,numberof(B); otherwise, the input basis vectors are B(..,i) for i=1,..,dimsof(B)(0). When called as a subroutine, the operation is done "in-place". | |
SEE ALSO: | SVdec |
pm
|
pm, a; Prints matrix A. Keyword FILE can be set to the name/stream of output file; the default is to use standard output. If keyword NDIGITS is set, floating-point values get printed with that number of significant digits; alternatively, keyword FORMAT may be set with the format for each element of A -- for complex numbers, the real and imaginary parts use the same format. | |
SEE ALSO: | write |
sv_dcmp
|
sv_dcmp(a) -or- sv_dcmp(a, full) Computes the singular value decomposition of matrix A. The result is an array of pointers: [&SIGMA, &U, &VT] where SIGMA is the vector of singular values of A, and the columns of U and the rows of VT are the left and right singular vectors (respectively). Using matrix notation: A = U.diag(SIGMA).V' = U.diag(SIGMA).VT and using Yorick notation: A = (U*SIGMA(-,))(,+)*VT(+,) = U(,+)*(SIGMA*VT)(+,) FULL has the same meaning as keyword FULL in SVdec (to see). | |
SEE ALSO: |
sv_dcmp,
sv_solve_trunc,
sv_solve_wiener,
SVdec |
sv_intro
|
sv_intro - introduction to SVD Yorick package Notes about matrix multiplication in Yorick: A.B = A(,+)*B(+,) A'.B = A(+,)*B(+,) // transpose of A times B A.B' = A(,+)*B(,+) // A times transpose of B diag(s).A = diag(s)(,+)*A(+,) = s*A = A*s diag(s).A' = diag(s)(,+)*A(,+) = transpose(s(-,)*A) = transpose(A*s(-,)) A.diag(s) = A(,+)*diag(s)(+,) = A*s(-,) = s(-,)*A A'.diag(s) = A(+,)*diag(s)(+,) = A(+,)*diag(s)(,+) = transpose(A*s) = transpose(s*A) Singular Value Decomposition: A = U.diag(SIGMA).V' = U.diag(SIGMA).VT = (U*SIGMA(-,))(,+)*VT(+,) = U(,+)*(SIGMA*VT)(+,) where: SIGMA = SVdec(A, U, VT) Columns of U and V form an orthonormal basis (i.e. U and V are column-orthonormal): U'.U = V'.V = I in Yorick notation: U(+,)*U(+,) = V(+,)*V(+,) = VT(,+)*VT(,+) = unit(N) Note (to be verified): if U and/or V are square, they are also row-orthonormal: U'.U = U.U' = I (if U is square) V'.V = V.V' = I (if V is square) Generalized-inverse of A: AP = V.diag(1/SIGMA).U' = VT'.diag(1/SIGMA).U' = VT(+,)*((1.0/SIGMA)(-,)*U)(,+) = ((1.0/SIGMA)*VT)(+,)*U(,+) Solving a linear problem: A.x ~ b with SVD (taking care of index ordering for faster matrix multiplication): X = V.diag(W).U'.B = (W*VT)(+,)*U(,+))(,+)*B(+,..) = (W*VT)(+,)*(U(+,)*B(+,..))(+,..) // sum over 1st indices: faster where W is an approximation of 1/SIGMA. | |
SEE ALSO: |
sv_dcmp,
sv_solve_trunc,
sv_solve_wiener,
SVdec, SVsolve |
sv_solve_trunc
|
sv_solve_trunc(a, b, eps) Solve linear problem A.x = b by truncated singular value method. A is either a matrix or the singular value value decomposition as returned by sv_dcmp (to see). B is the right hand side vector (or array to solve for several right hand sides at a time). EPS (in the range [0,1]) is the minimum relative singular value to keep, i.e. all singular values less than EPS*max(SIGMA) get discarded (SIGMA is the vector of singular values of A). The result is: (W*VT)(+,)*(U(+,)*B(+,..))(+,..) where SIGMA, U and VT are the components of the singular value decomposition of A and W is: W(i) = 1/SIGMA(i) if SIGMA(i) > EPS*max(SIGMA); 0 otherwise. | |
SEE ALSO: | sv_dcmp, sv_intro, sv_solve_wiener |
sv_solve_wiener
|
sv_solve_wiener(a, b, eps) Solve linear problem A.x = b by Wiener filtering of the singular values. A is either a matrix or the singular value value decomposition as returned by sv_dcmp (to see). B is the right hand side vector (or array to solve for several right hand sides at a time). EPS (in the range [0,1]) is the relative singular value filter level. The result is: (W*VT)(+,)*(U(+,)*B(+,..))(+,..) where SIGMA, U and VT are the components of the singular value decomposition of A and W is: W = SIGMA/(SIGMA^2 + (EPS*max(SIGMA))^2) | |
SEE ALSO: | sv_dcmp, sv_intro, sv_solve_trunc |
trace
|
trace(a) Returns the trace of matrix A. | |
SEE ALSO: | diag |