function x = xgrid(ax,bx,m,gridchoice)
%
% Specify grid points in space for solving a 2-point boundary value problem
% or time-dependent PDE in one space dimension.
%
% Grid has m-1 interior points on the interval [ax, bx].  
% gridchoice specifies the type of grid (see below).
% m+1 grid points (including boundaries) are returned in x.
%
% This file is a modified version of the original function taken from
% http://www.amath.washington.edu/~rjl/fdmbook/  (2007)

z = linspace(0, 1, m+1)';   % uniform grid in [0,1]

switch gridchoice

  case 'uniform'
     x = ax + (bx-ax)*z;

  case 'rtlayer'
     % Clustered near right boundary
     x =  ax + (bx-ax) * ((exp(-2*z) - 1)/(exp(-2)-1));
     
  case 'ltlayer'
     % Clustered near left boundary
     x =  ax + (bx-ax) * ((exp(2*z) - 1)/(exp(2)-1));
     
  case 'random'
     x = ax + (bx-ax)*sort(rand(m+1,1));  
     x(1) = ax;  
     x(m+1) = bx;

  case 'chebyshev'
     % Chebyshev extreme points
     x = ax + (bx-ax) * 0.5*(1 + cos(pi*(1-z)));
  end
