MAXCARDSEARCH Determine graph chordality and maximal cliques. [MAXCLIQUE,ISCHORDAL] = MAXCARDSEARCH(AADJ) Determine the maximal cliques for a chordal graph described by the adjacency matrix Aadj. Also test the graph for chordality. The algorithms for determining the maximal cliques and testing for chordality are described in [1]. Inputs: AADJ : The adjacency matrix of a graph. If the graph is chordal, the maximal cliques for this graph are calculated. Outputs: MAXCLIQUE : If the graph described by Aadj is chordal, maxclique returns a matrix describing the maximal cliques of the graph. Each column of maxclique represents a clique with nonzero entries indicating the nodes included in the maximal clique. If Aadj is not chordal, maxclique is set to NaN. ISCHORDAL : Returns 1 if the graph described by Aadj is chordal, otherwise returns 0. [1] Tarjan, R. E., and M. Yannakakis. "Simple Linear-Time Algorithms to Test Chordality of Graphs, Test Acyclicity of Hypergraphs, and Selectively Reduce Acyclic Hypergraphs." SIAM Journal on computing 13 (1984): 566.
0001 function [maxclique,ischordal] = maxCardSearch(Aadj) 0002 %MAXCARDSEARCH Determine graph chordality and maximal cliques. 0003 % [MAXCLIQUE,ISCHORDAL] = MAXCARDSEARCH(AADJ) 0004 % 0005 % Determine the maximal cliques for a chordal graph described by the 0006 % adjacency matrix Aadj. Also test the graph for chordality. The 0007 % algorithms for determining the maximal cliques and testing for 0008 % chordality are described in [1]. 0009 % 0010 % Inputs: 0011 % AADJ : The adjacency matrix of a graph. If the graph is chordal, 0012 % the maximal cliques for this graph are calculated. 0013 % 0014 % Outputs: 0015 % MAXCLIQUE : If the graph described by Aadj is chordal, maxclique 0016 % returns a matrix describing the maximal cliques of the graph. 0017 % Each column of maxclique represents a clique with nonzero 0018 % entries indicating the nodes included in the maximal clique. If 0019 % Aadj is not chordal, maxclique is set to NaN. 0020 % ISCHORDAL : Returns 1 if the graph described by Aadj is chordal, 0021 % otherwise returns 0. 0022 % 0023 % [1] Tarjan, R. E., and M. Yannakakis. "Simple Linear-Time Algorithms to 0024 % Test Chordality of Graphs, Test Acyclicity of Hypergraphs, and 0025 % Selectively Reduce Acyclic Hypergraphs." SIAM Journal on computing 13 0026 % (1984): 566. 0027 0028 % MATPOWER 0029 % Copyright (c) 2013-2019, Power Systems Engineering Research Center (PSERC) 0030 % by Daniel Molzahn, PSERC U of Wisc, Madison 0031 % 0032 % This file is part of MATPOWER/mx-sdp_pf. 0033 % Covered by the 3-clause BSD License (see LICENSE file for details). 0034 % See https://github.com/MATPOWER/mx-sdp_pf/ for more info. 0035 0036 %% Setup 0037 0038 nbus = size(Aadj,1); 0039 nline = nnz(Aadj)/2; 0040 0041 %% Create a perfect elimination ordering 0042 0043 % Store in s{i} all unnumbered vertices adjacent to exactly i numbered 0044 % vertices. 0045 s{nline} = []; 0046 s{1} = 1:nbus; % This corresponds to s{0} in the paper's notation 0047 0048 % Maintain the largest index j such that s{j} is nonempty 0049 jidx = 1; 0050 0051 % Candidate perfect elimination ordering, with vertex numbering in vnum. 0052 % Index of vnum corresponds to the original vertex number, and the value of 0053 % vnum indicates the new vertex number. 0054 alpha = zeros(nbus,1); 0055 0056 % v2s stores which set contains a given vertex? 0057 v2s = ones(nbus,1); 0058 0059 i = nbus; 0060 0061 while i >= 1 0062 0063 % To carry out a step of the search, we remove a vertex v from s(jidx) 0064 % and number it. 0065 v = s{jidx}(1); 0066 if length(s{jidx}(:)) > 1 0067 s{jidx} = s{jidx}(2:end); 0068 else 0069 s{jidx} = []; 0070 end 0071 alpha(v) = i; 0072 v2s(v) = 0; % This vertex is no longer in a set 0073 0074 % For each unnumbered vertex w adjacent to v, we move w from the set 0075 % containing it to the next set 0076 vadj = find(Aadj(:,v)); 0077 for k=1:length(vadj) 0078 % If this node isn't already numbered, remove it from the original set 0079 if v2s(vadj(k)) ~= 0 0080 s{v2s(vadj(k))} = s{v2s(vadj(k))}(s{v2s(vadj(k))}(:) ~= vadj(k)); 0081 0082 % Add it to the new set 0083 s{v2s(vadj(k))+1} = [s{v2s(vadj(k))+1}(:); vadj(k)]; 0084 0085 % Update v2s 0086 v2s(vadj(k)) = v2s(vadj(k)) + 1; 0087 end 0088 end 0089 0090 % Add one to jidx and then decrease jidx until reaching a non-empty s(jidx) 0091 jidx = jidx + 1; 0092 if jidx > length(s) 0093 jidx = jidx - 1; 0094 end 0095 while jidx >= 1 && isempty(s{jidx}) 0096 jidx = jidx - 1; 0097 end 0098 0099 i = i-1; 0100 0101 end 0102 0103 0104 %% Check for chordality 0105 0106 f = zeros(nbus,1); 0107 index = zeros(nbus,1); 0108 ischordal = 1; 0109 for i=1:nbus 0110 w = find(alpha == i); 0111 f(w) = w; 0112 index(w) = i; 0113 0114 valid_v = find(Aadj(:,w)); 0115 valid_v = valid_v(alpha(valid_v) < i); 0116 for vidx = 1:length(valid_v) 0117 v = valid_v(vidx); 0118 index(v) = i; 0119 if f(v) == v 0120 f(v) = w; 0121 end 0122 end 0123 0124 for vidx = 1:length(valid_v) 0125 v = valid_v(vidx); 0126 if index(f(v)) < i 0127 ischordal = 0; 0128 break; 0129 end 0130 end 0131 0132 if ~ischordal 0133 break; 0134 end 0135 end 0136 0137 0138 %% Determine maximal cliques 0139 0140 % According to https://en.wikipedia.org/wiki/Chordal_graph 0141 % "To list all maximal cliques of a chordal graph, simply find a perfect 0142 % elimination ordering, form a clique for each vertex v together with the 0143 % neighbors of v that are later than v in the perfect elimination ordering, 0144 % and test whether each of the resulting cliques is maximal." 0145 0146 % alpha is a candidate perfect elimination ordering. First form all 0147 % cliques. Then determine which cliques are maximal. Put maximal cliques in 0148 % the columns of the variable maxclique. 0149 0150 if ischordal 0151 % Form a matrix representation of the cliques 0152 clique = speye(nbus,nbus); 0153 for i=1:nbus 0154 % neighbors of node i that are later in the ordering 0155 neighbors = find(Aadj(i,:)); 0156 neighbors = neighbors(alpha(neighbors) > alpha(i)); 0157 0158 clique(i,neighbors) = 1; 0159 end 0160 0161 % Test whether each clique is maximal. A clique is not maximal if a node 0162 % neighboring any of the nodes in the clique, but not in the clique, is 0163 % connected to all nodes in the clique. 0164 i=0; 0165 nclique = size(clique,1); 0166 while i < nclique 0167 i = i+1; 0168 0169 neighbors = []; 0170 cliquei = find(clique(i,:)); 0171 cliquei_bool = clique(i,:).'; 0172 nnzi = sum(cliquei_bool); 0173 0174 for k = 1:length(cliquei) % Check all nodes adjacent to nodes in the clique, but not included in the clique 0175 neighbors = [neighbors find(Aadj(cliquei(k),:))]; 0176 end 0177 neighbors = unique(setdiff(neighbors,cliquei)); 0178 0179 not_maximal = 0; 0180 for m=1:length(neighbors) 0181 % If this neighbor is connected to all other nodes in the clique, 0182 % this is not a maximal clique 0183 if Aadj(neighbors(m),:) * cliquei_bool >= nnzi 0184 not_maximal = 1; 0185 break; 0186 end 0187 end 0188 0189 if not_maximal 0190 clique = clique([1:i-1 i+1:nclique],:); % Delete this clique 0191 i = i-1; % After deletion, the next clique is really numbered i 0192 nclique = nclique - 1; 0193 end 0194 end 0195 maxclique = clique.'; % put maximal cliques in the columns 0196 else 0197 maxclique = nan; % maximal clique algorithm only works for chordal graphs 0198 end