0001
0002
0003
0004
0005
0006
0007
0008 function [L] = DSA_BD_normOfOperator(blk,At,tol,m)
0009
0010
0011 if nargin < 4
0012 tol = 1e-6;
0013 end
0014 if nargin < 5
0015 m = size(At{1},2);
0016 end
0017
0018 fprintf('\n Calculating Norm of A');
0019 MAXITER = 2000;
0020 TOLERANCE = min(1e-6,tol);
0021 try
0022 AAt = sparse(m,m);
0023 for p = 1:size(blk,1)
0024 AAt = AAt + At{p,1}'*At{p,1};
0025 end
0026 pertdiag = 1e-13*ones(m,1);
0027 AAt = AAt + spdiags(pertdiag,0,m,m);
0028 options.disp = 0;
0029 options.maxit = MAXITER;
0030 options.tol = TOLERANCE;
0031 L = sqrt(svds(AAt,1,'L',options));
0032
0033 fprintf('\nNorm of A = %d\n',L);
0034 return;
0035 catch e
0036 fprintf('\nsvds(A,1,L,options) is not used.' );
0037 end
0038
0039 q = rand(size(At{1},2),1);
0040 q = q/norm(q);
0041 Atq = Atyfun1a(blk,At,q);
0042 L1 = 0;
0043 L2 = 1;
0044 q1 = q;
0045
0046 for i = 1:MAXITER
0047 z =DSA_BD_AXfun(blk,At,Atq);
0048 q = z/norm(z);
0049 Atq = DSA_BD_Atyfun(blk,At,q);
0050 L1 = q'*DSA_BD_AXfun(blk,At,Atq);
0051 error1 = abs(L1 - L2);
0052 error2 = norm(q - q1);
0053 if i > 1 && error1 < TOLERANCE && error2 < TOLERANCE
0054 break;
0055 end
0056 q1 = q;
0057 L2 = L1;
0058 end
0059 L =sqrt(L1);
0060 if (i == MAXITER)
0061 fprintf('WARNING: Power iteration did not converged. Errors: %f %f',error1,error2);
0062 end
0063 fprintf('\nPower iterations = %d',i);
0064 fprintf('\nNorm of A = %d\n',L);
0065 end