Flow around cylinders is observed in many engineering applications from laminar flow over sensors in petrochemical industries to turbulent flow in wind flow around buildings and electronic components. Two or more cylinders in proximity make the problem more complicated due to the interactions between shear layers, wakes and shedding vortices. The practical importance of cross flow over two proximate cylinder-like structures in nature and industry arises in flow around urban structures, high-rise buildings, towers, bridge piles, offshore structures, heat exchangers, electronic devices, etc. In fact, the staggered arrangement is the most common arrangement in cylinder pairs. The aims of this study are to perform the numerical investigations of three-dimensional flow properties and heat transfer around two square cylinders in staggered arrangement in the turbulent flow at Re=22000. In this study, the fluid flow and heat transfer on two cylinders with a square section in staggered arrangement for various spacing (P/D=1.5, 2, 3, 4, 5) are investigated with two and three-dimensional simulations. First the two-equation k-? SST model has been applied for numerical investigation of two-dimensional geometry and flow patterns, the mean and RMS of Nusselt, lift and drag coefficients has been determined. Secondly, the large eddy simulation has been applied for three-dimensional analysis. Grid and time step study for both two and three-dimensional simulations have been performed. Single body, periodic gap flow, periodic flow and synchronized vortex shedding patterns have been observed for different cylinder spacing (P/D). In low spacing, the lift and drag coefficients of cylinders differ from those of single cylinder values. In P/D=1.5, the drag coefficient was greater than that of single cylinder. In this case, the lift coefficient differs from zero value. It is observed that the value lift and drag coefficients approach to those of single cylinder values by increasing the distance between cylinders. In addition, the mean Nusselt number of each cylinder increases and approaches to those of single cylinder for higher spacing. In all distances except for P/D=2, the mean Nusselt number of downstream cylinder is greater than that of the upstream cylinder value. By comparing the results of the two and three-dimensional simulations, it was found that the two-dimensional simulations predict the large structures and flow patterns and trend of changes of mean parameters by changing the cylinder spacing relatively similar to those of three-dimensional results (LES). Two-dimensional simulations cannot predict the detail of real flow structures especially small scale ones. In addition, the predicted mean parameters such as drag coefficients for two- dimensional simulations are more than those of three-dimensional simulation. These differences are much higher for RMS of parameters, where it reaches to about twice as high as that of three-dimensional value. Key words : staggered square cylinders, numerical study, large eddy simulation, k-? SST model, flow structure, cylinder spacing, vortex shedding, heat transfer