function [zhat L ztilde Utilde] = sens_sel_approxnt(A, k) % Solves the problem % maximize log det (sum_{i=1}^m z_i a_i a_i^T) + kappa sum_{i=1}^m(log(z_i)+ log(1-z_i)) % subject to sum(z) = k % 0 <= z_i <= 1, i=1,..., m % variable z in R^m % problem parameters kappa (>0), a_1, ..., a_m in R^n % % see paper Sensor Selection via Convex Optimization % www.stanford.edu/~boyd/papers/sensor_selection.html % % Nov 2007 Siddharth Joshi & Stephen Boyd % Newton's method parameters MAXITER = 30; NT_TOL = 1e-3; GAP = 1.005; % Backtracking line search parameters alpha = 0.01; beta = 0.5; [m n] = size(A); z = ones(m,1)*(k/m); % initialize g = zeros(m,1); ones_m = ones(m,1); kappa = log(GAP)*n/m; % guarantees GM of lengths of semi-axes of ellipsoid corresponding to % ztilde <= 1.01 from optimal fprintf('\nIter. Step_size Newton_decr. Objective log_det\n'); fz = -log(det(A'*diag(z)*A)) - kappa*sum(log(z) + log(1-z)); fprintf(' 0\t -- \t -- %10.3f %10.3f\n', -fz, log(det(A'*diag(z)*A))); for iter=1:MAXITER W = inv(A'*diag(z)*A); V = A*W*A'; g = -diag(V)- kappa*(1./z - 1./(1-z)); H = V.*V + kappa*diag(1./(z.^2) + 1./((1-z).^2)); R = chol(H); Hinvg = (R\(R'\g)); Hinv1 = (R\(R'\ones_m)); dz = -Hinvg + ((ones_m'*Hinvg) / (ones_m'*Hinv1))*Hinv1; deczi = find(dz < 0); inczi = find(dz > 0); s = min([1; 0.99*[-z(deczi)./dz(deczi) ; (1-z(inczi))./dz(inczi)]]); while (1) zp = z + s*dz; fzp = -log(det(A'*diag(zp)*A)) - kappa*sum(log(zp) + log(1-zp)); if (fzp <= fz + alpha*s*g'*dz) break; end s = beta*s; end z = zp; fz = fzp; fprintf('%4d %10.3f %10.3f %10.3f %10.3f\n', iter, s, -g'*dz/2, -fz, log(det(A'*diag(z)*A))); if(-g'*dz/2 <= NT_TOL) break; end end zsort=sort(z); thres=zsort(m-k); zhat=(z>thres); L = log(det(A'*diag(zhat)*A)); ztilde = z; Utilde = log(det(A'*diag(z)*A)) + 2*m*kappa;