The meltinglike transition in unsupported Na N clusters (N 55, 92, 147, 181, 189, 215, 249, 271, 281 and 299) is studied by first-principles isokinetic molecular dynamics simulations. The irregular size dependence of the melting temperatures T m observed in the calorimetry experiments of Schmidt et al.[Nature (London) 393, 238 (1998)] is quantitatively reproduced. We demonstrate that structural effects alone can explain all broad features of experimental observations. Specifically, maxima in T m N correlate with high surface stability and with structural features such as a high compactness degree. DOI: 10.1103/PhysRevLett.94.233401 PACS numbers: 36.40.Ei, 31.15.Qg, 36.40.Sx, 64.70.Dv There is a fundamental interest in understanding the analog of the melting phase transition in small atomic clusters. Classical arguments predict that melting temperatures T m of small particles should show a monotonic decrease from the bulk limit as their sizes shrink. However, Jarrold and co-workers [1] have demonstrated that small gallium and tin clusters melt at temperatures higher than T bulk m . Also, calorimetry experiments by Haberland and co-workers [2] show that the size dependence of T m is not monotonic for Na N clusters in the size range N 50-350. The observed pattern of maxima and minima in the T m N curve cannot be fully explained either by electronic or geometric (icosahedral) shell closings. Despite much theoretical effort [3][4][5][6][7][8][9], no definite explanation for the calorimetry experiments has emerged yet. In this Letter, we show that molecular dynamics (MD) simulations, with an accurate (first-principles) representation of interatomic forces, reproduce and provide an explanation for the calorimetry results.We employ an orbital-free (OF) version of density-functional-theory (DFT) [10], which expresses the energy as an explicit functional of the valence electron density, and pseudopotentials to represent the ionic field acting on the valence electrons. The technical details of our OF-DFT implementation are explained in previous work [5,11]. Here we describe just the main improvements incorporated into our energy functional. Compared to Kohn-Sham (KS) [12] calculations, two new ingredients appear in OF-DFT: (i) an explicit expression for the functional dependence of the noninteracting electron kinetic energy on density, T s n , and (2) if pseudopotentials are employed, they must be local, that is, dependent only on spatial position and not on orbital angular momenta which are simply not defined in OF-DFT. In previous works [5], we employed the second-order gradient expansion for T s n [13] and the local pseudopotential of Fiolhais et al. [14]. Those pseudopotentials were adjusted to reproduce bulk properties within linear-response-theory and might not be optimal for cluster DFT studies. In this work, the T s n functional is given by [11] T s n T vW n T n ; (1)where T vW 1 8 R dr jrnj 2 n is the von Weizsäcker functional, kr 3 2 1=3 nr , k 0 f is the Fermi wave vector corresponding to a mean electron ...