function run_poisson(method)

% Solve the Poisson equation with Jacobi's method
% Without and with multigrid

N = 2048;
h = 1/N;
x = (0:h:1)';      % N+1 points
f = - 1 + 100*cos(24*pi*x) + 1000*cos(100*pi*x);

% Direct solver
e = ones(N+1,1);
A = spdiags([e -2*e e], -1:1, N-1, N-1)/h^2;
u = A \ f(2:(end-1));
u_exact = [0; u; 0];

if method == 1,
figure(1); clf; plot(x,u_exact,'b-');

elseif method == 2,
% Jacobi iterations
u = solve_poisson_jacobi(f, u_exact);

elseif method == 3,
% Multigrid Jacobi iterations
tic,
u = solve_poisson_multigrid(f);
toc
figure(1); clf;
plot(x,u_exact,'b-'); hold on;
plot(x,u,'rx');
end