function C_N=two_exp(sig, ell)
% C.E. Powell MATH64082
% Generate C_N matrix for Example 2.2 - assuming grid vertices are ordered
% left to right, bottom to top on the chosen grid on D=[0,1] x [0,1].
%
% inputs:
%        sig - standard deviation of random field
%        ell - correlation length of random field
% output:
%        C_N - N x N matrix whose (i,j)th entry is C(\vec{x}_i, \vec{x}_j)
%        where C(\vec{x},\vec{y}) is the two-dimensional exponential
%        covariance function.

n=100; h=1/n; x1=0:h:1; x2=x1;
[xx1,xx2]=meshgrid(x1,x2);
% xx1 contains first coordinates and xx2 contains second coordinates.
xx1=xx1'; xx1=xx1(:); xx2=xx2'; xx2=xx2(:);

bigN=(n+1)^2;
C_N=zeros(bigN,bigN);
% Double loop not efficient - but 'works'! 
for i=1:bigN
    for j=1:bigN
        C_N(i,j)=sig*sig*exp(-sqrt((((xx1(i)-xx1(j))^2 + (xx2(i)-xx2(j))^2)))/ell);
    end
end
