close all;

ns = [ 10 20 ];

for n=ns,
  x = linspace( 0, 1, 1000 );
  y = (n^2*x).^n .* exp( -n^2 * x );
  
  figure; plot( x, y ); title( 'function plot' );
  
  disp( 'Integral of our function is: ' );
  trapz( x, y )
end