  function [xMax,phiMax] = QuadForm(A)
% A is n-by-n and symmetric.
% xMax is a column vector that maximizes x'*A*x/x'*x  subject to the constraint
% that x(1) + ... + x(n) = 0.
% phiMax = xMax'*A*xMax/xMax'*xMax.

[n,n] = size(A);
e = ones(n,1);
[v,alfa] = House(e);
% Suppose Q = I-2vv' = [Q1 Q2] where Q1 is the first column of Q. 
% Note that Q2'*e = 0.
% phi(Q2*z) = (Q2*z)'*A*(Q2*z)/(Q2*z)'*(Q2*z) = z'*(Q2'*A*Q2)/(z'*z)
%
% Thus, the "best" z is the eigenvector associated with the most positive
% eigenvalue of the (n-1)-by-(n-1) matrix C = Q2'*A*Q2. Note that
% if Atilde = Q'*A*Q then C = Atilde(2:n,2:n).
%
% Atilde = (I - 2vv')A(I-2vv') = (I - 2vv')(A - 2Avv')
%        = A - 2Avv' - 2vv'A + 4vv'Avv'
%        = A - 2wv' - 2vw' + 4(v'w)vv'
 w = A*v;
 Atilde = A - (2*w)*v' - (2*v)*w' + (4*(v'*w))*v*v';
 [U,D] = schur(Atilde(2:n,2:n));
 % Find the most positive eigenvalue and extract the eigenvector
 d = diag(D);
 [m,iMax] = max(diag(D));
 zMax = U(:,iMax);
 
 %xMax = Q(:,2:n)*zMax = Q*[0;zMax] = (I-2vv')*z = z - (2*v'z)v where
 
 z = [0;zMax];
 
 xMax = z - (2*v'*z)*v; 
 phiMax = m;
 