Lisa, you are correcting above for the window, which is what I said initially. You're not working with power in the above example, which is why you are only using sum(w).
The factor of two gives you the correction for a single-sided spectral estimate for all frequencies except 0 and the Nyquist.
You should not multiply 0 and the Nyquist by two.
The mistake you have is that you are plotting the whole spectrum [-Nyquist,Nyquist). You only want to apply the factor of two if you are plotting just a single side, so:
fs=1000;
T=1/fs;
N=512;
n=0:1:512-1;
t=(0:1:N-1)*T;
xn=2*sin(2*pi*100/fs*n);
t=(0:1:N-1)*T;
w=flattopwin(N);
xn=xn(:);
xn1=xn.*w;
ws=sum(w);
Xk = fft(xn1);
Xk = Xk(1:length(Xk)/2+1);
Xk = Xk.*(1/ws);
Xk(2:end-1) = Xk(2:end-1)*2;
plot(abs(Xk));
