A high order adaptive method-of-lines package, BACOL, is developed for solving one dimensional parabolic partial differential equations. Collocation with a B-spline basis is used for the spatial discretization. An approximate solution is calculated in a piecewise polynomial subspace of degree p, and the spatial error estimate is obtained by using a second solution computed in a degree p + 1 piecewise polynomial subspace. BACOL controls both the spatial error and the time error. After each time step the spatial error is estimated, and if it is larger than the spatial error tolerance, an equidistribution principle is employed to refine the mesh. At the same time, the number of mesh points employed can be changed if necessary. The time integration is done by a differential-algebraic-eqation (DAE) solver, DASSL, which uses backward differentiation formulas. Modifications made to DASSL include replacing the original linear algebraic solver by the almost block diagonal system solver, COLROW and scaling the newton iteration matrix to avoid the large condition number generated by the index-1 DAE. Computational results, comparing with D03PPF, TOMS731, EPDCOL, indicate that BACOL is reliable and extremely efficient in dealing with problems having solutions with rapid variation.