FUNCTION SUMMARY: Subroutine PartialSymLU do partial symbolic LU factorization. [ERPU,CIndxU,ERPEQ,CIndxEQ] = PartialSymLU(CIndx,ERP,dim,Stop,BoundBus) INPUT DATA: CIndx - N*1 array containning the column index of the rows, N is the number of non-zeros elements in the input data ERP - (N_node+1)*1 array containning the end of row pointer data dim - scalar, dimension of the input matrix Stop - scalar, stop sign of the LU factorization (The LU factorization in the reduction process is not complete but partial) BoundBus- 1*n array, includes indices of boundary buses OUTPUT DATA: ERPU - N_dim*1 array containing end of row pointer of all rows except the last row which doesn have any off diagonal element, N_dim is the dimension of the input matrix A in the original Ax = b problem CIndxU - N*1 array containing the column index of off diagonal element in each row in the U matrix. The length N depends on the number of native plus filled non-zero elements in the off diagonal position in U matrix. CIndxU is unordered. ERPEQ, CIndxEQ, - both are N*1 arrays, together include the indices of equivalent branches Note: This subroutine will do: 1. Factorization of rows in the full model bus addmittance matrix corresponding to the external buses. This process will only generate the pointers of non-zeroes in L and U matrix; 2. Identify the equivalent branch indices. Identify the pointers (ERPEQ, CIndxEQ)of fills generated from step 1 in the rows and columns corresponding to the boundary buses. (The equivalent lines can only span the boundary buses).
0001 function [ERPU,CIndxU,ERPEQ,CIndxEQ] = PartialSymLU(CIndx,ERP,dim,Stop,BoundBus) 0002 %FUNCTION SUMMARY: 0003 % Subroutine PartialSymLU do partial symbolic LU factorization. 0004 % 0005 % [ERPU,CIndxU,ERPEQ,CIndxEQ] = PartialSymLU(CIndx,ERP,dim,Stop,BoundBus) 0006 % 0007 % INPUT DATA: 0008 % CIndx - N*1 array containning the column index of the rows, N is the 0009 % number of non-zeros elements in the input data 0010 % ERP - (N_node+1)*1 array containning the end of row pointer data 0011 % dim - scalar, dimension of the input matrix 0012 % Stop - scalar, stop sign of the LU factorization (The LU factorization 0013 % in the reduction process is not complete but partial) 0014 % BoundBus- 1*n array, includes indices of boundary buses 0015 % 0016 % 0017 % OUTPUT DATA: 0018 % ERPU - N_dim*1 array containing end of row pointer of all rows except 0019 % the last row which doesn have any off diagonal element, N_dim is 0020 % the dimension of the input matrix A in the original Ax = b 0021 % problem 0022 % CIndxU - N*1 array containing the column index of off diagonal element 0023 % in each row in the U matrix. The length N depends on the 0024 % number of native plus filled non-zero elements in the off 0025 % diagonal position in U matrix. CIndxU is unordered. 0026 % ERPEQ, CIndxEQ, - both are N*1 arrays, together include the indices of 0027 % equivalent branches 0028 % 0029 % 0030 % Note: This subroutine will do: 0031 % 1. Factorization of rows in the full model bus addmittance matrix 0032 % corresponding to the external buses. This process will only 0033 % generate the pointers of non-zeroes in L and U matrix; 0034 % 2. Identify the equivalent branch indices. Identify the pointers 0035 % (ERPEQ, CIndxEQ)of fills generated from step 1 in the rows and 0036 % columns corresponding to the boundary buses. (The equivalent lines 0037 % can only span the boundary buses). 0038 0039 % MATPOWER 0040 % Copyright (c) 2014-2016, Power Systems Engineering Research Center (PSERC) 0041 % by Yujia Zhu, PSERC ASU 0042 % 0043 % This file is part of MATPOWER/mx-reduction. 0044 % Covered by the 3-clause BSD License (see LICENSE file for details). 0045 % See https://github.com/MATPOWER/mx-reduction/ for more info. 0046 0047 %% Begin of the code 0048 numRow = dim; % number of rows of given data matrix 0049 % numCol = dim; % number of columns of given data matrix 0050 %% Initialization 0051 % the length of pos is equal to # of total non zero data - number of 0052 % diagonal elements - the number of non zero element in last row' 0053 %% step 0 0054 CIndxU = 0; 0055 ERPU= zeros(1,dim); % with aditional one digit for the 7th row 0056 Switch = zeros(1,dim); 0057 MinNod0 = inf; % This is a initiate large value of MinNod; not the MinNod used in building the symbolic struture 0058 MinNod = MinNod0; 0059 MinNod1=MinNod; 0060 Link = ERPU; 0061 CIndxEQ=0; 0062 ERPEQ=zeros(1,Stop+length(BoundBus)+1); 0063 Chain=zeros(1,Stop+length(BoundBus)); 0064 % preprocess the data by ordering the CIndx of every row in ascending order 0065 for i = 2:numRow+1 0066 RowColInd = ERP(i-1)+1:ERP(i); % for every row the pointer of the column index 0067 CIndx(RowColInd)=sort(CIndx(RowColInd)); 0068 end 0069 % initiate starting row 1 0070 RowIndex = 1; 0071 while RowIndex<=length(BoundBus)+Stop 0072 MinNod1=MinNod0; 0073 %% step 1 0074 if RowIndex<=Stop 0075 [CIndxU,ERPU,MinNod,Switch] = RODAssignment(ERP,CIndx,CIndxU,ERPU,MinNod,Switch,RowIndex,RowIndex); 0076 %% step 2 0077 % Check fill in ROD in current row 0078 SelfRef = RowIndex; % self referential link pointer to next availale link element 0079 while (Link(SelfRef)~=0) 0080 SelfRef = Link(SelfRef); % to refer for the next link element 0081 [CIndxU,ERPU,MinNod,Switch] = RODAssignment (ERPU,CIndxU,CIndxU,ERPU,MinNod,Switch,SelfRef,RowIndex); 0082 end 0083 Link(RowIndex) = 0; % zero the element in the Link list of current row 0084 %% step 3 0085 % update Link node (MinNod) to be current RowIndex 0086 if MinNod>RowIndex&&MinNod~=MinNod0 0087 SelfRef = MinNod; % start the assign value to Link 0088 while (Link(SelfRef)~=0) 0089 SelfRef = Link(SelfRef); 0090 end 0091 Link(SelfRef)=RowIndex; 0092 end 0093 %% step 4 0094 MinNod = MinNod0; % reset the MinNod value 0095 else 0096 if Link(RowIndex)~=0 0097 [LinkPos,LinkArray,Counter] = SelfLink(Link,RowIndex); 0098 for i = 1:Counter 0099 SelfRef=LinkArray(i); 0100 if SelfRef<=Stop 0101 if CIndxEQ==0 0102 StFlag=0; 0103 end 0104 [CIndxEQ,ERPEQ,MinNod,Switch,ChainFlag] = EQRODAssignment(ERPU,CIndxU,CIndxEQ,ERPEQ,MinNod,Switch,SelfRef,RowIndex,Chain,MinNod1); 0105 if StFlag==0&&any(CIndxEQ) 0106 StFlag=1; 0107 end 0108 if MinNod>RowIndex&&MinNod~=MinNod0||(MinNod1~=MinNod&&ChainFlag~=1) 0109 if MinNod1~=MinNod0 0110 if StFlag~=1&&MinNod1~=MinNod&&Chain(MinNod1)==0 0111 Link(SelfRef1)=0; 0112 % If the chain breaks, record where and which 0113 % row to be reconected 0114 Chain(MinNod1)=Link(MinNod1); 0115 end 0116 elseif StFlag~=1&&MinNod1~=MinNod&&ChainFlag==0 0117 Link(SelfRef1)=0; 0118 0119 end 0120 0121 if MinNod~=MinNod0 0122 SelfRef1=SelfRef; 0123 FillIn = MinNod; % start the assign value to Link 0124 while (Link(FillIn)~=0)&&Link(FillIn)~=FillIn; 0125 FillIn = Link(FillIn); 0126 end 0127 % if Link(FillIn)==0 0128 if FillIn~=SelfRef1 0129 Link(FillIn)=SelfRef1; 0130 else 0131 MinNod=MinNod; 0132 end 0133 Chain(MinNod)=0; 0134 % end 0135 end 0136 end 0137 %% step 4 0138 MinNod1=MinNod; 0139 0140 MinNod = MinNod0; % reset the MinNod value 0141 end 0142 end 0143 else 0144 ERPEQ(RowIndex+1)=ERPEQ(RowIndex); 0145 end 0146 Link(RowIndex) = 0; % zero the element in the Link list of current row 0147 end 0148 % ready for next loop 0149 RowIndex = RowIndex+1; 0150 end 0151 end