Merge pull request #122 from Marco-Congedo/dev

implemented `frf` and `invfrf` (whitening) functions

docs/src/linearAlgebra.md

@@ -141,6 +141,8 @@ congruence
 | Function   | Description |
 |:---------- |:----------- |
 | [`evd`](@ref) | Eigenvalue-Eigenvector decomposition of a matrix in ``UΛU'=P`` form|
+| [`frf`](@ref) | Full-rank factorization of an Hermitian matrix |
+| [`invfrf`](@ref) | Inverse of the full-rank factorization of an Hermitian matrix (whitening) |
 | [`spectralFunctions`](@ref) | Mother function for creating spectral functions of eigenvalues|
 | [`pow`](@ref)| Power of a positive matrix for any number of exponents in one pass|
 | [`invsqrt`](@ref)| Principal square root inverse (whitening) of a positive matrix|
@@ -151,6 +153,8 @@ congruence
 
 ```@docs
 evd
+frf
+invfrf
 spectralFunctions
 pow
 invsqrt

src/PosDefManifold.jl

@@ -137,6 +137,8 @@ export
     fVec,
     congruence, cong,
     evd,
+    frf,
+    invfrf,
     spectralFunctions,
     ispos,
     pow,

src/linearAlgebra.jl

@@ -1451,6 +1451,79 @@ function evd(S::Union{𝕄{T}, ℍ{T}}) where T<:RealOrComplex # return tuple (
 end
 
 
+
+"""
+```
+frf(P::ℍ{T}) where T<:RealOrComplex
+```
+Full-rank factorization of `Hermitian` matrix `P`.
+It is given by
+
+``F=UD^{1/2}``,
+
+where
+
+``\\textrm{EVD}(P)=UDU^{H}``
+
+is the eigenvalue-eigenvector decomposition of `P`. It verifies
+
+``FF^H=P``,
+
+thus ``F^{-1}`` is a whitening matrix.
+
+**See also**: [`invfrf`](@ref).
+
+## Examples
+```
+using LinearAlgebra, PosDefManifold
+P=randP(3)
+F = frf(P)
+F*F'≈P ? println(" ⭐ ") : println(" ⛔ ")
+```
+"""
+function frf(P::ℍ{T}) where T<:RealOrComplex
+   #size(P, 1)≠size(A, 2) && throw(ArgumentError(📌*", frf function: input matrix must be square"))
+   λ, U=eigen(P)
+   return U*Diagonal(sqrt.(λ))
+end
+
+
+"""
+```
+invfrf(P::ℍ{T}) where T<:RealOrComplex
+```
+Inverse of the full-rank factorization of `Hermitian` matrix `P`.
+It is given by
+
+``F=D^{-1/2}U^H``,
+
+where
+
+``\\textrm{EVD}(P)=UDU^{H}``
+
+is the eigenvalue-eigenvector decomposition of `P`. It verifies
+
+``FPF^H=I``,
+
+thus ``F`` is a whitening matrix.
+
+**See also**: [`frf`](@ref).
+
+## Examples
+```
+using LinearAlgebra, PosDefManifold
+P=randP(3)
+F = invfrf(P)
+F*P*F'≈I ? println(" ⭐ ") : println(" ⛔ ")
+```
+"""
+function invfrf(P::ℍ{T}) where T<:RealOrComplex
+   #size(P, 1)≠size(A, 2) && throw(ArgumentError(📌*", frf function: input matrix must be square"))
+   Λ, U=evd(P)
+   return invsqrt(Λ)*U'
+end
+
+
 """
     (1) spectralFunctions(P::ℍ{T}, func) where T<:RealOrComplex
     (2) spectralFunctions(D::𝔻{S}, func) where S<:Real

src/test.jl

@@ -280,6 +280,18 @@ function tests();
     U*Λ*U'≈PC ? OK() : OH(name*" Complex Input")
 
 
+    name="frf"; newTest(name)
+    F = frf(P)
+    F*F'≈P ? OK() : OH(name*" Real Input")
+    F = frf(PC)
+    F*F'≈PC ? OK() : OH(name*" Complex Input")
+
+    name="invfrf"; newTest(name)
+    F = invfrf(P)
+    F*P*F'≈I ? OK() : OH(name*" Real Input")
+    F = invfrf(PC)
+    F*PC*F'≈I ? OK() : OH(name*" Complex Input")
+
     name="spectralFunctions"; newTest(name);
     spectralFunctions(P, x->x+1); RUN()
     spectralFunctions(PC, abs2); RUN()