Is there an iron–iron bond in diiron nonacarbonyl?
tl;dr: There probably isn’t an Fe–Fe bond. We should probably stop drawing an Fe-Fe bond. Textbooks and lecturers should really stop teaching students that there is an Fe-Fe bond. But all that said, there isn’t one single satisfactory explanation for why there isn’t one…

Despite a lot of studies (and despite $\ce{[Fe2(CO)9]}$ being one of the earliest examples of a molecule with symmetrically bridging carbonyls) , the nature of bonding in $\ce{[Fe2(CO)9]}$ is still a little hard to get a clear picture of.
Whilst at a research level, most academics now accept that there is no significant bonding interaction between the two metal centres, this is not echoed at the undergraduate level where the Fe–Fe bonding model is still often taught.
Most textbooks either avoid the issue entirely (never providing a structure), or show the Fe–Fe bond (with some footnote that it might not be correct), and even those that do acknowledge the lack of the metal-metal bond (Elschenbroich being a good example) don’t really discuss how this works in practice, or how one would go about electron counting (a source of many a nasty tutorial problem).
Historical context
The simple answer is that there was never really any strong evidence for the presence of the Fe–Fe bond. As is often the case, it was adopted by several (well known) chemists at the time, and subsequent publications (and especially textbooks) propagated the idea, largely due to its convenience with respect to electron counting.
The earliest data that really looked at the nature of the Fe–Fe bond in $\ce{[Fe2(CO)9]}$  can be attributed to a crystal structure measured by Powell (2.46 Å)[1](later corrected by Cotton (2.53 Å)).[2]
Although these values didn’t necessarily prove the existence of an Fe–Fe bond (you can’t see a bond with X-ray, only the atoms, from which you can infer bonding), its presence was necessary to fulfil the 18-electron rule, and the values obtained were well in-keeping with other M–M bond lengths known.
The issue with the above argument (as convenient as it is) is that it’s impossible to tell whether the 2.53 Å represents a true bonding interaction or whether it’s simply a consequence of the bridging Fe–C(O)–Fe bonds pulling the metals close together. The iron itself is relatively contracted, so unlike in some other bridged systems, the steric problem of bringing the two metals together is relatively minimal.
The Beginning of Change
Although the Fe–Fe bond was fairly widely adopted, not all were in agreement.
Roald Hoffmann (of Woodward–Hoffmann fame) published several papers which appeared to contradict the presence of an Fe–Fe bond, showing that the interaction was in-fact “antibonding and repulsive”.[3]
Though not fully quantitative (simple Hückel calculations), Hoffmann was able to construct qualitative molecular orbital diagrams using molecular symmetry and isolobal/electronic fragments (from well behaved molecules) to piece together the likely bonding in $\ce{[Fe2(CO)9]}$. To date, his findings are largely still valid (though exact energies have now been calculated with far more accuracy).
Quantitative Computational Studies[#]
Modern computational studies agree that overall, the Fe–Fe interaction is net anti-bonding (though some bonding interactions have been found) with no genuine covalent Fe–Fe bond, as was initially proposed.
In Liddle’s Molecular Metal-Metal Bonds, McGrady provides a qualitative  molecular orbital diagram for the complex. The σ and σ* orbitals, as well as the δπ and δπ* orbitals (sandwiched between the former pair), are fully occupied and provide no net bonding character. However, the population of the Fe–Fe antibonding πδ* orbitals leads to a formal bond order of –2.[4]

This diagram correlates well to the results in one of the earliest computational studies on this system comparing various methods (DFT, MP2, SCF) by Reinhold.[5]
In the study, various properties of the molecules are calculated (M–M bond length, M–C bond lengths, bond angles etc.) and this data is compared to known data from crystal structures/other experimental work. They determine that only DFT gives accurate results (casting some doubt on earlier computational work using less rigorous methods).[##]
The molecular orbital diagram constructed (in agreement with the one from Ref [4] above) crucially shows a negative overlap of the Fe–Fe orbitals is observed, suggesting a direct repulsion. To quote directly from the paper:

in the iron systems, with more contracted metal orbitals, the direct metal–metal repulsion is relatively weak while the metal–bridge bonds are sufficiently strong

Other computational studies have looked to calculate properties of the complex, such as vibrational frequencies of any proposed Fe–Fe bond, or the NMR shieldings that would be present on the ligands with and without the presence of the Fe–Fe bond. This work is referenced in Comprehensive Organometallic Synthesis II (Vol. 6), but doesn’t seem to have lead to very useful results.
Conclusion
The take home message from this is that there really is unlikely to be an Fe–Fe bond in the way it is commonly taught to undergraduates, and certainly not for the purposes of electron counting.
A nice possible explanation is given by McGrady,[4] which provides some explanation for why the 18-electron rule isn’t important in this case (and at a level that an undergraduate could understand):

There are, in fact, many parallels between 34-electron $\ce{[Fe2(CO)9]}$ and square-planar 16-electron species, not least the tendency of both to react with nucleophiles. The breakdown of the 18-electron rule in square-planar complexes occurs because the vacant pz orbital fails to find an appropriate symmetry match among the linear combinations of ligand-based orbitals, and so is not used, either to form a bond or to accommodate metal-based lone pairs of electrons. [...]
In $\ce{[Fe2(CO)9]}$, there are two vacant $\mathrm{a_2''}$-symmetric metal-based combinations (primarily linear combinations of 4s and 4pz) but only one occupied ligand-based combination of the same symmetry. Thus, one of the two metal-based orbitals remains “unused” in the same sense that the pz orbital remains unused in square-planar species. This simple symmetry-based analysis suggests that there is no need to invoke an Fe–Fe bond in $\ce{[Fe2(CO)9]}$ simply because the metal centers do not conform to the 18-electron rule.


References &amp; Notes
[1]: Powell, H.M.; Ewens, R.V.G J. Chem. Soc 1939, 286
[2]: Cotton, F.A.; Troup, J.M. J. Chem. Soc. Dalton Trans. 1974, 800 
[3]: Lauher, J.W.; Elian, M.; Summerville, R.H.; Hoffman, R. J. Am. Chem. Soc., 1976, 98, 3219 
[#]: Full disclosure, I’ve read most of the papers that will follow at some point whilst trying to answer this in the past, but (big but) I’m not a computational chemist. I’m hoping someone else may be able to shed light on how valid these computational results actually are. 
[4]: McGrady, J.E. in Molecular Metal–Metal Bonds; Liddle, S., Eds.; Wiley-VCH: Weinheim, 2015 
[5]: Hunstock, E.; Mealli, C.; Jose-Calhorda, M.; Reinhold, J. Inorg. Chem. 1999, 38, 5053. Note that this isn’t a complete geometry optimisation, they start with some known parameters. New J. Chem. 1994, 18, 465 is a paper supposedly attempting the complete optimisation from scratch, but I can’t find it anywhere right now to access it
[##]: Even now there are a lot of issues with accurately calculating properties of organometallic complexes… it takes a vast vast amount of computational time due to the need for basis sets far expanded compared to those developed for simple organic systems