function output2 ( f, indx, n, nl, node, nprint, nu, problem, ul, ur, xl, xn, ...
  xr )

%% OUTPUT2 compares the computed and exact solutions.
%
%  Modified:
%
%    01 November 2006
%
%  Parameters:
%
%    Input, real F(NU), the solution of the linear equations.
%
%    Input, integer INDX(1:N+1).
%    For a node I, INDX(I) is the index of the unknown
%    associated with node I.
%    If INDX(I) is equal to -1, then no unknown is associated
%    with the node, because a boundary condition fixing the
%    value of U has been applied at the node instead.
%    Unknowns are numbered beginning with 1.
%    If IBC is 2 or 4, then there is an unknown value of U
%    at node 0, which will be unknown number 1.  Otherwise,
%    unknown number 1 will be associated with node 1.
%    If IBC is 1 or 4, then there is an unknown value of U
%    at node N, which will be unknown N or N+1,
%    depending on whether there was an unknown at node 0.
%
%    Input, integer N.
%    The number of subintervals into which the interval
%    [XL,XR] is broken.
%
%    Input, integer NL.
%    The number of basis functions used in a single
%    subinterval.  (NL-1) is the degree of the polynomials
%    used.  For this code, NL is fixed at 2, meaning that
%    piecewise linear functions are used as the basis.
%
%    Input, integer NODE(NL,N).
%    For each subinterval I:
%    NODE(1,I) is the number of the left node, and
%    NODE(2,I) is the number of the right node.
%
%    Input, integer NPRINT.
%    The number of points at which the computed solution
%    should be printed out when compared to the exact solution.
%
%    Input, integer NU.
%    NU is the number of unknowns in the linear system.
%    Depending on the value of IBC, there will be N-1,
%    N, or N+1 unknown values, which are the coefficients
%    of basis functions.
%
%    Input, integer PROBLEM, indicates which problem to be solved.
%    * 1, u = x, p=1, q=0, f=x, u(0)=0, u'(1)=1.
%    * 2, u = 2*(1-cos(0.5*pi*x))/pi, p=1, q=0,
%    f = -0.5*pi*cos(0.5*pi*x) + 2*sin(0.5*pi*x)*(1-cos(0.5*pi*x)/pi
%    u(0) = 0, u'(1)=1.
%
%    Input, real UL.
%    If IBC is 1 or 3, UL is the value that U is required
%    to have at X = XL.
%    If IBC is 2 or 4, UL is the value that U' is required
%    to have at X = XL.
%
%    Input, real UR.
%    If IBC is 2 or 3, UR is the value that U is required
%    to have at X = XR.
%    If IBC is 1 or 4, UR is the value that U' is required
%    to have at X = XR.
%
%    Input, real XL.
%    XL is the left endpoint of the interval over which the
%    differential equation is being solved.
%
%    Input, real XN(1:N+1).
%    XN(I) is the location of the I-th node.  XN(1) is XL,
%    and XN(N+1) is XR.
%
%    Input, real XR.
%    XR is the right endpoint of the interval over which the
%    differential equation is being solved.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'Compare computed and exact solutions:\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '                X             U(X)          U(exact)\n' );
  fprintf ( 1, '\n' );

  for i = 1 : nprint

    x = ( ( nprint - i     ) * xl   ...
        + (          i - 1 ) * xr ) ...
        / ( nprint     - 1 );

    ux = u_exact ( x, problem );

    for j = 1 : n

      xleft = xn(j);
      xrite = xn(j+1);
%
%  Search for the interval that X lies in.
%
      if ( xleft <= x & x <= xrite )

        u = 0.0;

        for k = 1 : nl

          ig = node(k,j);
          iu = indx(ig+1);
          [ phii, phiix ] = phi ( k, x, xleft, xrite );

          if ( iu <= 0 )
            if ( j == 1 & k == 1 )
              u = u + ul * phii;
            elseif ( j == n & k == nl )
              u = u + ur * phii;
            end
          else
            u = u + f(iu) * phii;
          end

        end

        break;

      end

    end

    fprintf ( 1, '          %12f  %12f  %12f\n', x, u, ux );

  end