function [tp,yp] = colloc4(N,Np);
%
%   function [tp,yp] = colloc4(N,Np);
%
% This function implements the collocation method
% for the problem y'' + y = -t subject to homogeneous
% Dirichlet boundary conditions. N denotes the number 
% of collocation points, and Np denotes the number of
% points in the output arrays for xp and up, for plotting
% the solution. The basis functions are (1-t)*t^k.
%

t  = linspace(1.0/(N+1.0), N/(N+1.0), N)';
tp = linspace(0,1,Np)';

% Create the coefficient matrix and the right-hand side

A = zeros(N,N);
b = -t;

for k=1:N
    A(:,k) = k*(k-1)*t.^(k-2) - (k+1)*k*t.^(k-1) + t.^k - t.^(k+1);
end;

% Solve for the coefficients

alpha = A \ b;

% Compute the function values

yp = alpha(N) * ones(Np,1);

for k=N-1:-1:1
    yp = yp .* tp + alpha(k);
end;

yp = yp .* tp .* (1.0 - tp);

maxer = max(abs(yp + tp - sin(tp) / sin(1.0)))