A mathematical model is developed and numerical modeling is performed to solve a scientific and industrial problem in the field of studying mass transfer processes in the “fracture set - matrix” system in a carbonate reservoir. The model is presented in a two-dimensional formulation, which is characterized by the presence of two subsystems with different sets of filtration parameters. The IMPES (implicit in pressure and explicit in saturation) and IMPIS (implicit in pressure and saturation) integration methods are used to numerically solve the system. An implicit finite-difference scheme is constructed for equations regarding pressure, when solving it, the saturations are set as fixed. Based on the solution of this scheme, the pressure in the fracture system and matrix is determined, and then the system is solved with respect to saturations. The explicit and implicit finite-difference schemes are considered to determine saturations. For a task on an industrial scale, significant resources are required to solve it. A parallel algorithm is proposed to speed up calculations with involves decomposing the computational domain into a lattice of domains of equal cardinality. Block-parallel processing is used when solving the implicit scheme. A series of computational experiments are performed to test the approach for both IMPES and IMPIS methods. Based on the obtained calculations, a comparative analysis of explicit and implicit methods for calculating equations for saturations is carried out. The calculation results confirm the effectiveness of the proposed fully implicit parallel algorithm.