Optimise creation of matern 3/2 covariance matrix

I'm making a function to create a Matern 3/2 covariance matrix (see below). The size of the matrix varies but is typically with n >= 1e5. My attempt currently takes very long.

How can I increase the performance of this code?

function K = matern32(n,tau)
% Function creates n x n matern covariance matrix

K = sparse(n,n);
eps = 1e-9;
for i = 1:n
    for j = 1:n
        K(i,j) = (1+sqrt(3)*abs(i-j)/tau)*exp(-sqrt(3)*abs(i-j)/tau);
        if K(i,j) < eps
            K(i,j) = 0;
        end 
    end 
end

end

Answer 1

[darova]

Try this

[X,Y] = meshgrid(1:n);
K = (1+sqrt(3)*abs(X-Y)/tau).*exp(-sqrt(3)*abs(X-Y)/tau);
K = (K>1e-9).*K;

read this