Uncertainty quantification of practical engineering applications using the intrusive spectral stochastic finite element methods (SSFEM) may involve solving a system of linear equations in the order of billions of unknowns. Therefore, in this thesis the intrusive polynomial chaos expansion (PCE) based two-level domain decomposition (DD) algorithms for stochastic partial differential equations (PDEs) are extended to handle high resolution numerical models using an in-house scalable parallel solvers toolkit. First, the attention is given to facilitate the numerical simulation of the stochastic PDEs with a large number of random variables to address the so-called curse of dimensionality issue. Second, for three-dimensional coupled stochastic PDE systems such as equations of linear elasticity, the extended wirebasket-based coarse grid is developed to improve the performance and overcome the scalability issues of the DD based iterative solvers with a vertex-based coarse grid. Third, the developed DD solvers for the SSFEM are coupled with FEniCS deterministic finite element assembly routines in order to reduce the coding required for the implementation and generalize the application of these solvers to a variety of PDEs using FEniCS. Fourth, the intrusive SSFEM with scalable DD solver is shown to outperform the non-intrusive SSFEM with the sparse grid quadrature for a stochastic PDE with the non-Gaussian random variables. This highlights the advantages of intrusive approach and demonstrates the necessity of the scalable parallel solvers for uncertainty quantification. This thesis also elaborates on the HPC implementational aspects of the DD solvers for SSFEM. Three-level nested sparse iterative solvers, which employ an efficient DD based preconditioner at each level, are used to simulate two and three-dimensional scalar and vector-valued stochastic PDEs. The random system parameters and the solution process are modeled as a non-Gaussian stochastic process characterized by using up to 25 random variables with a third-order expansion resulting in 3276 PCE terms. Numerical and parallel scalabilities of these algorithms are investigated employing large-scale high-performance computing clusters with MPI, PETSc and FEniCS libraries. For all these scalability measurements, we report the results for different numbers of random variables and different orders of PCE to measure the sensitivity of the algorithms to the stochastic parameter.