function [Y,S]=URFrameworkTest(P,E) % Test if a bar-framework with configuration/position P and graph % specified by edge set E is universal rigid (UR); % coded on 2/6/2011. Updated 2/10/2011 to output the configuration % matrix % Input: % P: d x n matrix where jth column is the position of jth point in R^d % E: 2 x |E| matrix where each column specifies two point-idexes of % an edge % Output: % Y: a max-rank configuration matrix satisfying the edge distance % constraint approximately % % rank(Y)=d+1 if and only if the framework is UR, if the graph contains % a (d+1) clique and the (d+1) points are in general position; see % So and Ye 2005. % % S: a max-rank PSD (dual) stress matrix approximately orthogonal to Y % % rank(S)=n-d-1 implies that the framework is UR, if all points are in % general position; see Alfakih and Ye 2010, and Connelly 1999 and % Gortler/Thurston 2009 for points in generic position. % % Check the dimension of the configuration and the number of edges % [d,n]=size(P); [m,mtotal]=size(E); toler = 1.e-8; % % Compute edge distances specified by E and set up the SDP inputs % b=[]; A=[]; for k =1:mtotal i = E(1,k); j = E(2,k); a=sparse(zeros(n,1)); a(i)=1; a(j)=-1; A=[A vec(a*a')]; dd= P(:,i)-P(:,j); rr= dd'*dd; b=[b; rr]; end; C=vec(sparse(n,n)); % % Set up algorithm parameter to call the SDP solver sedumi % K.s=[n]; pars.eps=toler; pars.stepdif=0; [Y,y,info] = sedumi(A,b,C,K,pars); Y = mat(Y); S = zeros(n,n); for k =1:mtotal S=S-y(k)*mat(A(:,k)); end; % % Prepare the outputs % YR = chol(Y); rankY = rank(Y,sqrt(toler)) r = rankY; Y =YR(1:r,:); r = - min([min(eig(S)),0]); S = (r+toler^2)*eye(n,n)+S; SR = chol(S); rankS = rank(S,sqrt(toler)) r = rankS; SR=SR(1:r,:); S = SR'*SR; if (r>0) & (r==n-d-1), YRP = SR\(-SR*Y'); Y=Y+YRP'; end % end of the function