This dissertation focuses on the formulation and analysis of various fitted numerical methods for the solution of singularly perturbed parabolic diffusion-type problems characterized by boundary layers, both in the presence and absence of small temporal delays and spatial shift parameters. These types of problems arise frequently across a wide range of modeling contexts, including biological and chemical reaction processes, heat and mass transfer, control theory, neuronal variability, and population dynamics and epidemiology. Because the highest-order derivative term in the problems under consideration is multiplied by a small parameter e, 0 < e 1, the solutions exhibit multi-scale behavior, including boundary layers, interior layers, oscillations, and sharp gradients confined to narrow regions of the spatial domain. These features pose significant numerical challenges. Consequently, standard numerical methods such as finite difference, finite element, and finite volume methods applied on uniform meshes generally fail to produce accurate approximations unless prohibitively fine discretizations are employed. Such refinements, however, significantly increase computational cost and may also amplify round-off errors. To overcome these computational challenges, a Taylor series approximation is employed to transform the delay and advance terms, thereby reducing the original time-delay and spatial-shift problem into an asymptotically equivalent singularly perturbed problem. This reformulation places the problem within a standard singular perturbation framework while preserving its essential asymptotic structure. Based on this reformulated framework, we construct and rigorously analyze several fitted numerical methods for the classes of problems under consideration. The proposed approaches employ uniform meshes for temporal discretization and both uniform and piecewise-uniform Shishkin meshes for spatial discretization, depending on the nature of the problem. Specifically, the dissertation develops and analyzes secondand fourth-order fitted spatial discretization methods, including nonstandard finite difference schemes, fitted linear multistep methods, exponentially fitted upwind schemes, midpoint-upwind schemes, and hybridized numerical schemes, together with Crank–Nicolson and implicit Euler time-stepping strategies. The validity of the proposed methods is established through rigorous stability and parameter-uniform convergence analyses and is further supported by extensive numerical experiments. The maximum absolute errors and corresponding orders of convergence are computed and presented in tabular and graphical forms for clarity and comparison. The results demonstrate excellent agreement between theoretical predictions and numerical experiments, confirming that the proposed methods exhibit superior accuracy and convergence behavior compared with several existing approaches in the literature.