Home > K25104 > LinearSystems > deflation_alt_2.m

deflation_alt_2

PURPOSE ^

Given an nxn matrix A which has a 2-dimensional eigenspace

SYNOPSIS ^

function [ D ] = deflation_alt_2( A,v1,v2 )

DESCRIPTION ^

 Given an nxn matrix A which has a 2-dimensional eigenspace 
   spanned by v1 and v2, the function performs deflation creating 
   an nxn matrix which is similar to A. That is it has the same 
   eigenvalues as A. Uses an alternative algorithm.
 Input arguments:
   A, square matrix
   v1,v2, two vectors which span a known eigenspace of A
 Output arguments:
   D, deflated matrix

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [ D ] = deflation_alt_2( A,v1,v2 )
0002 % Given an nxn matrix A which has a 2-dimensional eigenspace
0003 %   spanned by v1 and v2, the function performs deflation creating
0004 %   an nxn matrix which is similar to A. That is it has the same
0005 %   eigenvalues as A. Uses an alternative algorithm.
0006 % Input arguments:
0007 %   A, square matrix
0008 %   v1,v2, two vectors which span a known eigenspace of A
0009 % Output arguments:
0010 %   D, deflated matrix
0011 
0012 [n,m]=size(A);      % finding the size of A
0013 [p,q]=size(v1);     % finding the size of v1
0014 [r,s]=size(v2);     % finding the size of v2
0015 if n~= m;
0016     error('input is not a square matrix');
0017 elseif p~=n || q~=1 || r~=n || s~=1;
0018     error('input dimensions do not agree');
0019 end
0020 
0021 % find index of first nonzero entry in v1
0022 k1=1;
0023 while v1(k1)==0
0024     k1=k1+1;
0025 end
0026 if k1>n
0027     error('zero vector is invalid input');
0028 end
0029 
0030 % find index of first nonzero entry in v2
0031 k2=1;
0032 while v2(k2)==0
0033     k2=k2+1;
0034 end
0035 if k2>n
0036     error('zero vector is invalid input');
0037 end
0038 
0039 if k1>1
0040     if k2>1         % first entry in v1 and v2 is zero
0041         P=eye(n);   % generate permutation matrix such that
0042         P(1,1)=0;   % first entry in P*v1 is non-zero
0043         P(1,k1)=1;
0044         P(k1,1)=1;
0045         P(k1,k1)=0;
0046         v1 = P*v1;
0047         v2 = P*v2;
0048         A = P*A*P'; % adjust A
0049     else            % first entry of v2 is non-zero, swap v1 and v2
0050         temp = v1;
0051         v1 = v2;
0052         v2 = temp;
0053     end
0054 end
0055 
0056 S1=eye(n);  % initializing S1, to the nxn identity
0057 for i=2:n;                 
0058     S1(i,1)=-v1(i)/v1(1);   % adjusting S1 so that the bottom n-1 entries
0059                             % of S1*v1 are zero
0060 end
0061 S1_inv=2*eye(n)-S1; % calculating the inverse of S1
0062 D = S1*A*S1_inv;
0063 
0064 v_hat=S1*v2;        % initializing v_hat to S1*v2
0065 % find index >= 2 of first nonzero entry in v_hat
0066 k=2;
0067 while v_hat(k)==0
0068     k=k+1;
0069 end
0070 if k>n
0071     error('zero vector is invalid input');
0072 end
0073 if k>2
0074     P=eye(n);   % generate permutation matrix such that
0075     P(2,2)=0;   % first entry in P*v is non-zero
0076     P(2,k)=1;
0077     P(k,2)=1;
0078     P(k,k)=0;
0079     v_hat = P*v_hat;
0080     D = P*D*P'; % adjust D
0081 end
0082 S_hat=eye(n-1); % initializing S_hat to the (n-1)x(n-1) identity
0083 v_hat=v_hat(2:n);   % letting v_hat be the lower n-1 entries of S1*v2
0084 for i=2:n-1;                                                         
0085     S_hat(i,1)=-v_hat(i)/v_hat(1);  % adjusting S_hat such that the bottom
0086                                     % n-2 entries of S_hat*v_hat are zero
0087 end
0088 S_hat_inv=2*eye(n-1)-S_hat; % calculating inverse of S_hat
0089 
0090 S2=eye(n);S2_inv=eye(n);    % initializing S2, S2_inv to the nxn identity
0091 
0092 S2(2:n,2:n)=S_hat;  % setting the bottom-right (n-1)x(n-1) submatrix of S2
0093                     % to S_hat, so that S2*(S1*v2) is a linear combination
0094                     % of the first two standard basis vectors.
0095                     
0096 S2_inv(2:n,2:n)=S_hat_inv;  % calculating inverse of S2
0097 
0098 D=S2*D*S2_inv;    % generating deflated matrix
0099 
0100 end
0101

Generated on Mon 18-Jan-2016 10:25:49 by m2html © 2005