Difference between revisions of "SingularValueDecomp"

m (adding doc status category)
(11 intermediate revisions by 3 users not shown)
Line 1: Line 1:
 
[[category:Matrix Functions]]
 
[[category:Matrix Functions]]
[[Category:Doc Status C]] <!-- For Lumina use, do not change -->
 
 
Computes the singular value decomposition of a matrix.
 
  
{{stub}}
+
== SingularValueDecomp(a, i, j, j2) ==
 +
 
 +
[[SingularValueDecomp]] computes the singular value decomposition of a matrix. Singular value decomposition is often used with sets of equations or matrices that are singular or ill-conditioned (that is, very close to singular). It factors a matrix «a», indexed by «i» and «j», with ''[[IndexLength]](i) >= [[IndexLength]](i)'', into three matrices, ''U'', ''W'', and ''V'', such that:
 +
:a = U . W . V<sup>T</sup>
 +
 
 +
where ''U'' and ''V'' are orthogonal matrices and ''W'' is a diagonal matrix. ''U'' is dimensioned by «i» and «j», ''W'' by «j» and «j2», and ''V'' by «j» and «j2». In Analytica notation:
 +
 
 +
:<code>Variable A := [[Sum]]([[Sum]](U*W, J)*[[Transpose]](V, J, J2), J2)</code>
 +
 
 +
The index «j2» must be the same size as «j» and is used to index the resulting ''W'' and ''V'' arrays.  [[SingularValueDecomp]] returns an array of three elements indexed by a special system index named <code>SvdIndex</code> with each element, ''U'', ''W'', and ''V'', being a reference to the corresponding array.
 +
 
 +
Use the [[Using References|# (dereference) operator]] to obtain the matrix value from each reference, as in:
 +
 
 +
:<code>Index J2 := [[CopyIndex]](J)</code>
 +
:<code>Variable SvdResult := [[SingularValueDecomp]](A, I, J, J2)</code>
 +
:<code>Variable U := #SvdResult[SvdIndex = 'U']</code>
 +
:<code>Variable W := #SvdResult[SvdIndex = 'W']</code>
 +
:<code>Variable V := #SvdResult[SvdIndex = 'V']</code>
 +
 
 +
== Matrix inverse ==
 +
The inverse of a square matrix A, in Analytica syntax, is
 +
<code>
 +
:[[Local]] Winv := [[If]] J=J2 [[Then]] 1/W [[Else]] 0;
 +
:[[Transpose]]([[Sum]]([[Sum]](U*Winv, J)*[[Transpose]](V, J, J2), J2),I,J)
 +
</code>
 +
 
 +
Singular value decomposition can be used for matrix inverse when the matrix A is ill-conditioned, in which case the [[Invert]] function may encounter numeric instabilities. When the matrix is ill-conditioned (the [[Determinant]] is very close to zero), then some of the elements of the diagonal of <code>W</code> will be very close to zero. To avoid the numerical instabilities, the diagonal entries corresponding to the very small <code>W</code> can be replaced with 0 in <code>Winv</code>:
 +
 
 +
<code>
 +
:[[Local]] Winv := [[If]] J=J2 [[And]] [[Abs]](W)>1e-4 [[Then]] 1/W [[Else]] 0;
 +
:[[Transpose]]([[Sum]]([[Sum]](U*Winv, J)*[[Transpose]](V, J, J2), J2),I,J)
 +
</code>
 +
 
 +
To make this convenient to use, you can introduce a new [[User-Defined Function]] as follows:
 +
 
 +
:'''Function''' <code>MatInvert( A : [I,J] ; I,J : Index )</code>
 +
:'''Definition:'''<code>
 +
::Index J2 := J;
 +
::Local svd := SingularValueDecomp(A,I,J,J2);
 +
::Local U := #svd[SvdIndex='U'];
 +
::Local W := #svd[SvdIndex='W'];
 +
::Local Winv := if J=J2 And W>1e-5 then 1/W else 0;
 +
::Local V := #svd[SvdIndex='V'];
 +
::Transpose(Sum(Sum(U*Winv, J)*Transpose(V, J, J2), J2),I,J)
 +
</code>
 +
 
 +
You can then use <code>MatInvert(A,I,J)</code> in place of <code>[[Invert]](A,I,J)</code>.
 +
 
 +
== See Also ==
 +
* [[EigenDecomp]]
 +
* [[Decompose]]
 +
* [[Transpose]]
 +
* [[MatrixMultiply]]
 +
* [[Matrix functions]]
 +
* [[:Category:Matrix Functions]]

Revision as of 18:47, 24 May 2021


SingularValueDecomp(a, i, j, j2)

SingularValueDecomp computes the singular value decomposition of a matrix. Singular value decomposition is often used with sets of equations or matrices that are singular or ill-conditioned (that is, very close to singular). It factors a matrix «a», indexed by «i» and «j», with IndexLength(i) >= IndexLength(i), into three matrices, U, W, and V, such that:

a = U . W . VT

where U and V are orthogonal matrices and W is a diagonal matrix. U is dimensioned by «i» and «j», W by «j» and «j2», and V by «j» and «j2». In Analytica notation:

Variable A := Sum(Sum(U*W, J)*Transpose(V, J, J2), J2)

The index «j2» must be the same size as «j» and is used to index the resulting W and V arrays. SingularValueDecomp returns an array of three elements indexed by a special system index named SvdIndex with each element, U, W, and V, being a reference to the corresponding array.

Use the # (dereference) operator to obtain the matrix value from each reference, as in:

Index J2 := CopyIndex(J)
Variable SvdResult := SingularValueDecomp(A, I, J, J2)
Variable U := #SvdResult[SvdIndex = 'U']
Variable W := #SvdResult[SvdIndex = 'W']
Variable V := #SvdResult[SvdIndex = 'V']

Matrix inverse

The inverse of a square matrix A, in Analytica syntax, is

Local Winv := If J=J2 Then 1/W Else 0;
Transpose(Sum(Sum(U*Winv, J)*Transpose(V, J, J2), J2),I,J)

Singular value decomposition can be used for matrix inverse when the matrix A is ill-conditioned, in which case the Invert function may encounter numeric instabilities. When the matrix is ill-conditioned (the Determinant is very close to zero), then some of the elements of the diagonal of W will be very close to zero. To avoid the numerical instabilities, the diagonal entries corresponding to the very small W can be replaced with 0 in Winv:

Local Winv := If J=J2 And Abs(W)>1e-4 Then 1/W Else 0;
Transpose(Sum(Sum(U*Winv, J)*Transpose(V, J, J2), J2),I,J)

To make this convenient to use, you can introduce a new User-Defined Function as follows:

Function MatInvert( A : [I,J] ; I,J : Index )
Definition:
Index J2 := J;
Local svd := SingularValueDecomp(A,I,J,J2);
Local U := #svd[SvdIndex='U'];
Local W := #svd[SvdIndex='W'];
Local Winv := if J=J2 And W>1e-5 then 1/W else 0;
Local V := #svd[SvdIndex='V'];
Transpose(Sum(Sum(U*Winv, J)*Transpose(V, J, J2), J2),I,J)

You can then use MatInvert(A,I,J) in place of Invert(A,I,J).

See Also

Comments


You are not allowed to post comments.