function gibbs(N,x0)
%This script allows us to study the Gibbs phenomenon
%We can either let N be a list of positive integers, or a single
%integer.
%If N is a single integer, then we compute the sum_k=1^N sin(kx)/k for x in
%the interval [-x0,x0].
%
%If N is an arithmetic sequence, then we compute sum_k=1^N(j) sin(kx)/k for x in
%the interval [-x0/j,x0/j].
%
if nargin == 1
    x0=pi;
end
lN=length(N);
figure
for j=1:lN
x=linspace(-x0/j,x0/j,2000);
XX=x'*(1:N(j));
XX=sin(XX);
coef=(1./(1:N(j)))';
SN=XX*coef;
MX=max(SN)/(pi/2);
mn=min(SN)/(pi/2);
subplot(1,lN,j)
plot(x,SN)
xlim([-x0/j,x0/j])
title(sprintf('N = %g, max = %g, min = %g',N(j), MX,mn))
hold on
end
