The thesis deals with the numerical modeling of nonreactive and nonlinearly sorbing solutes in groundwater and analysis of the effect of heterogeneity resulting from spatial variation of physical and chemical parameters on the transport of solutes. For this purpose, a hyƯbrid flux corrected transport (FCT) and central difference method based on operator-split approach is developed for advection-dispersion solute transport equation. The advective transport is solved using the FCT technique, while the dispersive transport is solved using a conventional, fully implicit, finite difference scheme. Three FCT methods are developed and extension to multidimensional cases are discussed. The FCT models developed are anlaysed using test problems possessing analytƯical solutions for one and two dimensional cases, while analysing advection and dispersion dominated transport situations. Different initial and boundary conditions, which mimic the laboratory and field situations are analysed in order to study numerical dispersion, peak cliping and grid orientation. The developed models are tested to study their relative merits and weaknesses for various grid Peclet and Courant numbers. It is observed from the one dimensional results that all the FCT models perform well for continuous solute sources under varying degrees of Courant number restriction. For sharp solute pulses FCT1 and FCT3 methods fail to simulate the fronts for advection dominated situations even for modƯerate Courant numbers. All the FCT models can be extended to multidimensions using a dimensional-split approach while FCT3 can be made fully multidimensional. It is observed that a dimensional-split approach allows use of higher Courant numbers while tracking the fronts accurately for the cases studied. The capability of the FCT2 model is demonstrated in handling situations where the flow is not aligned along the grid direction. It is observed that FCT2 method is devoid of grid orientation error, which is a common pr.