In spite of many advantages of using nanoparticles in convective heat transfer, there are still some hidden aspects of nanofluids regarding simulations. Pool boiling flows are complicated regarding analytical aspects, and the presence of particles can noticeably extend this complexity. One of the most important aspects of nanofluid pool boiling is concerned with the changes in nucleation site and bubble diameter due to particles deposition and surface roughness. To include these effects, new correlations are implemented as a user-defined function for nucleation site density and bubble departure diameter. On the other hand, the particles are introduced and tracked everywhere in the domain in the Lagrangian frame by using discrete model. As an application, a tube bundle with four tubes is considered with different orientation angles concerning each other and different pitch distances. Unsteady Eulerian two-fluid model in ANSYS-Fluent is employed to simulate the liquid and vapour flows in the computational domain. In this work, pool boiling nanofluid flow is numerically solved around a two-dimensional horizontal cylinder with 20 mm diameter and compared with experimental data. The nanofluid is consist of distilled water and aluminium oxide with a particle size of 38 nm. The good agreement is found, and further discussion regarding particles migration and deposition are presented. It is found that the percentage of deposition is dependent on heat flux and particle concentration. Also, heat transfer coefficient increases with expanding the horizontal distance between cylinders and then decreases to a fixed value.