Numerical simulations of laminar forced convection heat transfer for reciprocating, two-dimensional channel flows are performed as a function of the penetration length, Womersley (alpha) and Prandtl (Pr) numbers. The numerical algorithm is based on a spectral element formulation, which enables high-order spatial resolution with exponential decay of discretization errors, and second-order time-accuracy. Uniform heat flux and constant temperature boundary conditions are imposed on certain regions of the top surface, while the bottom surface is kept insulated. Periodicity of velocity and temperature fields is imposed on the side boundaries, while the flow is driven by an oscillating pressure gradient. These sets of boundary conditions enable time-periodic solution of the problem. Instantaneous and time-averaged surface and bulk temperature distributions, and Nusselt number variations are presented. For high a flows, the temperature field is significantly affected by the Richardson's annular effect. Overall, forced convection increases by increasing the penetration length, alpha and Pr Corresponding steady-flow simulations are performed by matching the volumetric flowrate. For the limited parameter space investigated in this paper steady unidirectional forced convection is more effective than the reciprocating flow forced convection.