The energy balance of a fixed bed system can be modelled by a Ist or 2nd order evolutionary differential equation. In the first order case, the equation exhibits shock behaviour. Standard numerical techniques fail to resolve the shock adequately without a large number of finely spaced points in the spatial dimension. This leads to the use of adaptive grid methods to attempt to reduce the computational requirements for fixed bed system simulation.This paper describes an adaptive grid algorithm and the derivation of the nonuniform numerical approximations for the derivative terms in the model equation and the boundary conditions. Two time stepping algorithms are compared : an implicit scheme based on the Crank-Nicholson approximation coupled with a Newton-type solver and the explicit Runge-Kutte-Gill 4th order method. Results are presented for the fixed bed system with appropriate initial and boundary value conditions. A discussion of the effects of the interpolation method used (required for generation of new grids adaptively) is presented, indicating the most suitable choice for a problem of this type.