3D multiphase flow was analyzed in the area of Rijeka bay in the Adriatic Sea. The necessary morphology data in the range of interest of the coastal bottom area were collected and the spatial surface was created. The functionality of the 3D model was studied in the large area of the realistic stochastic structure of the bottom and the shore. The obtained results give a detailed view of the velocity fields in the horizontal plane of different depths. Numerical simulation was performed in the open source program OpenFOAM with Volume of Fluid (VOF) method using the Eulerian approach. For solving this problem interFOAM solver for two incompressible, isothermal, immiscible fluids was used. The resulting simulations showed dominant flow from the western coast of the Krk island to the eastern coast of the Istrian peninsula. Seawater enters the bay through the Middle gate and Tihi channel and exits the bay through the Great gate. The movement of pollutants from the underwater marine outfalls Delta at different depths and the spread of pollution in the case of a tanker disaster and pollution from the Urinj refinery are presented. This research has shown that using a VOF method can be successfully implemented for describing fluid motion in large areas such as bays and oceans.