Three-dimensional numerical simulation of Taylor gas bubbles as primary unites of slug flow patterns rising in non-Newtonian environments is performed in the context of Direct Numerical Simulation (DNS) of the governing equations, where the whole physics of fluid motions will be taken into account. State-of-the-art numerical tools are proposed to tackle the numerical challenges in the DNS study of this problem. E.g. a coupled level-set volume-of-fluid (CLSVOF) interface capturing method is used to solve the topological changes of the interface. Physical formulations are integrated with moving-mesh (MM) technique to decrease the computational cost of 3D simulations and adaptivemesh-refinement (AMR) technique to increase the local accuracy around the interface. The governing equations are solved using High-Performance Computing (HPC) parallel approaches. To the best of the authors’ knowledge, this is the first work dealing with three-dimensional direct numerical simulation of Taylor bubbles rising in non-Newtonian environments.