Symmetrizer Matrix
The symmetrizer matrix $\mathbf{S}$ is a matrix that symmetrizes tensors or higher-order Kronecker products. It ensures that the resulting product is symmetric with respect to its indices. In the context of polynomial dynamical systems and Kronecker products, the symmetrizer matrix is crucial for handling symmetric terms efficiently and accurately.
Definition
The Symmetrizer Matrix $\mathbf{S}_{n,p} \in \mathbb{R}^{n^p \times n^p}$ of order $p$ for a vector of length $n$ is defined such that it symmetrizes a tensor obtained from the $p$-th Kronecker power of a vector $\mathbf{x} \in \mathbb{R}^n$.
Given the $p$-th Kronecker power:
\[\mathbf{x}^{[p]} = \underbrace{\mathbf{x} \otimes \mathbf{x} \otimes \cdots \otimes \mathbf{x}}_{p \text{ times}} \in \mathbb{R}^{n^p}\]
The symmetrizer matrix $\mathbf{S}_{n,p}$ symmetrizes $\mathbf{x}^{[p]}$ to produce a vector where all permutations of indices are averaged. The symmetrized vector is:
\[\mathbf{x}^{\{p\}} = \mathbf{S}_{n,p} \mathbf{x}^{[p]}\]
Properties of the Symmetrizer Matrix:
- Idempotent: $\mathbf{S}_{n,p}^2 = \mathbf{S}_{n,p}$.
- Symmetric: $\mathbf{S}_{n,p}^\top = \mathbf{S}_{n,p}$.
- Projection Matrix: It projects any tensor onto the space of symmetric tensors.
For the quadratic case ($p = 2$), the symmetrizer matrix can be explicitly defined using the commutation matrix $\mathbf{K}_{n,n}$:
\[\mathbf{S}_{n,2} = \frac{1}{2} \left( \mathbf{I}_{n^2} + \mathbf{K}_{n,n} \right)\]
where $\mathbf{I}_{n^2}$ is the identity matrix of size $n^2 \times n^2$.
UniqueKronecker.symmtzrmat
— Functionsymmtzrmat(n::Int, p::Int) -> Sp::SparseMatrixCSC{Float64}
Create a symmetrizer matrix of order p
for a vector of length n
Magnus and Neudecker [1].
Arguments
n::Int
: The length of the vector.p::Int
: The order of the symmetrizer matrix, e.g.,p = 2
for x ⊗ x.
Output
Sp::SparseMatrixCSC{Float64}
: The symmetrizer matrix of orderp
.
Example
julia> symmtzrmat(2,2)
4×4 SparseMatrixCSC{Float64, Int64} with 6 stored entries:
1.0 ⋅ ⋅ ⋅
⋅ 0.5 0.5 ⋅
⋅ 0.5 0.5 ⋅
⋅ ⋅ ⋅ 1.0
Use Cases
The symmetrizer matrix is particularly useful in the following contexts:
Symmetrizing Kronecker Products:
When dealing with polynomial terms in dynamical systems, it's often necessary to ensure that the terms are symmetric. The symmetrizer matrix averages all permutations of indices in the Kronecker product to produce a symmetric tensor.
Example:
For $\mathbf{x} \in \mathbb{R}^n$ and $p = 2$:
\[\mathbf{x}^{\{2\}} = \mathbf{S}_{n,2} (\mathbf{x} \otimes \mathbf{x})\]
This results in a vector containing terms like $x_i x_j$ averaged over all permutations.
Ensuring Symmetric Operators:
In polynomial dynamical systems, operators (matrices) that act on Kronecker powers of state vectors should often be symmetric to preserve certain properties like stability. The symmetrizer matrix can be used to enforce symmetry in these operators.
Duplicating Symmetric Coefficients:
When converting between unique and redundant representations of polynomial operators, the symmetrizer matrix ensures that coefficients are duplicated symmetrically, maintaining the symmetry of the operator.
UniqueKronecker.duplicate_symmetric
— Functionduplicate_symmetric(A::AbstractArray, p::Int)
Duplicate the redundant polynomial coefficients in the matrix A
with a unique set of coefficients and return the matrix with redundant coefficients which are duplicated symmetrically. This guarantees that the operator is symmetric. The difference from duplicate
is that we use the elimination matrix Lp
and the symmetric commutation matrix Sp
to multiply the A
matrix.
Arguments
A::AbstractArray
: A matrixp::Int
: The order of the polynomial, e.g.,p = 2
for x ⊗ x.
Returns
- matrix with redundant coefficients duplicated symmetrically
Example
julia> n = 2; P = rand(n,n); P *= P'; p = vec(P)
4-element Vector{Float64}:
0.5085988756090203
0.7704767970682769
0.7704767970682769
1.310279680309927
julia> Q = rand(n,n); Q *= Q'; q = vec(Q)
4-element Vector{Float64}:
0.40940214810208353
0.2295272821417254
0.2295272821417254
0.25503767587483905
julia> A2 = [p'; q']
2×4 Matrix{Float64}:
0.257622 0.44721 0.0202203 0.247649
1.55077 0.871029 0.958499 0.650717
julia> D2 = dupmat(2,2)
4×3 SparseMatrixCSC{Int64, Int64} with 4 stored entries:
1 ⋅ ⋅
⋅ 1 ⋅
⋅ 1 ⋅
⋅ ⋅ 1
julia> A2 * D2
2×3 Matrix{Float64}:
0.508599 1.54095 1.31028
0.409402 0.459055 0.255038
julia> duplicate_symmetric(A2 * D2, 2)
2×4 Matrix{Float64}:
0.508599 0.770477 0.770477 1.31028
0.409402 0.229527 0.229527 0.255038