Simulations of the quantum dynamics of the HF molecule immersed in a field of five overlapping, intense, linearly polarized, infrared laser pulses of subpicosecond duration are performed. The HF molecule, initially in its ground state, is modeled as a rotating oscillator interacting with a classical laser field via electric dipole interaction. Realistic potential and dipole functions are used. Optimal overlaps of the five laser pulses, as well as the optimal carrier frequencies of the laser pulses, are found which maximize the HF dissociation yield. A maximal yield of 45% in a single combined pulse is achieved using the best available potential and dipole moment functions. The optimal infrared multiphoton dissociation pathway for the HF molecule includes a series of the Delta v=1 vibrational-rotational transitions followed by a series of Delta v greater than or equal to 2 vibrational-rotational transitions. The latter is necessary as a consequence of the vanishing Delta v=1 transition moment around v=12. In the Delta v=1 regime, both P and R branch transitions are found to be important. The angular distribution of the dissociative flux is computed. Robustness of the results with respect to changes in the interatomic potentials, dipole functions and reduced mass, as well as to changes in laser pulse parameters (carrier frequencies, timings, phases, field amplitudes, and pulse durations) is investigated.