n = 10;
pn = 0;
format long
disp(['       p(n)                   pi'])
disp('   --------------------------------------')
for k = 0:n
   pn = pn+ (1/16^k)*(4/(8*k+1)-2/(8*k+4)-1/(8*k+5)-1/(8*k+6)) ;   
   disp([pn pi])
end
disp('Lopullinen:')
disp(['p(n)      pi'])
disp([pn pi])
format short