In past decades, high performance computing has became a valuable tool in many fields of environmental science and technology to utilize computational power for better characterization of the complexity of environmental systems as well as predicting their evolution in time. In this work, a parallel computing technique is presented for the numerical simulation of two-phase flow processes in porous media. The Galerkin finite element method (FEM) is used to solve the initial boundary value problem arising from the underlying mathematical model. The PETSc package is utilized for parallelization of the computational task in both the global assembly of the system of linear equations and the linear solver. In order to parallelize the global assembly of the linear equation system, the overlapping domain decomposition method is used. The preset parallel FEM approach is realized within the framework of OpenGeoSys, an open source C++ finite element code for numerical simulation of thermal, hydraulic, mechanical and chemical processes in fractured porous media. The computational efficiency of the approach has been tested with three examples of increasing complexity, the five spot benchmark, dense non-aquaeous phase liquid infiltration into a inhomogeneous porous medium and a real-world application to the CO2 storage research site: Ketzin, in Germany.