Numerical simulation of fully developed laminar longitudinal forced convection heat transfer between triangular/square arrays of parallel cylinders immersed in a porous medium

This study numerically investigates fully developed laminar forced convection through triangular and square cylinder arrays embedded in a porous medium, introducing a novel analysis of submerged cylinders beyond conventional fluid-flow studies, with relevance to heat exchangers and energy storage applications. The modified Brinkman–Forchheimer–extended Darcy model is utilized to describe fluid flow in the porous medium. In the analysis of forced convection within flow channels, it is assumed that the fluid and solid phases are in local thermal equilibrium throughout the porous medium. The cylinder walls were subjected to two classical thermal boundary conditions, including constant wall temperature (CWT) and constant wall heat flux (CWHF). Taking advantage of the flow field’s symmetry, the simulation was performed on a reduced computational domain limited to the symmetry planes and the cylinder walls to minimize computational effort. The finite element approach was applied to discretize the governing equations subject to the given boundary conditions, and the solution was obtained using COMSOL Multiphysics software. Extensive comparisons with prior solutions verified the correctness of the solution method and yielded excellent consistency. The results are expressed in terms of how parameters such as porosity, ratio of the thermal conductivity of solid particles to that of the fluid ( Kr ), and the longitudinal pitch-to-radius ratio of the cylinders ( Spr ) influence the average Nusselt number. For every scenario analyzed, the triangular arrangement of cylinders produced a higher average Nusselt number than the square configuration. For instance, at ε = 0.6, Kr = 10 , and Spr = 1.6 , the triangular arrangement yields an average Nusselt number approximately 5.7 units greater than that of the square arrangement under the CWHF condition, while the corresponding increase under the CWT condition is only about 1.8. For ε = 0.5, Kr = 75 , and Spr = 1.6 , the corresponding values are 14.6 and 13.0 for the triangular and square arrangements, respectively. The results further indicate that, under the conditions investigated, the CWHF boundary condition consistently yields a higher average Nusselt number than the CWT boundary condition. With the CWHF boundary condition imposed, the mean Nu is roughly twice that found under the CWT scenario. The Nusselt number is found to decrease with increases in porosity and longitudinal pitch, whereas it increases with increasing Kr .