The main objective of this thesis is to utilize the powerful approximation properties of spline basis functions for numerical solutions of engineering problems that arise in the field of fluid mechanics. Special types of spline functions, the so-called Fup basis functions, are used as representative members of the spline family. However, the techniques developed in this work are quite general with respect to the choice of different spline functions. The application of this work is twofold. The first practical goal is the development of a novel numerical model for groundwater flow in karst aquifers. The concept of isogeometric analysis (IGA) is presented as a unified framework for multiscale representation of the geometry, material heterogeneity and solution. Moreover, this fundamentally higher-order approach enables the description of all fields as continuous and smooth functions by using a linear combination of spline basis functions. Since classical IGA uses the Galerkin and collocation approach, in this thesis, a third concept, in the form of control volume isogeometric analysis (CV-IGA), is developed and set as the foundation for the development of a karst flow numerical model. A discrete-continuum (hybrid) approach is used, in which a three-dimensional laminar matrix flow is coupled with a one-dimensional turbulent conduit flow. The model is capable of describing variably saturated conditions in both flow domains. Since realistic verification of karst flow models is an extremely difficult task, the particular contribution of this work is the construction of a specially designed 3D physical model (dimensions: 5.66 x 2.95 x 2.00 m) to verify the developed numerical model under controlled laboratory conditions. As a second application, this thesis presents the development of a full space-time adaptive collocation algorithm with particular application to advection-dominated problems. Since these problems are usually characterized by numerical instabilities, the novel adaptive algorithm accurately resolves small-scale features while controlling the numerical error and spurious numerical oscillations without need for any special stabilization technique. The previously developed spatial adaptive strategy dynamically changes the computational grid at each global time step, while the novel adaptive temporal strategy uses different local time steps for different collocation points based on the estimation of the temporal discretization error. Thus, in parts of the domain where temporal changes are demanding, the algorithm uses smaller local time steps, while in other parts, larger local time steps can be used without affecting the overall solution accuracy and stability. In contrast to existing local time stepping methods, the developed method is applicable to implicit discretization and resolves all temporal scales independently of the spatial scales. The efficiency and accuracy of the full space-time adaptive algorithm is verified with some classic 1D and 2D advection-diffusion benchmark test cases.