Representations

Representation in RepLAB, whether they are exact or inexact, reducible or irreducible, are described by the Rep base class.

  • Rep is a real or complex representation of a compact group.

Some representations encode a division algebra using real or complex coefficients.

  • DivisionAlgebra provides information about the division algebra encoding handled by RepLAB.

A particular type of representation is a representation of a finite group given by generator images.

  • RepByImages is a real or complex representation described by the (matrix) images of the group generators.

Once a Rep has been constructed, it can be decomposed into irreducible representations using the rep.decomposition method.

  • SubRep describes a subrepresentation of an existing representation. The instance can remember which representation it splits (SubRep.parent) and the change of basis maps (SubRep.injection, SubRep.projection).

  • Isotypic are isotypic components, which group equivalent irreps present in a representation.

  • Irreducible represents the decomposition of a representation into isotypic components.

Rep

class replab.Rep(group, field, dimension, varargin)

Bases: replab.Obj

Describes a finite dimensional representation of a compact group

This class has properties that correspond to information that can be computed and cached after the replab.Rep instance is constructed, for example isUnitary or trivialDimension.

Notes

While we do not expect users to implement their own subclass of Rep, to do so only the map from group elements to matrix images need to be implemented.

Either the override the image_double_sparse method if the user implementation provides floating-point images, or both the isExact and image_exact methods.

The method computeErrorBound must also be overriden (or the value cached at construction).

This class implements also extra methods about the action of the representation on matrices, methods which can be overloaded for performance.

If this class has representations as “parts”, override the decomposeTerm and composeTerm methods.

Own members

dimension

Representation dimension

Type

integer

divisionAlgebraName

Name of the division algebra encoded by this representation (see replab.DivisionAlgebra)

Type

‘R->C’, ‘C->R’, ‘H->C’, ‘H->R:rep’, ‘’

field

Vector space defined on real (R) or complex (C) field

Type

{‘R’, ‘C’}

group

Group being represented

Type

replab.CompactGroup

isUnitary

Whether the representation is unitary

Type

logical

Rep(group, field, dimension, varargin)

Constructs a finite dimension representation

Parameters
  • group (replab.CompactGroup) – Group being represented

  • field ({'R', 'C'}) – Whether the representation is real (R) or complex (C)

  • dimension (integer) – Representation dimension

  • isUnitary – Whether the representation is unitary, default: false

  • 'H->R (divisionAlgebraName % ('R->C', 'C->R', 'H->C',) – rep’, ‘’, optional): Name of the division algebra encoding this representation respects, default ''

  • isIrreducible – Whether this representation is irreducible, default []

  • trivialDimension – Dimension of the trivial subrepresentation, default []

  • frobeniusSchurIndicator – Exact value of the Frobenius-Schur indicator, default []

Kwtype isUnitary

logical, optional

Kwtype isIrreducible

logical or [], optional

Kwtype trivialDimension

integer or [], optional

Kwtype frobeniusSchurIndicator

integer or [], optional

alternatingSquare(self)

Returns the alternating square of this representation

Let \(V\) be the vector space corresponding to this representation, and \(V \otimes V\) be its second tensor power. Let \(T\) be the linear map such that \(T(v_1 \otimes v_2) = v_2 \otimes v_1\).

Then the alternating square \(A\) is the subspace \(A subset V \otimes V\) such that \(T(a) = -a\) for all \(a \in A\).

Returns

A subrepresentation of the second tensor power of this representation

Return type

SubRep

antilinearInvariant(self, type)

Returns the equivariant space of matrices representing equivariant antilinear maps

Let F be an antilinear map such that F(alpha * x) = conj(alpha) * F(x) for any scalar alpha and vector x. We describe F using a matrix J such that F(x) = J * conj(x).

The following conditions are equivalent.

  • F is an equivariant map: F(rho.image(g) * x) == rho.image(g) * F(x)

  • J is an equivariant matrix: rho.image(g) * J == J * conj(rho.image(g))

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The space of matrices describing equivariant antilinear maps

Return type

replab.Equivariant

bilinearInvariant(self, type)

Returns the equivariant space of matrices representing equivariant bilinear maps

Let F be an equivariant bilinear map, i.e. F(x,y) = F(rho(g) x, rho(g) y). We describe this map by a matrix X such that F(x,y) = x.' * X * y. When the map F is equivariant, x.' * X * y = x' * rho(g).' * X * rho(g) * y for all vectors, so that rho(g).' * X = X * rho(g).

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The space of matrices describing equivariant bilinear maps

Return type

replab.Equivariant

blkdiag(varargin)

Direct sum of representations

See replab.CompactGroup.directSumRep

blockSubRep(self, block)

Constructs a subrepresentation from a subset of Euclidean coordinates

The block argument can be obtained from invariantBlocks, and this method constructs the sparse injection and projection maps with which to call subRep.

Parameters

block (integer(1,*)) – Euclidean coordinates of the block

Returns

Subrepresentation

Return type

SubRep

character(self)

Returns the character corresponding to this representation, provided the group represented is finite

If this representation is approximate, a trick from Dixon is used to compute the corresponding exact character.

Returns

Character corresponding to this representation

Return type

Character

commutant(self, type)

Returns the commutant of this representation

This is the algebra of matrices that commute with the representation, i.e. the vector space isomorphism to the equivariant space from this rep to this rep.

For any g in G, we have rho(g) * X = X * rho(g).

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The commutant algebra represented as an equivariant space

Return type

replab.Equivariant

complexification(self)

Returns the complexification of a real representation

Returns

The complexification of this representation

Return type

replab.Rep

Raises

An error if this representation is already complex.

compose(self, applyFirst)

Composition of a representation with a morphism, with the morphism applied first

Parameters

applyFirst (Morphism) – Morphism from a finite group

Returns

Representation

Return type

Rep

conditionNumberEstimate(self)

Returns an estimation of the condition number of the change of basis that makes this representation unitary

Returns

Condition number of the change of basis matrix

Return type

double

conj(self)

Returns the complex conjugate representation of this representation

See https://en.wikipedia.org/wiki/Complex_conjugate_representation

It obeys rep.conj.image(g) = conj(rep.image(g))

If this representation is real, it is returned unchanged.

Returns

The complex conjugate of this representation

Return type

replab.Rep

contramap(self, morphism)

Returns the representation composed with the given morphism applied first

Parameters

morphism (replab.FiniteMorphism) – Morphism of finite groups such that morphism.target == self.group

Returns

Representation on the finite group morphism.source

Return type

replab.Rep

decomposition(self, type)

Returns the irreducible decomposition of this representation

The type 'double/sparse' is the same as double.

Exact decomposition depends on the availability of an algorithm: for example, we decompose finite groups with an available character table using the Serre projection formulas.

The approximate decomposition (type = ‘double’) is not deterministic in all aspects as it uses random samples from the commutant equivariant space. But there is also a catch with the exact decomposition algorithm: it retrieves the character table of the group from the atlas, by finding an isomorphism to this representation’s group. But this isomorphism is unique only up to automorphism of group.

Parameters

type ('double', 'double/sparse', 'exact', optional) – Decomposition type, default: double

Returns

The irreducible decomposition

Return type

replab.Irreducible

directSumOfCopies(self, n)

Returns a direct sum of copies of this representation

Parameters

n (integer) – Number of copies

Returns

The direct sum representation

Return type

replab.Rep

dual(self)

Returns the dual representation of this representation

See https://en.wikipedia.org/wiki/Dual_representation

It obeys rep.dual.image(g) = rep.inverseImage(g).'

Returns

The dual representation

Return type

replab.Rep

equivariantFrom(self, repC, varargin)

Returns the space of equivariant linear maps from another rep to this rep

The equivariant vector space contains the matrices X such that

self.image(g) * X == X * repC.image(g)

Parameters
  • repC (replab.Rep) – Representation on the source/column space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

Kwtype special

charstring or ‘’, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Returns

The equivariant space

Return type

replab.Equivariant

equivariantTo(self, repR, varargin)

Returns the space of equivariant linear maps from this rep to another rep

The equivariant vector space contains the matrices X such that

repR.image(g) * X == X * self.image(g)

Parameters
  • repR (replab.Rep) – Representation on the target/row space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

Kwtype special

charstring or ‘’, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Returns

The equivariant space

Return type

replab.Equivariant

errorBound(self)

Returns a bound on the approximation error of this representation

Returns

There exists an exact representation repE such that ||rep(g) - repE(g)||fro <= errorBound

Return type

double

frobeniusSchurIndicator(self)

Returns the Frobenius-Schur indicator of this representation

It is an integer corresponding to the value \(\iota = \int_{g \in G} tr[\rho_g^2] d \mu\) or \(\iota = \frac{1}{|G|} \sum_{g \in G} tr[\rho_g^2]\).

For real irreducible representations, the Frobenius-Schur indicator can take values:

  • 1 if the representation is of real-type; its complexification is then also irreducible

  • 0 if the representation is of complex-type; it decomposes into two conjugate irreducible representations over the complex numbers

  • -2 if the representation is of quaternion-type.

For complex irreducible representations, it can take the values:

  • 1 if the representation is of real-type,

  • 0 if the representation is of complex-type, i.e. is not equivalent to its conjugate,

  • -1 if the representation is of quaternion-type.

Returns

Value of the indicator

Return type

integer

hasTorusImage(self)

Returns whether a description of representation of the maximal torus subgroup is available

This description is used internally in RepLAB to speed up the group averaging process when computing equivariants.

Returns

True if the call to torusImage succeeds

Return type

logical

hermitianInvariant(self, type)

Returns the equivariant space of matrices representing equivariant Hermitian sesquilinear maps

Similar sesquilinearInvariant with the additional constraint that the matrix X == X'.

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The space of matrices describing equivariant Hermitian sesquilinear maps

Return type

replab.Equivariant

identifyIrrep(self)

Identifies the irrep type and returns an equivalent representation with explicit structure

This representation must be irreducible.

Example

>>> Q = replab.QuaternionGroup();
>>> rep = Q.naturalRep.complexification;
>>> irreps = rep.split;
>>> assert(length(irreps) == 2); % splits in two equivalent irreps
>>> sub = irreps{1};
>>> sub1 = sub.identifyIrrep;
>>> strcmp(sub1.divisionAlgebraName, 'H->C')
    1
Returns

Similar representation with the canonical division algebra encoding

Return type

replab.SubRep

image(self, g, type)

Returns the image of a group element

Raises

An error if type is 'exact' and isExact is false.

Parameters
  • g (element of group) – Element being represented

  • type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

Image of the given element for this representation

Return type

double(*,*) or cyclotomic(*,*)

imap(self, f)

Maps the representation under an isomorphism

Parameters

f (FiniteIsomorphism) – Isomorphism with self.group.isSubgroupOf(f.source)

Returns

Representation satisfying newRep.group.isSubgroup(f.target).

Return type

Rep

invariantBlocks(self)

Returns a partition of the set 1:self.dimension such that the subsets of coordinates correspond to invariant spaces

This method does not guarantee that the finest partition is returned.

Example

>>> G = replab.SignedPermutationGroup.of([1 4 7 2 5 8 3 6 9], [1 -2 -3 -4 5 6 -7 8 9], [1 3 2 4 6 5 -7 -9 -8]);
>>> rep = G.naturalRep;
>>> p = rep.invariantBlocks;
>>> isequal(p.blocks, {[1] [2 3 4 7] [5 6 8 9]})
    1
Returns

Partition of coarse invariant blocks

Return type

Partition

inverseImage(self, g, type)

Returns the image of the inverse of a group element

Parameters
  • g (element of group) – Element of which the inverse is represented

  • type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

Image of the inverse of the given element for this representation

Return type

double(*,*) or cyclotomic(*,*)

isExact(self)

Returns whether this representation can compute exact images

Returns

True if the call self.image(g, 'exact') always succeeds

Return type

logical

isIrreducible(self)

Returns whether this representation is irreducible over its current field

Note that a real representation may become reducible once it is computed over the complex numbers.

Returns

True if this representation is irreducible, false if it has a nontrivial subrepresentation

Return type

logical

isIrreducibleAndCanonical(self)

Returns whether this representation is irreducible, and has its division algebra in the canonical encoding

This is always true for complex irreps (though we may decide to canonicalize quaternion-type complex representations later). For real irreps, it is true if the Frobenius-Schur indicator is known and matches the encoded divisionAlgebraName.

Returns

True if the representation is irreducible and has its division algebra in the canonical encoding

Return type

logical

kernel(self)

Returns the kernel of the given representation

Only works if group is a finite group.

Raises

An error if this representation is not precise enough to compute the kernel.

Example

>>> S3 = replab.S(3);
>>> K = S3.signRep.kernel;
>>> K == replab.PermutationGroup.alternating(3)
    1
Returns

The group K such that rho.image(k) == id for all k in K

Return type

replab.FiniteGroup

knownIrreducible(self)

Returns whether this representation is known to be irreducible; only a true result is significant

Returns

True if isIrreducible is known and is true

Return type

logical

knownProperties(self, keys)

Returns the known properties among the set of given keys as a key-value cell array

Parameters

cell (1,*) – Property names

Returns

Key-value pairs

Return type

cell(1,*)

knownReducible(self)

Returns whether this representation is known to be reducible; only a true result is significant

Returns

True if isIrreducible is known and is false

Return type

logical

kron(varargin)

Tensor product of representations

See replab.CompactGroup.tensorRep

maschke(self, sub1)

Given a subrepresentation, returns the complementary subrepresentation

Note that we optimize special cases when the representation and/or the injection is unitary

Parameters

sub1 (SubRep) – Subrepresentation of this representation

Returns

Complementary subrepresentation

Return type

SubRep

matrixColAction(self, g, M, type)

Computes the matrix-representation product

We multiply by the inverse of the image, so this stays a left action.

self.matrixColAction(g, self.matrixColAction(h, M))

is the same as

self.matrixColAction(self.group.compose(g, h), M)

Parameters
  • g (group element) – Group element acting

  • M (double(*,*)) – Matrix acted upon

  • type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’.

Returns

The matrix M * self.image(inverse(g))

Return type

double(*,*)

matrixRowAction(self, g, M, type)

Computes the representation-matrix product

This is a left action:

self.matrixRowAction(g, self.matrixRowAction(h, M))

is the same as

self.matrixRowAction(self.group.compose(g, h), M)

Parameters
  • g (group element) – Group element acting

  • M (double(*,*) or cyclotomic (*,*)) – Matrix acted upon

  • type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’.

Returns

The matrix self.image(g) * M

Return type

double(*,*) or cyclotomic (*,*)

overC(self)

Returns whether this representation is defined over the complex field

Returns

True if this representation is defined over the complex field

Return type

logical

overR(self)

Returns whether this representation is defined over the real field

Returns

True if this representation is defined over the real field

Return type

logical

permutedDirectSumOfCopies(self, mu)

Returns a direct sum of copies of this representation

Example

>>> U2 = replab.U(2);
>>> S3 = replab.S(3);
>>> G = U2.directProduct(S3);
>>> GtoU2 = G.projection(1); % morphism from G to U2
>>> GtoS3 = G.projection(2); % morphism from G to S3
>>> repU2 = GtoU2.andThen(U2.definingRep);
>>> rep = repU2.permutedDirectSumOfCopies(GtoS3);
>>> rep.dimension
    6
Parameters

mu (Morphism) – Morphism from group to a PermutationGroup

Returns

The permuted tensor power representation

Return type

Rep

permutedTensorPower(self, mu)

Returns a tensor power of this representation, whose factors are permuted by the group

Example

>>> U2 = replab.U(2);
>>> S3 = replab.S(3);
>>> G = U2.directProduct(S3);
>>> GtoU2 = G.projection(1); % morphism from G to U2
>>> GtoS3 = G.projection(2); % morphism from G to S3
>>> repU2 = GtoU2.andThen(U2.definingRep);
>>> rep = repU2.permutedTensorPower(GtoS3);
>>> rep.dimension
    8
Parameters

mu (Morphism) – Morphism from group to a PermutationGroup

Returns

The permuted tensor power representation

Return type

Rep

restrictedTo(self, subgroup)

Returns the restricted representation to the given subgroup

Parameters

subgroup (CompactGroup) – Subgroup to restrict the representation to, must be a subgroup of group

sample(self, type)

Samples an element from the group and returns its image under the representation

Raises

An error if type is 'exact' and isExact is false.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Optionally, the inverse of the image can be returned as well.

Returns

  • rho (double(*,*)) – Image of the random group element

  • rhoInverse (double(*,*)) – Inverse of image rho

sesquilinearInvariant(self, type)

Returns the space of matrices representing equivariant sesquilinear maps

Let F be an equivariant sesquilinear map, i.e. F(x,y) = F(rho(g) x, rho(g) y). We describe this map by a matrix X such that F(x,y) = x' * X * y. When the map F is equivariant, x' * X * y = x' * rho(g)' * X * rho(g) * y for all vectors, so that rho(g)' * X = X * rho(g).

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The equivariant space

Return type

replab.Equivariant

similarRep(self, A, varargin)

Returns a similar representation under a change of basis

It returns a representation rep1 such that

rep1.image(g) = A * self.image(g) * inv(A).

Additional keyword arguments will be passed on to the Rep constructor. The argument isUnitary must be provided, except when the inverse is provided and both this representation and the change of basis is unitary.

Parameters

A (double(*,*) or cyclotomic (*,*)) – Change of basis matrix

Keywords Args:

inverse (double(*,*) or cyclotomic (*,*)): Inverse of the change of basis matrix

Returns

The similar representation

Return type

replab.SubRep

simplify(self, varargin)

Returns a representation identical to this, but possibly with its structure simplified

It pushes SubRep to the outer level, and replab.rep.DerivedRep to the inner levels; expands tensor products.

Additional keyword arguments can be provided as key-value pairs.

Parameters
  • dense – Whether to allow multiplication of non-sparse matrices

  • approximate – Whether to allow multiplication of non-exact matrices, implies dense = true

Kwtype dense

logical, optional

Kwtype approximate

logical, optional

Returns

A possibly simplified representation

Return type

replab.Rep

split(self, varargin)

Splits this representation into irreducible subrepresentations

Parameters

forceNonUnitaryAlgorithms – Whether to force the use of algorithms for not necessarily unitary representations, default: false

Kwtype forceNonUnitaryAlgorithms

logical, optional

subRep(self, injection, varargin)

Constructs a subrepresentation of this representation

A subrepresentation is specified by an invariant subspace of the current representation.

For simplicity, this subspace can be provided as column vectors in a matrix passed as the first argument to this function. This matrix is called injection for reasons that are clarified below.

The user should also provide the keyword argument isUnitary, except when this representation is unitary, and the injection map corresponds to an isometry; then RepLAB deduces that the subrepresentation is unitary as well.

More precisely, the subrepresentation is defined by maps between two vector spaces:

  • the vector space \(V\) corresponding to this representation, of dimension \(D\),

  • the vector space \(W\) corresponding to the created subrepresentation, of dimension \(d\).

The two maps are:

  • the injection \(I: W \rightarrow V\) that takes a vector in \(W\) and injects it in the parent representation,

  • the projection \(P: V \rightarrow W\) that takes a vector in \(V\) and extracts/returns its component in \(W\).

Correspondingly, the map \(I\) is given as a \(D x d\) matrix, while the map \(P\) is given as a \(d x D\) matrix.

If only the injection map \(I\) is given as an argument, a projection \(P\) is computed; it is not necessarily unique. In particular, if the range of \(I\) spans irreducible representations with non-trivial multiplicities, the recovered projection is chosen arbitrarily.

However, all choices of the projection map provide the same results when computing images of the subrepresentation.

If a floating-point approximation \(\tilde{I}\) of the injection map is given instead of the exact map \(I\), an upper bound on the error can be provided; otherwise, it will be automatically estimated. In that case, the projection map \(\tilde{P}\) is also considered as approximate, regardless of whether it is user-provided or computed.

The error bound is mapErrorBound, and corresponds to an upper bound on both \(|| I \tilde{P} - id ||_F\) and \(|| \tilde{I} P - id ||_F\); in the expression, we assume that \(I\) and \(P\) are the exact injections/projections closest to \(\tilde{I}\) and \(\tilde{P}\).

Parameters
  • injection (double(D,d) or cyclotomic (D,d), may be sparse) – Basis / Injection map

  • projection – Projection map

  • mapErrorBound – Upper bound as described above

  • mapConditionNumberEstimate – Upper bound on the condition number of both \(P\) and \(I\)

  • isUnitary – Whether the resulting representation is unitary, may be omitted

  • largeScale – Whether to use the large-scale version of the algorithm, default automatic selection

  • tolerances – Termination criteria

  • nSamples – Number of samples to use in the large-scale version of the algorithm, default 5

  • forceReal – Whether to force the computed projection map to be real, default false

Kwtype projection

double(D,d) or cyclotomic (D,d), may be sparse, optional

Kwtype mapErrorBound

double, optional

Kwtype mapConditionNumberEstimate

double, optional

Kwtype isUnitary

logical, optional

Kwtype largeScale

logical, optional

Kwtype tolerances

replab.rep.Tolerances

Kwtype nSamples

integer, optional

Kwtype forceReal

logical, optional

Returns

Subrepresentation

Return type

replab.SubRep

symmetricInvariant(self, type)

Returns the equivariant space of matrices representing equivariant symmetric bilinear maps

Similar bilinearInvariant with the additional constraint that the matrix X == X.'.

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The space of matrices describing equivariant symmetric bilinear maps

Return type

replab.Equivariant

symmetricSquare(self)

Returns the symmetric square of this representation

Let \(V\) be the vector space corresponding to this representation, and \(V \otimes V\) be its second tensor power. Let \(T\) be the linear map such that \(T(v_1 \otimes v_2) = v_2 \otimes v_1\).

Then the symmetric square \(S\) is the subspace \(S subset V \otimes V\) such that \(T(s) = s\) for all \(s \in S\).

Returns

A subrepresentation of the second tensor power of this representation

Return type

SubRep

tensorPower(self, n)

Returns a tensor power of this representation, with the blocks permuted by the group

Parameters

n (integer) – Exponent of the tensor power

Returns

The tensor power representation

Return type

Rep

torusImage(self)

Returns a simple description of the representation of the group maximal torus

Let: * t be an element of the maximal torus T of group, * mu be the morphism from T to group, * torusMap be the integer matrix representing the morphism from T to T1 (a torus of dimension dimension).

We have equality between * self.image(mu.imageElement(t)) and * torusInjection * T1.definingRep.image(torusMap * t) * torusProjection

Returns

  • torusMap (integer(d, r)) – Exponents, with d the representation dimension dimension and r the torus rank

  • torusInjection (double(d, d), may be sparse) – Injection from the torus representation to this representation

  • torusProjection (double(d, d) may be sparse) – Projection from this representation to the torus representation

trivialColSpace(self, type)

Returns an equivariant space to this representation from a trivial representation

The trivial representation has the same dimension as this representation

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The equivariant space

Return type

replab.Equivariant

trivialComponent(self, type)

Returns the trivial isotypic component present in this representation

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The trivial isotypic component

Return type

Isotypic

trivialDimension(self)

Returns the dimension of the trivial subrepresentation in this representation

Returns

Dimension

Return type

integer

trivialProjector(self, type)

Returns the projector into the trivial component present in this representation

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

  • proj (double(*,*) or cyclotomic (*,*)) – Projector

  • err (double) – Estimated error in Frobenius norm

trivialRowSpace(self, type)

Returns an equivariant space to a trivial representation from this representation

The trivial representation has the same dimension as this representation

The computation is cached.

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The equivariant space

Return type

replab.Equivariant

unitarize(self)

Returns a unitary representation equivalent to this representation

The returned representation is of type SubRep, from which the change of basis matrix can be obtained.

If the representation is already unitary, the returned SubRep has the identity matrix as injection/projection maps;

Example

>>> S3 = replab.S(3);
>>> defRep = S3.naturalRep.complexification;
>>> C = randn(3,3) + 1i * rand(3,3);
>>> nonUnitaryRep = defRep.similarRep(C);
>>> unitaryRep = nonUnitaryRep.unitarize;
>>> U = unitaryRep.sample;
>>> norm(U*U' - eye(3)) < 1e-10
    1
Returns

Unitary similar representation

Return type

replab.SubRep

DivisionAlgebra

class replab.DivisionAlgebra(name)

Bases: replab.domain.VectorSpace

Describes an encoding of a division algebra using real or complex matrices

The three division algebras supported are the real numbers (R), complex numbers (C) and quaternions (H). They can be encoded using real matrix blocks (R) or complex matrix blocks (C).

A complex scalar is written s = a + i b where a and b are reals and i is the imaginary unit. A quaternion scalar is written s = a + i b + j c + k d, where i, j , k are the fundamental quaternion units, obeying the relations i^2 = j^2 = k^2 = i j k = -1.

It supports the following encodings:

  • 'R->C': A real scalar is encoded as one complex scalar.

  • 'C->R': A complex scalar is encoded as a 2x2 real matrix block.

  • 'H->C': A quaternion scalar is encoded as a 2x2 complex matrix block.

  • 'H->R:rep' : A quaternion scalar is encoded as a 4x4 real matrix block, using the rep encoding.

  • 'H->R:equivariant' : A quaternion scalar is encoded as a 4x4 real matrix block, using the equivariant encoding.

The 'C->R' encoding represents a complex number \(a + i b\) as:

[ a -b
  b  a ]

It has dimension 2 and matrixSize 2 as well. Note that another encoding could be [a b; -b a]; the encoding we choose seems more widespread but is otherwise arbitrary.

For the quaternions, we need three matrices \(I, J, K\) obeying those relations. If we further require the matrices to be sparse, there are 48 different possible encodings, see R. W. Farebrother, J. Groß, and S.-O. Troschke, “Matrix representation of quaternions”, Linear Algebra and its Applications, vol. 362, pp. 251–255, Mar. 2003 https://doi.org/10.1016/S0024-3795(02)00535-9

Before introducing our choice, let us remark that the quaternion \(q\) can be encoded used the complex matrix:

Quaternion a+i*b+j*c+k*d, encoding using a 2x2 complex matrix (H->C)

  [ a+i*b c+i*d
   -c+i*d a-i*b ]

Now, expanding each complex coefficient using the complex algebra real encoding described above, we obtain:

Quaternion equivariant encoding (H->R:equivariant)

  [  a -b  c -d
     b  a  d  c
    -c -d  a  b
     d -c -b  a  ]

This is the encoding that will be looked for when sampling from quaternionic-type commutant matrices, thus we named it that way.

Now, quaternionic-type representations will have another encoding of the quaternion algebra that commutes with the equivariant encoding. It is:

Quaternion representation encoding (H->R:rep)

  [  a -b -c -d
     b  a -d  c
     c  d  a -b
     d -c  b  a  ]

Own members

basis

Division algebra basis

Type

double(matrixSize,matrixSize,dimension)

dualBasis

Dual abasis under the product (D,B) = trace(D'*B)

Type

double(matrixSize,matrixSize,dimension)

indexMatrix

Index matrix

Type

double(matrixSize,matrix)

name

Name of the division algebra

Type

‘complex’, ‘quaternion.rep’, ‘quaternion.equivariant’

decode(self, X)

Decodes the division algebra elements present in a matrix encoding

This method accepts a matrix of size (ms*d1) x (ms*d2) as argument, where ms is matrixSize and d1, d2 are integers.

Parameters

X (double(ms*d1,ms*d2)) – (Block) matrix to decode

Returns

Decoded elements where d is dimension

Return type

double(d1,d2,d)

encode(self, X)

Takes a matrix of division algebra elements and encodes them in a matrix

This method accepts a tensor of size d1 x d2 x d as argument, where d is dimension and d1, d2 are integers.

Parameters

X (double(d1,d2,d)) – Tensor encoding the matrix

Returns

Encoded matrix

Return type

double(ms*d1, ms*d2)

project(self, X)

Projects the given matrix in this division algebra

This method accepts a matrix of size (ms*d1) x (ms*d2) as argument, where ms is matrixSize and d1, d2 are integers. The projection is then done on each ms x ms block.

Parameters

X (double(ms*d1,ms*d2)) – (Block) matrix to project

Returns

Projected matrix

Return type

double(ms*d1, ms*d2)

RepByImages

class replab.RepByImages(group, field, dimension, preimages, images, imagesErrorBound, varargin)

Bases: replab.Rep

A finite dimensional representation of a finite group

Own members

images

Images

Type

cell(1,n) of double(d,d) or cyclotomic(d,d) or intval(d,d), may be sparse

imagesErrorBound

Error bound on the given images

Type

double(1,n)

preimages

Preimages

Type

cell(1,n) of group elements

static fromExactRep(rep)

Constructs a replab.RepByImages from an existing representation with exact images

rewriteTerm_hasBlockStructure(self, options)

Rewrite rule: rewrite this representation as a direct sum if it has a block structure

rewriteTerm_isTrivial(self, options)

Rewrite rule: replace this representation by a trivial representation if it is a representation of the trivial group

SubRep

class replab.SubRep(parent, injection_internal, projection_internal, varargin)

Bases: replab.Rep

Describes a subrepresentation of a finite dimensional representation

A subrepresentation corresponds to a subspace \(W\) of a parent representation associated with the vector space \(V\).

For simplicity, we say that the subspace \(W\) is specified using a set of basis vectors, amalgamated as column vectors in a matrix injection.

More precisely, the subrepresentation is defined by two maps, an injection and a projection.

  • The injection map \(I: W \rightarrow V\) expresses a vector in \(V\) from a vector in \(W\).

  • The projection map \(P: V \rightarrow W\) extracts the \(W\) part from a vector in the parent space \(V\).

Those maps obey \(P I = id\).

If one wants an operator \(\Pi: V \rightarrow V\) acting as a projector (\(\Pi^2 = \Pi\)), then one simply defines \(\Pi = I P\).

Note that often, only approximations \(\tilde{P}\) and \(\tilde{I}\) are known by RepLAB. We thus store a measure of the error associated with those maps in projectorErrorBound.

Note

The terminology of “injection” and “projection” linear maps was inspired by https://arxiv.org/pdf/0812.4969.pdf

Own members

injection_internal

Injection map

Type

double(D,d) or cyclotomic (D,d), may be sparse

isSimilarRep

True if this SubRep encodes a similarity transformation

Type

logical

mapsAreAdjoint

True if parent is unitary (so the Hermitian adjoint makes sense) and injection is the conjugate transpose of projection

Type

logical

parent

Parent representation of dimension \(D\)

Type

replab.Rep

projection_internal

Projection map

Type

double(d,D) or cyclotomic (d,D), may be sparse

SubRep(parent, injection_internal, projection_internal, varargin)

Constructs a subrepresentation of a parent representation

Additional keyword arguments are passed to the Rep constructor.

If both parent is unitary, and injection_internal is adjoint to projection_internal, then the isUnitary argument may be omitted as it can be deduced to be true.

Parameters
  • parent (replab.Rep) – Parent representation of dimension \(D\)

  • injection_internal (double(D,d) or cyclotomic (D,d), may be sparse) – Injection map \(I\)

  • projection_internal (double(d,D) or cyclotomic (d,D), may be sparse) – Projection map \(P\)

  • isUnitary – Whether the resulting representation is unitary, may be omitted (see above)

  • projectorErrorBound – Upper bound on || I P - tilde{I} tilde{P} ||_F

  • injectionConditionNumberEstimate – Upper bound of the condition number of \(\tilde{I}\) (and thus \(\tilde{P}\))

Kwtype isUnitary

logical, optional

Kwtype projectorErrorBound

double, optional

Kwtype injectionConditionNumberEstimate

double, optional

basis(self, indices, type)

Returns (part of) the basis of the subrepresentation in the parent representation

It is an alias for I(:, indices) where I = self.injection(type).

Parameters
  • indices (integer(1,*) or []) – Indices of the basis elements to return

  • type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The basis of the subrepresentation given as column vectors

Return type

double(*,*) or cyclotomic (*,*)

biorthogonalityErrorBound(self)

Returns an upper bound on the biorthogonality of the injection/projection maps

This returns an upper bound on norm(self.projection*self.injection - eye(d), 'fro'), where d == self.dimension.

collapse(self)

Simplifies a SubRep of a SubRep

Raises

An error if parent is not of type SubRep

Returns

A subrepresentation of .parent.parent

Return type

SubRep

static directSumBiorthogonal(parent, subReps)

Computes the direct sum of subrepresentations of the same parent representation

The subrepresentations must have biorthogonal injection projection maps: subReps{i}.projection * subReps{i}.injection ~= eye(subReps{i}.dimension), but subReps{i}.projection * subReps{j}.injection ~= 0 for i ~= j.

Parameters
  • parent (replab.Rep) – Parent representation

  • subReps (cell(1,*) of replab.SubRep) – Subrepresentations of the parent representation

Returns

A block-diagonal subrepresentation composed of the given subrepresentations

Return type

replab.SubRep

static identical(parent)

Creates a full subrepresentation of the given representation, with identity projection/injection maps

Parameters

parent (replab.Rep) – Representaiton

Returns

Subrepresentation identical to parent

Return type

replab.SubRep

injection(self, type)

Returns the injection from the subrepresentation to the parent representation

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The injection map

Return type

double(*,*) or cyclotomic (*,*)

injectionConditionNumberEstimate(self)

Returns an upper bound on the condition number of both the injection and the projection map

Returns

Upper bound on the condition number

Return type

double

injectionIsExact(self)

Returns whether the injection map is written either with cyclotomics or integer coefficients

Returns

True if that is the case

Return type

logical

mapsAreExact(self)

Returns whether both the injection map and the projection map are written either with cyclotomics or integer coefficients

Returns

True if both injection and projection are cyclotomic matrices or have Gaussian integer entries

Return type

logical

mapsAreIntegerValued(self)

Returns whether both the injection map and the projection map are expressed with Gaussian integer coefficients

Returns

True if both injection and projection have Gaussian integer entries

Return type

logical

nice(self)

Returns a representation similar to the current subrepresentation, with a nicer basis (EXPERIMENTAL)

The “niceness” of the basis is implementation dependent. As of the first implementation of this feature, RepLAB tries to make the basis real, and then with small integer coefficients.

The returned subrepresentation is not necessarily unitary.

In the case no improvement could be made, the original subrepresentation is returned.

Returns

  • res (SubRep) – A subrepresentation of self.parent

  • better (logical) – Whether a nicer subrepresentation has been found

projection(self, type)

Returns the projection from the parent representation to the subrepresentation

Parameters

type ('double', 'double/sparse' or 'exact', optional) – Type of the returned value, default: ‘double’

Returns

The projection map

Return type

double(*,*) or cyclotomic (*,*)

projector(self, type)

Returns the projector on this subrepresentation

Note that the projector is in general not uniquely defined by the injection map; it is defined by the injection map only when the subrepresentation is a direct sum of isotypic components.

Returns

Projector matrix on the subrepresentation

Return type

double(*,*)

projectorErrorBound(self)

Returns an upper bound on the quality of the approximation of the injection/projection maps

Let \(I\) and \(P\) be the respective exact injection and projection maps. Let \(\tilde{I}\) and \(\tilde{P}\) be their approximate counterparts.

This returns an upper bound on the Frobenius norm of \(I P - \tilde{I} \tilde{P}\).

Returns

Upper bound on the Frobenius distance of the approximate projector to the commutant space

Return type

double

refine(self, varargin)

Refines this subrepresentation

Assumes that the subrepresentation is already close to an exact subrepresentation, and refines its subspace by an iterative procedure applied on its injection and projection maps.

Parameters
  • tolerances – Termination criteria

  • largeScale – Whether to use the large-scale version of the algorithm, default automatic choice

  • nSamples – Number of samples to use in the large-scale version of the algorithm, default 5

  • nInnerIterations – Number of inner iterations in the medium-scale version of the algorithm, default 3

  • injectionBiortho – Injection map of known multiplicity space to remove from this subrepresentation

  • projectionBiortho – Projection map of known multiplicity space to remove from this subrepresentation

Kwtype tolerances

replab.rep.Tolerances

Kwtype largeScale

logical, optional

Kwtype nSamples

integer, optional

Kwtype nInnerIterations

integer, optional

Kwtype injectionBiortho

double(*,*), may be sparse

Kwtype projectionBiortho

double(*,*), may be sparse

Returns

sub1 – Subrepresentation with refined subspace (injection/projection maps)

Return type

SubRep

subEquivariantFrom(self, repC, varargin)

Returns the space of equivariant linear maps from another subrepresentation to this subrepresentation

The equivariant vector space contains the matrices X such that

X * repC.image(g) = self.image(g) * X

Example

>>> S3 = replab.S(3);
>>> rep = S3.naturalRep;
>>> Xrep = randn(3, 3);
>>> triv = rep.subRep([1;1;1]);
>>> std = rep.maschke(triv);
>>> E = triv.subEquivariantFrom(std);
>>> Xsub = E.projectFromParent(Xrep);
>>> isequal(size(Xsub), [1 2]) % map to trivial (dim=1) from std (dim=2)
    1
>>> norm(Xsub) <= 1e-15
    1
Parameters
  • repC (replab.SubRep) – Subrepresentation on the source/column space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

  • parent – Equivariant space of self.parent and repC.parent, default: []

Kwtype special

charstring, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Kwtype parent

Equivariant, optional

Returns

The equivariant vector space

Return type

replab.SubEquivariant

subEquivariantTo(self, repR, varargin)

Returns the space of equivariant linear maps from this subrepresentation to another subrepresentation

The equivariant vector space contains the matrices X such that

X * self.image(g) = repR.image(g) * X

Example

>>> S3 = replab.S(3);
>>> rep = S3.naturalRep;
>>> Xrep = randn(3, 3);
>>> triv = rep.subRep([1;1;1]);
>>> std = rep.maschke(triv);
>>> E = triv.subEquivariantTo(std);
>>> Xsub = E.projectFromParent(Xrep);
>>> isequal(size(Xsub), [2 1]) % map to std (dim=1) from triv (dim=1)
    1
>>> norm(Xsub) <= 1e-15
    1
Parameters
  • repR (replab.SubRep) – Subrepresentation on the target/row space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

  • parent – Equivariant space of repR.parent and self.parent, default: []

Kwtype special

charstring, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Kwtype parent

Equivariant, optional

Returns

The equivariant vector space

Return type

replab.SubEquivariant

withNoise(self, injectionMapNoise, projectionMapNoise)

Adds Gaussian noise to the injection/projection maps

This method has two call conventions.

If sub.withNoise(sigma) is called, then the injection map is changed to injection1 = injection + Delta, where Delta is a matrix with normally distributed entries of standard deviation sigma. Delta is real (resp. complex) if the representation is real (resp. complex). Then the original projection map is ignored and projection1 is recovered as in replab.Rep.subRep.

If sub.withNoise(sigmaI, sigmaP) is called, then noise of magnitude sigmaI is added to the injection map and noise of magnitude sigmaP is added to the projection map. Then the new noisy maps are corrected so that injection1 * projection1 is still a projector.

Parameters
  • injectionMapNoise (double) – Magnitude of the noise to add to the injection map

  • projectionMapNoise (double, optional) – Magnitude of the noise to add to projection map

Returns

A noisy subrepresentation

Return type

SubRep

withUpdatedMaps(self, injection, projection, varargin)

Returns a copy of this subrepresentation with updated injection and projection maps

The following properties are copied from the original subrepresentation:

Additional keywords arguments are passed to the subRep method of parent.

Parameters
  • injection (double(D,d) or cyclotomic (D,d), may be sparse) – Injection map

  • projection (double(d,D) or cyclotomic (d,D), may be sparse) – Projection map

Returns

The updated subrepresentation

Return type

SubRep

Isotypic

class replab.Isotypic(parent, irreps, modelIrrep, projection)

Bases: replab.SubRep

Describes an isotypic component in the decomposition of a representation

An isotypic component regroups equivalent irreducible representations expressed in the same basis. Note that if the multiplicity is not one, there is a degeneracy in the picking of a basis of of the multiplicity space, and the basis is not necessarily chosen in a deterministic way by RepLAB.

However the subspace spanned by an isotypic component as a whole is unique.

The isotypic component is expressed as a subrepresentation of the representation being decomposed. However key methods are implemented more efficiently as more structure is available. In particular the computation of images is done in a way that minimizes numerical error and returns truly block diagonal matrices.

To cater for empty isotypic components, the isotypic component stores a “model” of the irreducible representation.

The isotypic stores its own projection map, as the following equation is not necessarily satisfied for individual irreps:

`` irreps{i}.projection * irrep{j}.injection = (i == j) * eye(d) ``

Own members

irreps

Equivalent irreducible subrepresentations in this isotypic component

Type

cell(1,*) of SubRep

modelIrrep

Irreducible representation identical to the subrepresentations present in this isotypic component

Type

Rep

Isotypic(parent, irreps, modelIrrep, projection)

Constructs an isotypic component of a parent representation

The injection map of the isotypic component comes from the sum of the injection maps of the irreps, while its projection map is supplied as an argument.

The static method fromIrreps computes this projection map if necessary.

Additional keyword arguments are passed to the Rep constructor.

Parameters
  • parent (replab.Rep) – Parent representation of which we construct a subrepresentation

  • irreps (cell(1,*) of replab.SubRep) – Irreducible subrepresentations

  • modelIrrep (replab.Rep) – Canonical model of the irreducible representation contained in this component

  • projection (double(*,*), may be sparse) – Embedding map of the isotypic component

static fromIrreps(parent, irreps, modelIrrep, varargin)

Creates an isotypic component from linearly independent equivalent irreducible subrepresentations

Parameters
  • parent (Rep) – Parent representation

  • irreps (cell(1,*) of SubRep) – Irreducible equivalent subrepresentations of parent

  • modelIrrep (Rep or []) – Model of the irreducible representation in this component (can be empty, and is then set to irreps{1})

  • irrepsAreBiorthogonal – True if the irreps injection/projection are biorthogonal, default: false

  • irrepsAreHarmonized – True if the irreps are in the same basis as modelIrrep, default: false

Kwtype irrepsAreBiorthogonal

logical, optional

Kwtype irrepsAreHarmonized

logical, optional

Returns

Isotypic component

Return type

Isotypic

static fromTrivialSubRep(trivial)

Builds an isotypic component from the trivial subrepresentation

The trivial subrepresentation must be complete, i.e. represent the space of all vectors that are fixed by the action of the representation.

Parameters

trivial (SubRep) – Maximal trivial subrepresentation of parent

Returns

The corresponding trivial isotypic component

Return type

Isotypic

irrep(self, i)

Returns the i-th copy of the irreducible representation

Returns

Irreducible subrepresentation of parent

Return type

SubRep

irrepDimension(self)

Returns the dimension of a single irreducible representation contained in this component

Returns

Irrep dimension

Return type

integer

irrepRange(self, i)

Returns the coefficient range where the i-th irrep lies

Parameters

i (integer) – Irreducible representation index

Returns

Range of the rows/columns block of the i-th irrep

Return type

integer(1,*)

isotypicEquivariantFrom(self, repC, varargin)

Returns the space of equivariant linear maps from another isotypic component to this isotypic component

The equivariant vector space contains the matrices X such that

X * repC.image(g) = self.image(g) * X

Both isotypic components must be harmonized.

Parameters
  • repC (replab.Isotypic) – Isotypic component, representation on the source/column space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

Kwtype special

charstring, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Returns

The equivariant vector space, or [] if the space has dimension zero or contains only the zero matrix

Return type

replab.IsotypicEquivariant or []

isotypicEquivariantTo(self, repR, varargin)

Returns the space of equivariant linear maps from this isotypic component to another isotypic component

The equivariant vector space contains the matrices X such that

X * self.image(g) = repR.image(g) * X

Both isotypic components must be harmonized.

Parameters
  • repR (replab.Isotypic) – Isotypic component, representation on the target/row space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

Kwtype special

charstring, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Returns

The equivariant vector space, or [] if the space has dimension zero or contains only the zero matrix

Return type

replab.IsotypicEquivariant or []

multiplicity(self)

Number of equivalent irreducible representations in this isotypic component

Returns

Multiplicity

Return type

integer

nIrreps(self)

Returns the number of irreps in this isotypic component, which is their multiplicity

Returns

Number of irreducible representations

Return type

integer

refine(self, varargin)

Refines this isotypic component

Assumes that the isotypic component is already close to an exact isotypic component, and refines its subspace by an iterative procedure applied on its injection and projection maps.

Parameters
  • numNonImproving – Number of non-improving steps before stopping the large-scale algorithm, default 20

  • nSamples – Number of samples to use in the large-scale version of the algorithm, default 5

  • maxIterations – Maximum number of (outer) iterations, default 1000

Kwtype numNonImproving

integer, optional

Kwtype nSamples

integer, optional

Kwtype maxIterations

integer, optional

Returns

Isotypic component with refined subspace (injection/projection maps)

Return type

Isotypic

Irreducible

class replab.Irreducible(parent, components)

Bases: replab.SubRep

Describes the irreducible decomposition of a representation

The mathematical background is based on Section 2.6 of Jean-Pierre Serre, “Linear representations of finite groups”.

An irreducible decomposition contains a sequence of isotypic components. When the irreducible decomposition is obtained through replab.Rep.decomposition with an exact argument, the order of isotypic components matches the real or complex character table of the represented group, and thus the irreducible decomposition can contain empty isotypic components. Those empty components can be stripped using the squeeze method.

As the decomposition of a representation is an expensive operation in RepLAB, Irreducible provides a pair of methods export and import that return a plain struct that can be saved into a .mat file under both MATLAB and Octave.

Own members

components

Isotypic components

Type

cell(1,*) of replab.Isotypic

component(self, i)

Returns a particular isotypic component in the decomposition

Parameters

i (logical) – Index of the isotypic component

Returns

The i-th isotypic component

Return type

replab.Isotypic

export(self)

Returns a plain data structure with the information about this decomposition

The returned data can be saved in a .mat file in both MATLAB and Octave and passed to the static import static method to reconstruct this irreducible decomposition.

Returns

Plain data structure

Return type

struct

static import(parent, data)

Reconstruct an irreducible decomposition from exported data

Parameters
  • parent (Rep) – Representation to decompose

  • data (struct) – Plain data resulting from export

irreducibleEquivariantFrom(self, repC, varargin)

Returns the space of equivariant linear maps from another irreducible decomposition to this irreducible decomposition

The equivariant vector space contains the matrices X such that

X * repC.image(g) = self.image(g) * X

Both irreducible decompositions must be harmonized.

Parameters
  • repC (replab.Irreducible) – Irreducible decomposition, representation on the source/column space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

Kwtype special

charstring, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Returns

The equivariant vector space, or [] if the space has dimension zero or contains only the zero matrix

Return type

replab.IrreducibleEquivariant or []

irreducibleEquivariantTo(self, repR, varargin)

Returns the space of equivariant linear maps from this irreducible decomposition to another irreducible decomposition

The equivariant vector space contains the matrices X such that

X * self.image(g) = repR.image(g) * X

Both irreducible decompositions must be harmonized.

Parameters
  • repR (replab.Irreducible) – Irreducible decomposition, representation on the target/row space

  • special – Special structure if applicable, see Equivariant, default: ‘’

  • type – Whether to obtain an exact equivariant space, default ‘double’ (‘double’ and ‘double/sparse’ are equivalent)

Kwtype special

‘commutant’, ‘hermitian’, ‘trivialRows’, ‘trivialCols’ or ‘’, optional

Kwtype type

‘exact’, ‘double’ or ‘double/sparse’, optional

Returns

The equivariant vector space, or [] if the space has dimension zero or contains only the zero matrix

Return type

replab.IrreducibleEquivariant or []

irrep(self, i, j)

Returns a subrepresentation in the irreducible decomposition

Parameters
  • i (integer) – Index of the isotypic component

  • j (integer, optional) – Index of the copy in the i-th isotypic component Default value is 1.

Returns

An irreducible subrepresentation

Return type

replab.SubRep

nComponents(self)

Returns the number of isotypic components in the decomposition

Returns

Number of isotypoic components

Return type

integer

squeeze(self)

Returns an irreducible decomposition where empty isotypic components have been removed

Useful to clean up the result of the exact decomposition, as the exact decomposition contains isotypic components for all irreducible representations, even if they are not present in the decomposed representation.

Example

>>> S3 = replab.S(3);
>>> rep = S3.naturalRep;
>>> dec = rep.decomposition('exact'); % doctest: +cyclotomic
>>> dec.nComponents
    3
>>> dec.squeeze.nComponents
    2
Returns

  • I (Irreducible) – Irreducible decomposition with empty isotypic components removed

  • ind (integer(1,*)) – Indices of non-empty components

toSubRep(self)

Returns the block-diagonal similar representation corresponding to this decomposition

Both this Irreducible object and the returned representation will be instances of SubRep. This method helps us test the implementation of Irreducible as the returned SubRep loses its particular structure.

The returned representation is the parent representation with an explicit change of basis, so it does not look as clean as this Irreducible.

Returns

The block-diagonal representation as a representation similar to parent

Return type

replab.SubRep