Direct numerical simulations (DNSs) have been carried out for single and multiple square jets issuing normally into a cross-flow, with the primary aim of studying the flow structures and interaction mechanisms associated with the jet in cross-flow (JICF) problems. The single JICF configuration follows a similar study previously done by Sau et al. (2004, Phys. Rev. E, 69, p. 066302) and the multiple JICF configurations are arranged side-by-side in the spanwise direction with a jet-to-jet adjacent edge distance (H) for the twin-jet case and an additional third jet downstream along the centerline with a jet-to-jet adjacent edge distance (L) for the triple-jet case. Simulations are performed for two twin-jet cases with H=1D,2D, respectively, and for one triple-jet case with H=1D, L=2D, where D is the jet exit width. Flow conditions similar to Sau et al. are considered, i.e., the jet to the cross-flow velocity ratio R=2.5 and the Reynolds number 225, based on the freestream velocity and the jet exit width. For the single jet in cross-flow, the vortical structures from our DNS are in good qualitative agreement with the findings of Sau et al. For the side-by-side twin-jet configuration, results have shown that the merging process of the two initially separated counter-rotating vortex pairs (CRVPs) from each jet hole exit is strongly dependent on the jet-to-jet adjacent edge distance H with earlier merging observed for the case H=1D. Downstream, the flow is dominated by a larger CRVP structure, accompanied by a smaller inner vortex pair. The inner vortex pair is found not to survive in the far-field as it rapidly dissipates before exiting the computational domain. These observations are in good agreement with the experimental findings in the literature. Simulations of the triple-jet in cross-flow case have shown some complicated jet-jet and jet-cross-flow interactions with three vortex pairs observed downstream, significantly different from that seen in the twin-jet cases. The evidence of these flow structures and interaction characteristics could provide a valuable reference database for future in-depth flow physics studies of laboratory experimental and numerical investigations.