function A = myhilb(n,method)
% A = myhilb(n,method) computes the Hilbert matrix of order n using a
% method specified by the variable method, which can be an integer from 1
% to 6.

% Initialise the matrix as an n x n matrix with all entries equal to 0
% It would be also *possible* to skip this line and sipmly build up the
% matrix one row at a time. In this case, however, Matlab does not know
% from the start which size the whole matrix will end up and thus has to
% assign the matrix a new place in the memory during each iteration. This
% can greatly slow down the execution and should therefore be avoided
% whenever possible.
A = zeros(n,n);

% set the default method to be the first in case the variable method has
% not been specified by the user
if nargin < 2
    method = 1;
end

switch method
    case 1
        % the following lines contain a straightforward implementation,
        % filling up the rows one after each other:
        for i=1:n
            A(i,:) = 1./(i:i+n-1);
        end

    case 2
        % the following lines contain a straightforward implementation,
        % filling up the columns one after each other. Because of the way
        % matrices are stored in Matlab, this code turns out to be faster
        % than filling up the rows.
        for i=1:n
            A(:,i) = 1./(i:i+n-1)';
        end
        
    case 3
        % Here the division of the entries is performed only after the
        % matrix is already built. For small matrices this method can be
        % faster than the first method, but it usually becomes slower for
        % larger matrices.
        for i=1:n
            A(i,:) = i:i+n-1;
        end
        A = 1./A;

    case 4
        % The same as method 3 but filling up the matrix by its columns.
        % This method is surprisingly competitive and often faster than the
        % first 3 methods for large matrices.
        for i=1:n
            A(:,i) = (i:i+n-1)';
        end
        A = 1./A;

    case 5
        % Yhe idea of the following lines is to perform all necessary
        % arithmetic operations before building up A. This leads to a much
        % smaller number of operations that have to be performed (2n-1
        % versus n^2).
        entries = 1./(1:2*n-1);
        for i= 1:n
            A(i,:) = entries(i:i+n-1);
        end
        
    case 6
        % The same idea as method 4, but building up A by its columns. This
        % is by far the fastest method (and for large matrices can be even
        % faster than the built in function hilb).
        entries = 1./(1:2*n-1)';
        for i= 1:n
            A(:,i) = entries(i:i+n-1);
        end
        
    otherwise
        error('Undefined method for computing the Hilbert matrix');
end