Turbulence in viscoelastic flows is a fascinating phenomenon with important technological implications, e.g. drag reduction at high Reynolds numbers and increased mixing efficiencies at low Reynolds numbers. The dynamics of these flows have been extensively studied experimentally over the last seventy years and more recently, in direct numerical simulations (DNS). However, theoretical progress in viscoelastic turbulence has been hindered by the fundamental challenges posed by the need to account for both the velocity as well as the elastic deformation history, encapsulated in the positive-definite conformation tensor. Due to the positivity constraint, the latter tensor is not a vector space quantity and thus classical approaches used to quantitatively analyze turbulence in Newtonian flows cannot be directly extended to viscoelastic flows. This fundamental issue is addressed in the present thesis in two parts. Firstly, we develop a decomposition of the conformation tensor about a given base-state that respects the mathematical and physical nature of this quantity. Scalar measures to quantify the resulting fluctuating conformation tensor are developed based on the non-Euclidean Riemannian geometry of the set of positive-definite tensors. The three measures are (a) the logarithmic volume ratio of the conformation tensor with respect to the base--state conformation tensor (b) the squared geodesic distance of the conformation tensor from the base--state, (c) the geodesic distance of the fluctuating conformation tensor from the closest isotropic tensor. Secondly, we develop an approach to perturb the conformation tensor in a physically consistent manner. This approach is an alternative to the classical weakly nonlinear expansion of vector space quantities, and is thus termed the weakly nonlinear deformation. When specialized to linear perturbations, this approach reveals the correct Hilbert space structure for the linearized problem. Viscoelastic (FENE-P) channel flow DNS are developed and used to illustrate the theoretical framework: fully turbulent flow is used for the first part, and the nonlinear evolution of Tollmien-Schlichting waves are considered for the second part. Several important insights are gleaned from these simulations, demonstrating the efficacy of the proposed approach. The fundamental contributions in the present thesis pave the road for theoretical modelling and analysis of viscoelastic turbulence.