theta0 = -10:0.1:10; % parameter will vary from -10 to 10 by 0.1
theta1 = theta0; % the same

% create all possible combinations:
[theta0m,theta1m] = meshgrid(theta0, theta1);

x = 1:3; % training instances
y = x;  % this is just y = 1*x + 0
m = length(x);

h = zeros(size(theta0m)(1), size(theta0m)(2), m);
J = zeros(size(theta0m));

for i = 1:m
      h(:,:,i) = theta0m .+ theta1m * x(i);
          J = J + (h(:,:,i) - y(i)) .^ 2;
endfor

J = J/(2*m);

contour(theta0, theta1, J, 1:5);
xlabel('\theta_0');
ylabel('\theta_1');
title('Contour plot for cost function J');
print -dpng contour.png