ConvergenceDemo(SourceSwitch);
ConvergenceDemo(SourceSwitch) compares the polynomial EVD approximations by
the functions SBR2() and SMD() for different types of parahermitian matrices.
The demo provides the reduction in off-diagonal energy in increments of 5
iteration steps.
Theinput parameter SourceSwitch can be selected as:
ComparisonDemo('Random') uses an arbitrary complex valued 5x5 polynomial
parahermitian matrix of order 10.
ComparisonDemo('Simple1') compares decompositions of a 3x3 real valued
polynomial parahermitian matrix R(z) of order 4:
| 1 -.4*z 0 |
R(z) = | -.4*z 1 .5*z^{-2} |
| 0 .5*z^2 3 |
ComparisonDemo('Simple2') compares the decompositions of a 5x5 complex valued
polynomial parahermitian matrix R(z) of order 4:
| 1 .5*z^2 -.4j*z^{-1} 0 .2j*z^2 |
| .5*z^{-2} 1 0 0 0 |
R(z) = | .4j*z 0 3 .3z^{-1} 0 |
| 0 0 .3*z .5 0 |
|-.2j*z^{-2} 0 0 0 .25 |
This matrix admits a PEVD with a non-polynomial Gamma(z).
The function generates a plot showing the remaining off-diagonal power,
normalised w.r.t. the total power of the matrix; this indicates how well the
two algorithms perform their diagonalisation task. A second figure compares
the off-diagonal power in relation to the order of the paraunitary matrices
required in order to accomplish this decomposition.
Input parameter:
SourceSwitch selects parahermitian matrx to be decomposed
default: 'Random'
Output parameter: none
0001 function ConvergenceDemo(SourceSwitch); 0002 %ConvergenceDemo(SourceSwitch); 0003 % 0004 % ConvergenceDemo(SourceSwitch) compares the polynomial EVD approximations by 0005 % the functions SBR2() and SMD() for different types of parahermitian matrices. 0006 % The demo provides the reduction in off-diagonal energy in increments of 5 0007 % iteration steps. 0008 % 0009 % Theinput parameter SourceSwitch can be selected as: 0010 % 0011 % ComparisonDemo('Random') uses an arbitrary complex valued 5x5 polynomial 0012 % parahermitian matrix of order 10. 0013 % 0014 % ComparisonDemo('Simple1') compares decompositions of a 3x3 real valued 0015 % polynomial parahermitian matrix R(z) of order 4: 0016 % | 1 -.4*z 0 | 0017 % R(z) = | -.4*z 1 .5*z^{-2} | 0018 % | 0 .5*z^2 3 | 0019 % 0020 % ComparisonDemo('Simple2') compares the decompositions of a 5x5 complex valued 0021 % polynomial parahermitian matrix R(z) of order 4: 0022 % | 1 .5*z^2 -.4j*z^{-1} 0 .2j*z^2 | 0023 % | .5*z^{-2} 1 0 0 0 | 0024 % R(z) = | .4j*z 0 3 .3z^{-1} 0 | 0025 % | 0 0 .3*z .5 0 | 0026 % |-.2j*z^{-2} 0 0 0 .25 | 0027 % This matrix admits a PEVD with a non-polynomial Gamma(z). 0028 % 0029 % The function generates a plot showing the remaining off-diagonal power, 0030 % normalised w.r.t. the total power of the matrix; this indicates how well the 0031 % two algorithms perform their diagonalisation task. A second figure compares 0032 % the off-diagonal power in relation to the order of the paraunitary matrices 0033 % required in order to accomplish this decomposition. 0034 % 0035 % Input parameter: 0036 % SourceSwitch selects parahermitian matrx to be decomposed 0037 % default: 'Random' 0038 % Output parameter: none 0039 0040 % S. Weiss, University of Southampton, 8/10/2014 0041 0042 if nargin==0, 0043 SourceSwitch='Random'; 0044 end; 0045 0046 disp('SBR2 and SMD Convergence Comparison Demo'); 0047 disp('--------------------------------------------'); 0048 0049 %--------------------------------------- 0050 % Define Scenario 0051 %--------------------------------------- 0052 % create parahermitian matrix to be decomposed 0053 if strcmp(SourceSwitch,'Random')==1, 0054 A = randn(5,5,6) + sqrt(-1)*randn(5,5,6);; 0055 R = PolyMatConv(A,ParaHerm(A)); 0056 elseif strcmp(SourceSwitch,'Simple1')==1, 0057 R = zeros(3,3,5); 0058 R(:,:,3) = diag([1 1 3]); 0059 R(1,2,1) = .5; 0060 R(2,1,5) = .5; 0061 R(1,2,2) = -.4; 0062 R(2,1,4) = -.4; 0063 elseif strcmp(SourceSwitch,'Simple2')==1, 0064 R = zeros(5,5,5); 0065 R(:,:,3) = diag([1 1 3 .5 .25]); 0066 R(1,2,1) = .5; 0067 R(2,1,5) = .5; 0068 R(1,3,4) = -.4*sqrt(-1); 0069 R(3,1,2) = .4*sqrt(-1); 0070 R(1,5,1) = .2*sqrt(-1); 0071 R(5,1,5) = -.2*sqrt(-1); 0072 R(4,3,2) = .3; 0073 R(3,4,4) = .3; 0074 else 0075 error('option for input parameter SourceSwitch not implemented'); 0076 end; 0077 0078 %--------------------------------------- 0079 % Set incremental parameters 0080 %--------------------------------------- 0081 % parameters 0082 maxiter = 5; % max number of iterations per update 0083 blocks = 20; 0084 epsilon = 10^(-10); % alternative stopping criterion 0085 % ensure that this is not invoked 0086 mu = 0; % truncate true zeroes 0087 GammaSMD = R; % initialisations 0088 GammaSBR2 = R; 0089 N1 = PolyMatNorm(R); 0090 N2 = PolyMatNorm(R,'OffDiag'); 0091 OffDiagNormSMD(1) = N2/N1; % normalised remaining off-diag. energy 0092 OffDiagNormSBR2(1) = N2/N1; 0093 OrderSMD(1)=0; 0094 OrderSBR2(1)=0; 0095 0096 %--------------------------------------- 0097 % Recursively decompose in blocks of 5 iterations 0098 %--------------------------------------- 0099 for n = 2:blocks, 0100 % SMD iteration 0101 [H,GammaSMD] = SMD(GammaSMD,maxiter,epsilon,mu,'SMD'); 0102 OffDiagNormSMD(n) = PolyMatNorm(GammaSMD,'OffDiag')/N1; 0103 OrderSMD(n) = OrderSMD(n-1)+size(H,3)-1; 0104 % SBR2 iteration 0105 [H,GammaSBR2] = SBR2(GammaSBR2,maxiter,epsilon,mu,'SBR2'); 0106 OffDiagNormSBR2(n) = PolyMatNorm(GammaSBR2,'OffDiag')/N1; 0107 OrderSBR2(n) = OrderSBR2(n-1)+size(H,3)-1; 0108 end; 0109 0110 %--------------------------------------- 0111 % Display resulting matrices 0112 %--------------------------------------- 0113 % display normalised remaining off-diagonal power vs iterations 0114 figure(1); clf; 0115 plot(1:maxiter:maxiter*blocks,5*log10(OffDiagNormSBR2),'bo'); 0116 hold on; 0117 plot(1:maxiter:maxiter*blocks,5*log10(OffDiagNormSMD),'r*'); 0118 legend('SBR2','SMD'); 0119 xlabel('iterations'); 0120 ylabel('norm. remaining off-diagonal power / [dB]'); 0121 disp('Figure 1: norm. remaining off-diagonal power vs iterations'); 0122 0123 % display normalised remaining off-diagonal power vs paraunitary order 0124 figure(2); clf; 0125 plot(OrderSBR2,5*log10(OffDiagNormSBR2),'bo'); 0126 hold on; 0127 plot(OrderSMD,5*log10(OffDiagNormSMD),'r*'); 0128 legend('SBR2','SMD'); 0129 xlabel('paraunitary order (without truncation)'); 0130 ylabel('norm. remaining off-diagonal power / [dB]'); 0131 disp('Figure 2: norm. remaining off-diagonal power vs paraunitary order');