A mathematical model describing the behavior of elution peaks in preparative liquid chromatography for multicomponent separation has been solved numerically with the application of the perturbation and modified collocation methods. With application of such a collocation method on the time and axial domains, the governing equations are transformed into a set of recursive, algebraic equations that can be solved efficiently. The simulation results excellently match the reported experimental elution profiles resulting from the separation of a binary mixture of dipeptides on a reversed phase column in the literature. The results from the binary system can be extrapolated to multicomponent separation. For a given system of stationary and mobile phases, the model simulation predicts that the sample volume, concentrations of solutes in the sample, flow velocity, and particle diameter influence the performance of the elution peaks considerably. When samples having a given composition are chromatographied on a prepacked column, the throughputs and optimal injection conditions for achieving the touching-band separation are predictable. The simulation results show that the longer the column length, the larger the production rate, and the column length required,to achieve the touching-band separation increases with sample volume. The production rates for every component also increase with sample volume injected for the preparative proposes.