Theory of melting lines

Our understanding of the three basic states of matter (solids, liquids and gases) is based on temperature and pressure phase diagrams with three phase transition lines: solid-gas, liquid-gas and solid-liquid lines. There are analytical expressions $P(T)$ for the first two lines derived on a purely general-theoretical thermodynamic basis. In contrast, there exists no similar function for the third, melting, line (ML). Here, we develop a general two-phase theory of MLs and their analytical form. This theory predicts the parabolic form of the MLs for normal melting, relates the MLs to thermal and elastic properties of liquid and solid phases and quantitatively agrees with experimental MLs in different system types. We show that the parameters of the ML parabola are governed by fundamental physical constants. In this sense, parabolic MLs possess universality across different systems.


INTRODUCTION
In statistical physics and thermodynamics, coexistence of different phases is set by the equality of thermodynamic potentials.This results in the Clausius-Clapeyron (CC) equation describing phase transition lines on the phase diagram: where dP dT is the slope of the transition line and S and V are entropy and volume of any two phases [1][2][3].
Eq. ( 1) readily gives important results for two transition lines, the sublimation and boiling lines [1,2].This involves two steps.First, the gas volume is much larger than volume of solids and liquids so that V 1 −V 2 becomes the gas volume obeying the gas equation of state (EoS).Second, assuming the constant latent heat q = T (s 1 −s 2 ) in Eq. (1) (s i are entropies per particle) and integrating gives a purely theoretical line of saturated vapor pressure along the sublimation and boiling lines as Eq. ( 2) is a textbook function for the sublimation and boiling lines [1,2].It is broadly considered to be in agreement with experiments [2,4], although cases of disagreement are discussed too [5].
In contrast, there is no analytical function P (T ) derived for the third, melting, line of a theoretical generality comparable to sublimation and boiling lines.This generality includes not relying, as was done in the past, on phenomenological fitting functions, arbitrary series truncations, assumptions of a particular system-specific melting mechanism, using computer modelling to calculate the melting line (ML), its parts or related data, choosing specific liquid structure or interactions and so on [3,[6][7][8].
We note that the ML is not limited by triple and critical points from above as the other two transition lines are.Continuing up to very high temperature where sys-tem properties change due to, e.g., ionisation, the ML extends over a much larger range of the phase diagram.
Understanding mechanisms of melting has a long history [3,6,7,9].Starting from 1913, Sommerfeld and Brillouin aimed to develop a thermodynamic theory of liquids based on phonons [10][11][12].This was around the time when phonon-based theories of Einstein and Debye were published [13,14] laying the foundations of the solid state theory.Sommerfeld and Brillouin considered that, differently from solids, liquids do not support shear stress and shear phonons [10][11][12].Brillouin [15] and independently Born [16] then defined melting as the point where solid's rigidity disappears.This elicited an objection from Frenkel [17]: liquids too support shear stress, albeit at high frequency.This is a valid point to which we return to later.Other mechanisms include semiphenomenological Lindemann-Gilvarry criteria [18,19] based on threshold atomic displacements a solid can tolerate before melting.These criteria continue to remain of interest [20], although they are known to agree with experiments in some systems only but not others [6,7].Other approaches include defect-mediated mechanisms with system-specific defects such as dislocations, surfaces and premelting [3,6,7].It is therefore unlikely that a microscopic theory based on a melting mechanism exists that is generally applicable.We also note that the earlier criteria above are one-phase approaches and focus on the solid phase only, whereas both solid and liquid phases are equally relevant for understanding melting which is set by the equality of thermodynamic potentials of the two phases as reflected in the CC equation (1) [1,3,6].
Faced with these historical issues, we can nevertheless aim for a general thermodynamic theory of the MLs which is independent of a particular melting mechanism and where system-specific effects contribute to thermodynamic potentials.This is of importance for general theory as well as applications: predicting phase diagrams is an active research extending to materials science, with the need to understand liquids emphasised [21,22].As discussed in the next section, developing a general theory of MLs faced, until recently, fundamental problems of liquid description.
Here, we develop a general two-phase theory of MLs and their analytical form.This theory predicts the parabolic form of the MLs for normal melting, relates the MLs to thermal and elastic properties of solid and liquid phases and agrees with experiments in different system types including noble, molecular, network and metallic.We show that the parameters of the ML parabola are governed by fundamental physical constants.In this sense, parabolic melting lines possess universality across different systems.

RESULTS AND DISCUSSION
A. The main idea: liquid energy The main reason for the absence of a general analytical function for the ML is related to fundamental problems involved in liquid theory.In this section, we recall the origin of these problems.This is followed by the summary of how these problems are overcome in the phonon theory of liquid thermodynamics.This theory gives an equation for the liquid energy which will be used in the next section to derive the differential equation for the ML.
Deriving the ML theoretically is faced with several fundamental problems, explaining why it was not achieved before.First, q for the sublimation and boiling lines is related to transitions between a condensed solid or liquid state and the uncondensed gas state.For this reason, q includes the cohesive energy of the condensed state and can be considered roughly constant.This does not apply to the ML because it separates the states which are both condensed.Hence one needs to explicitly know the liquid entropy in order to solve the CC equation (1) for the ML.Differently from solids and gases, liquid thermodynamic properties including entropy are not generally known from theory [1].Second, one needs to explicitly know the liquid V as a function of P and T and hence the liquid EoS in order to solve the CC equation (1).This EoS is not generally known from theory [1], differently from the gas EoS used to derive the boiling and sublimation lines from Eq. (1).Finally, V 1 − V 2 in Eq. ( 1) can't be approximated by V of either phase as is done when deriving the sublimation and boiling lines because solids and liquids are both condensed states with relative differences of V of about 10% only [6].
These problems, the unavailability of a general theory of liquid thermodynamic properties including entropy and EoS, are part of a more general fundamental problem of liquid theory.As discussed by Landau and Lifshitz (LL) [1], the simplifying features of solids and gases are small oscillations and weak interactions, respectively.Liquids have neither because atomic displacements are large due to flow and interactions are as strong as in solids.This fundamental problem is summarised to say that liquid theory has no small parameter [23].Because the interactions are strong, all thermodynamic properties strongly depend on the type of interactions and therefore are strongly system-dependent.This, conclude LL [1], rules out a generally applicable theory of liquid thermodynamic properties and their temperature dependence (this assertion was consistent with Landau's own view according to Peierls [24]).These fundamental problems extend to all theories based on explicitly considering liquid interactions and structure, the mainstream approach to liquid physics in the last century [25][26][27][28][29][30][31][32][33][34][35][36][37][38].This includes expanding interactions into short-range repulsive reference and attractive terms in simple models (see, e.g., [35][36][37][38]): the expansions and their coefficients are systemdependent and the results are not generally applicable.Apart from research, these problems impacted teaching [39].Importantly, these issues do not emerge in the solid state theory: this theory operates in terms of excitations in the system, phonons, rather than considers interactions and structure explicitly [9].As a result, the derived thermodynamic properties are applicable to all solids [1].
The liquid entropy in the CC equation ( 1) is particularly hard to approach head on: it includes vibrational and configurational terms, and its unclear how to balance the two [6,40,41].This entropy is not available, apart from simple models such as hard-sphere or van der Waals systems used to discuss liquids [42,43].These models give the specific heat 1,9,44,45] and therefore describe the ideal gas thermodynamically.This is far removed from real liquids where c v ≈ 3k B close to the ML [9,41,[44][45][46] which we are interested in and where c v is as large as in solids.
Our main idea to overcome the problem of solving the CC equation ( 1) for the ML is to write this equation in terms of energy instead of entropy.The liquid energy is known on the basis of the recent phonon theory of liquids.This theory [8,9,41,47,48] is based on the main premise of statistical physics: properties of an interacting system are governed by its excitations [1,49].In solids, these are collective excitations, phonons.The theory further focuses on energy E as the primary property in statistical physics [1], and this makes the liquid problem tractable (entropy can be found by integrating 1 T dE dT if needed).In view of long-lasting fundamental problems of liquid theory above, it is useful to briefly summarise the theory and its main result which we will use.
The key to addressing the fundamental small parameter problem stated by Landau, Lifshitz and Pitaevskii is that this parameter does exist in liquids and is the same as solids, small phonon displacements (this is grounded in microscopic liquid dynamics: between diffusive jumps enabling liquid flow, liquid particles vibrate around quasiequilibrium positions [17,40,41]).In an important difference to solids where the phase space available to phonons is fixed, this phase space in liquids is variable [9,41].
Experiments have shown that propagating phonons in liquids and solids are remarkably similar in terms of dispersion curves [9,41,[50][51][52][53][54][55][56][57][58][59].The phonon theory of liquid thermodynamics [8,9,41,47,60] further observes that transverse phonons in liquids exist above the gap in the reciprocal space k g = 1 cτ [41,[61][62][63] only, where τ is the liquid relaxation time identified with the average hopping time of molecules from one quasi-equilibrium point to the next and c is the speed of sound [9,17,40].This is equivalent to propagating phonons (in a sense of frequency exceeding the decay rate) above the hopping frequency ω F = 1 τ [9,64] as envisaged by Frenkel originally [17].This gives the variability of the phonon phase space in liquids.There are several ways to measure the hopping frequency ω F , including using the approximate relation ω F = G η , where η is viscosity and G is the high-frequency shear modulus [17,40].
The energy of phonons excitations in the liquid is the sum of the longitudinal phonon energy and the energy of propagating gapped transverse phonons [9,41,64].Differently from solids, liquids have another type of excitations: these are diffusing atoms which enable liquid flow.Adding the energy of this motion to the phonon energy gives the liquid energy as [8,9,41,47,60] where ω D is the maximal Debye frequency that sets the high-temperature limit for ω F .Here and below, k B = 1.Eq. ( 3) predicts that liquid specific heat at the melting point where ωF ωD is small (see next section for a detailed discussion).Therefore, Eq. ( 3) predicts that liquid c v is close to that in solids in this regime.
As ω F increases with temperature, E l changes from 3N T to 2N T when ω F reaches ω D at the Frenkel line where transverse phonons disappear from the liquid spectrum [65] (recall that the absence of transverse phonons was seen as a general property of liquids by Sommerfeld and Brillouin [10][11][12]).This change of E l gives the decrease of c v from 3 to 2. This decrease is universally seen in liquids and is related to the reduction of the phonon phase space [8,9,41,46,47].
Eq. ( 3) and its extensions to quantum effects [9, 41] have undergone independent "detailed and rigorous" tests against experiments for different types of liquids in a wide range of T and P [8,47].Noting that the theory is falsifiable and pursuing this falsifiability, these tests concluded that the theory is supported by experiments in a number of ways.

B. Differential equation for the melting line
We write the CC equation ( 1) as where indexes l and s refer to liquids and solids from now on, differentiate (5) and multiply the result by T : The derivatives dV l,s dT in Eq. ( 6) are along the ML.Below we will relate the ML and its properties to thermal and elastic characteristics of the solid and liquid phases.This includes the thermal expansion coefficient α and the bulk modulus B calculated from separate experimental data.This data enables us to calculate dV dT in Eq. ( 6) along the ML directly.However, it is useful to relate  6) gives The right-hand side is the difference between the constant-pressure heat capacities C p of solids and liquids across the ML.C p is related to the constant-volume heat capacity C v = ∂E ∂T , where E is energy and the derivative is at constant volume, as [1] Using this in Eq. ( 7) gives the equation depending on energies E l and E s of both phases as We now recall the energies E s and E l : E s = 3N T and E l in Eq. ( 3).Using E s and E l in Eq. ( 9) gives All terms in Eq. ( 10) scale with N and V .However, the first term on the right-hand side is much smaller than all other terms and can be dropped because ωF ωD ≪ 1 at the ML.Recalling ω F = G η gives ωF ωD at the ML as where T m is the melting temperature and η 0 is the hightemperature limiting value of η.
The smallness of ωF ωD in liquids at the ML is in perfect agreement with the experimental c v of many liquids.According to Eq. ( 3), ωF ωD ≪ 1 gives the liquid c v close to its solid value, c v = 3.This c v is seen in liquids close to the ML in a wide range of temperature and pressure and includes liquids with different structure and bonding types: noble, metallic, semiconducting and molecular [8,9,41,[44][45][46][47].This experimental fact supports the smallness of the energy difference term E l − E s in Eq. ( 9) and of ωF ωD in Eq. (10).Physically, this smallness of the energy difference E l − E s at the ML is due to solids and liquids being both condensed states of matter where the characteristic scales of energy and interatomic separation are set by the Rydberg energy and Bohr radius involving fundamental physical constants [71].The closeness of most important sys-tem properties then follows [71].This does not apply to phases separated by the other two transition lines, sublimation and boiling lines, because they involve the non-condensed gas phase where the Rydberg energy and Bohr radius do not operate.
This last point makes an interesting connection to our discussion in the previous section where we noted one of the problems of solving the CC equation ( 1) for the ML: it can not be assumed that the latent heat q includes a large contribution from the cohesive energy of solids (liquids) as is the case for sublimation (boiling) lines.Indeed, solids and liquids are both condensed phases, and q is related to the system-specific difference of cohesive energies of the two states.In our theory here, the condensed nature of solids and liquids is turned into an advantage: the closeness of energies of these states in Eq. ( 9), E l and E s , the consequence of both states being condensed, simplifies the theory.
For ωF ωD ≪ 1, Eq. ( 10) reads This second-order equation does not suggest a link to a second-order transition [1].We are considering the firstorder melting transition with parameters which enter in Eq. ( 12) through the entropy derivative in Eq. (8).
As set out in the previous section, Eq. ( 12) eliminates the problem of knowing the liquid entropy and depends on thermal and elastic properties of solids and liquids only.At the same time, Eq. ( 12) is quite general and follows directly from the CC equation ( 1) using one assumption only: the equality of c v (or energy) of solids and liquids across the ML.As discussed earlier, this equality is backed up by theory and numerous experiments.
Eq. ( 12) predicts the function P (T ) of the ML.The parameters of this function are set by the change of thermal and elastic properties of solids and liquids on crossing the ML.This naturally makes our theory a two-phase theory and different to earlier one-phase theories which considered how the properties of solids change on melting and left liquids out from consideration.These theories include those of Born, Brillouin, Lindemann-Gilvarry and so on [3,[6][7][8][15][16][17][18][19].These one-phase approaches fo-cusing solely on solids came about because solids were understood whereas liquids weren't as discussed in the previous section.However, there is no physical reason to single out solids to understand the transition between solids and liquids: both phases are equally relevant to understand this transition which is set by the equality of thermodynamic potentials of the two phases [1,6].Recall that the two-phase approach is implicit in deriving the sublimation and boiling lines where the latent heat q = T ∆s reflects the entropy change between two phases and enters the analytical form (2) as a parameter.
Eq. ( 12) provides a relation between different properties: thermal and elastic properties of solids and liquids on one hand and melting lines on the other hand.Providing such relations is viewed as an essence of a physical theory [73].
We note in that Eq. ( 3) assumes that the energy of each contributing phonon is given by T in the harmonic approximation.Accounting for phonon anharmonicity modifies E l and E s as [9,41]: where setting the hopping frequency ω F = 0 in the liquid corresponds to the solid and E a l = E s l as expected.Using E a l − E a s in the energy difference in Eq. ( 9) and ωF ωD ≪ 1 as before contributes the extra term 3N (α l − α s ) to the right-hand side of Eq. (12).Using the Grüneisen parameter γ = V αB Cv ≈ V αB 3N [7], this term is This approximately doubles the term on the right-hand side of Eq. (12).

C. Solution, its properties and the parabolic form of the melting line
Written in terms of functions a, b and d, Eq. ( 12) is where r is the ratio of the liquid and solid volumes and is typically about 1.1 [6] (the relative volume increase at melting is typically 10-15% for noble and molecular systems, 1-5% in monoatomic metals and semiconductors including negative values and 10-20% for binary salts [6]).
In general, a, b and d depend on P and T because α, B and r do.In this case, the ML comes from solving Eq. ( 14) numerically.a(P, T ), b(P, T ) and d(P, T ) can be taken from (a) experiments or modelling data including the EoS or (b) combining this data at low pressure with models predicting the variation of α and B with P and T [7].In the latter case, analytical solutions of Eq. ( 14) might be possible depending on the model.
To analytically solve Eq. ( 14) in order to check theory predictions, we make two approximations.These approximations are independent and sequentially interchangeable.First, we assume that α, B and r do not substantially change for sufficiently small variations of P and T and can be approximately considered constant [7].This assumption is analogous to considering the constant latent heat q in the derivation of the boiling and sublimation lines from the CC equation [1,2].In liquids, α and B depend on both viscous and elastic components of motion [17,74].The viscous component is governed by the activation energy that needs to be surmounted by the jumping molecule, U .At small pressure, U = 8πGr∆r 2 , where G is the liquid high-frequency shear modulus, r is the cage radius and ∆r is the temperature-induced fluctuating increase of the cage size that enables the molecule to escape the cage [17,40].Pressure typically increases U (unless there is an anomalous structural change from, for example, covalent to molecular liquid).This increase is the extra work needed to expand the cage from radius r to r + ∆r, 4πr 2 ∆rP .Then, U becomes U = 4πr∆r(G∆r + P r) The internal elastic and external pressure effects become comparable when P = G ∆r r .Taking ∆r ≲ r gives P on the order of GPa in noble and molecular systems and 10 GPa in stiffer systems such as Fe.A similar estimation of the variation of B due to the elastic component can be done by writing B = B 0 + ∂B ∂P P , where B 0 is the bulk modulus at zero pressure and ∂B ∂P is typically 4-6 for different minerals and oxides [7].In these systems, high-temperature α typically increases by about 3% over 100 K at ambient pressure.
We add two remarks regarding variation of α, B and r in Eq. ( 14).First, pressure and temperature have competing effects on α and B (pressure usually increases B and decreases α, whereas temperature does the opposite).This gives slower variation of α and B along the ML if pressure along the line increases with temperature as it usually does.Second, Eq. ( 14) contains V l and V s in addition to α and B, however V l and V s enter Eq. ( 14) as the ratio r = V l Vs which varies slower than the volumes themselves.For example, this ratio is constant along the ML for the inverse power potentials [75] (see Section 2 F).
Considering constant α, B and r gives constant a, b and d in Eq. ( 14).The solution of Eq. ( 14) is then where P 0 is the integration constant which can be used to fix the location of the triple point for a particular system and we set the second integration constant in the argument of cos to zero.We consider the most general case taking phonon anharmonicity into account.As discussed in the previous section, this anharmonicity doubles d in Eq. ( 14).Using a, b and 2d in Eq. ( 14) gives ∆ = 2bd − a 2 4 as In common non-anomalous cases where α l > α s , B s > B l (the reason these inequalities often apply is discussed later) and with typical r, √ ∆ is real.This is confirmed by using the actual values of α and B for liquids and solids [74] and a range of r in Eq. (17).
Eq. ( 16) makes several predictions regarding the ML.First, below we will consider normal cases where α l > α s , B s > B l and r > 1.In this case, the coefficient in front of the linear term in Eq. ( 16), a 2b , is positive according to Eq. ( 14).Second, the last term in Eq. ( 16) with real √ ∆ in normal systems gives the increase of the slope of the P (T ) melting line, dP dT , with temperature.We will discuss experimental data showing this increase below.However, we note that a combination of anomalous r, α and B in Eq. ( 17) can give imaginary √ ∆.In this case, cos √ ∆T in Eq. ( 16) becomes cosh |∆|T , resulting in the decrease of dP dT with temperature as seen in anomalous cases such as Ce [8].dP dT in Eq. ( 16) can also become negative (with a 2b being either positive or negative if α l < 0 in Eq. ( 14)) as seen in anomalous systems such as water or in systems with re-entrant melting [3,76].dP dT can also change its sign, corresponding to a non-monotonic behavior of the ML [3].These different scenarios of anomalous melting are described by the theory, however the detailed discussion of these scenarios is outside the scope of this paper.Here, we discuss commonly encountered normal melting where dP dT is positive and increases with temperature.
For small √ ∆T , Eq. ( 16) reads In the temperature range set by the approximation above, Eq. ( 18) predicts that the melting line is a parabola.This is an attractive function in view of its simplicity and common occurrence in physics.For a given pressure, Eq. ( 18) predicts the temperature of melting.
Fourth, the second term in Eqs. ( 16) and ( 18) predicts the value of the linear slope contributing to the overall dP dT of the ML and the last term predicts the deviation from linearity.Further insight into these values follows from making the second approximation to a, b and d in Eq. ( 14).α l is often substantially larger than α s and B l is substantially smaller than B s .This is because the temperature and elastic response includes both elastic and viscous components in liquids and only an elastic component in solids [17,74] (for example, typical ranges of α l and α s are (10 −5 − 10 −3 ) K −1 and (10 −6 − 10 −5 ) K −1 , respectively).This gives a = 2α l r r−1 , b = r B l 1 r−1 in Eq. ( 14), a 2b = α l B l and ∆ 2b = r 2(r−1) α 2 l B l .We can check these approximations to a 2b and ∆ 2b by using α and B in liquids and solids [74] and the typical value r = 1.1 [6].We find that α l B l underestimates a 2b by a factor of 1.4 on average and r 2(r−1) α 2 l B l overestimates ∆ 2b by 1.7 on average while the order of magnitude of a 2b and ∆ 2b is not changed as a result of these approximations.Then, Eq. ( 18) reads The ratio of quadratic and linear terms contains α l T which might suggest that the quadratic term is small compared to the linear term.However, Vs is close to 1. Therefore, the quadratic term does not need to be small and should not be viewed as a next-order correction.This is consistent with the CC Eq. ( 1): V l ≈ V s gives large dP dT .In our theory, this large dP dT is described by large quadratic term in Eq. ( 19).As the last remark in this section, we observe that equating the liquid and solid energies across the ML in Eq. ( 9) (setting the small first term on the right-hand side in Eq. ( 10) to zero) meant that T cancelled out from all terms in Eq. (10), resulting in Eq. ( 14).If instead we chose not to use the equality of E l and E s and consider the right-hand side of Eq. ( 7) as a new function ∆C p , T would not cancel in Eq. ( 14), and d would become d = ∆Cp T (V l −Vs) .This would alter Eq. ( 14) in two respects.First, this d would include an additional and unknown thermodynamic property, ∆C p , whereas the current Eq. ( 14) sets the ML in terms of thermal and elastic properties only.Second, even if we assumed ∆C p and V l − V s to be constant alongside constant a and b in Eq. ( 14), d ∝ 1 T would mean that Eq. ( 14) could not be integrated in a closed form (the solution would be an integral of a combination of special functions).Hence, using a physical insight regarding the equality of E l and E s significantly simplifies mathematical description and gives the parabolic form of the ML (18).

D. Further properties of the melting line and its relation to fundamental physical constants
The linear term in Eq. ( 19), α l B l T , gives the following physical interpretation of the ML in terms of generic thermal and elastic effects.Lets consider a liquid just below the ML at temperature T .A small temperature increment ∆T gives the relative volume increase This combination of thermal and elastic effects increases the melting pressure at temperature T + ∆T by the amount ∆P in comparison to what this pressure would have been in the absence of these effects.This corresponds to the slope ∆P ∆T = α l B l from Eq. ( 20) and the linear term in Eq. (19).
The same effect operates in the solid phase (our second approximation dropped this effect): negative temperature increment ∆T in the solid above the ML gives the relative volume decrease ∆V V = α s ∆T , generating extra pressure ∆P = α s B s ∆T and corresponding to the slope α s B s .The actual linear slope of the ML and dP dT are set by the interplay of thermal and elastic effects in both liquid and solid phases according to Eq. ( 12) or (14).
The coefficients in front of the linear and quadratic terms of the ML parabola (19) are set by α l and B l and, in a more general case of Eq. ( 14), by α and B of both solids and liquids.Whereas α and B vary in different systems, their characteristic magnitude is set by fundamental physical constants.Indeed, B is set by the density of electromagnetic or bonding energy , where a 0 is the interatomic separation [71].A shown by Weisskopf and Bernstein [77], the magnitude of α is set by the inverse of The scale of E b and a 0 is set by the Rydberg energy and the Bohr radius which, in turn, are fixed by the fundamental physical constants including the Planck constant and the electron mass m e and charge e.As a result, we have [71]: Eqs. ( 21), together with Eq. ( 19) or (14), imply that (a) the parameters of the ML parabolas are governed by fundamental constants and (b) variation of ML parameters in different systems is due to system-specific proportionality coefficients in Eq. ( 21).Hence, parabolic MLs (19) can not be too different in different systems and in this sense possess universality.In the next section, we will see that parameters of experimental MLs in different systems indeed fall in the characteristic range set by typical α and B. We note that pressure typically increases B and decreases α, whereas temperature has the opposite effect.As a result, the product of B and high-temperature α varies slower with pressure or temperature than α or B. In liquids, α l B l is in the range (3-16) MPa K as we will later see in Table 1.In solid minerals, αB at ambient pressure and high T is in a similar range (4-7) MPa K [7].Slow variation of αB might be useful for developing models of temperature and pressure variations of a, b and d in Eq. ( 14) and for using Eq. ( 19) where αB features.
Eqs. (21) help explain this slower variation of αB: B = E b a 3 and α = 1 E b imply that E b cancels in the product αB.As a result, this product becomes slower-varying with P and T and in different systems.

E. Comparison to experiments
To compare our theory to experiments and test its general applicability, we consider systems with different structure and bonding types: noble Ar and He, molecular H 2 , network H 2 O and metallic Fe and In.Fig. 1 shows the MLs for these systems.Consistent with experiments, Eqs. ( 16) and ( 18) predict the increase of the slope of the ML with temperature.We note that low-pressure water and its ML have anomalies affected by coordination changes [78].The ML in Fig. 1 is at GPa pressures where anomalies disappear and water properties become similar to those in simple systems such as Ar [46].
We now address the quantitative agreement between theory and experiment.In view of approximations made, this agreement is expected to be in terms of the right magnitude of ML properties rather than a close match.
We first check whether the first linear term in Eq. ( 19) contributing α l B l to the slope of the ML gives a sensible agreement with experiments.Recall that our first approximation to derive Eq. ( 19) is related to the low pressure range where a, b and d were considered constant.Hence we calculate the experimental slope of MLs, ∆P ∆T , in the low-pressure range 0-3 GPa as discussed earlier and shown in Fig. 1b.In much stiffer Fe, this range is below 100 GPa as discussed earlier and is considered in the range 0-30 GPa.Some MLs in Fig. 1b are close to linear and others are less so.In the latter case, the calculated slope is an average with a characteristic value.We compare the experimental ∆P ∆T to 2α l B l in Table 1.The ratios between experimental and theoretical slopes of the ML in Table 1 are about 2-3 for Ar, He, H 2 O and Fe and 5-6 for H 2 and In.We note that the experimental The ML for In is from an empirical EoS with experimental parameters [80].The ML for Fe is from Ref. [81] shown in the inset due to a different temperature range and up to 60 GPa where static experiments exist.(b): the same MLs in the low-pressure and low-temperature range (circles).The lines show the fits to Eq. ( 19) as discussed in text.
values in the second line in Table 1 systematically exceed the theoretical linear term α l B l .This is expected since, as discussed earlier, α l B l in Eq. ( 19) underestimate the linear slope in Eq. ( 16).This underestimation brings theoretical and experimental values closer.
The overall agreement indicates that the theory predicts the right order of parameters of experimental MLs and is quantitatively sensible given approximations made.In particular, the theory predicts that the slope of MLs is set by thermal and elastic properties of the system and is on the order of α l B l or about 1-10  Theoretical linear slopes of MLs, α l B l , evaluated using α l and B l calculated from the data in Ref. [80] for In and Ref. [82] for Fe in the middle of the pressure range in Fig. 1b to get characteristic values across the ML.For Ar, He, H2 and H2O, α l and B l are calculated from the NIST data around 1 GPa, the highest pressure at which the data is available [46].Experimental ∆P ∆T are evaluated as average slopes of MLs in Fig. 1b.The last entry shows α l B l from the fits of experimental MLs in Fig. 1b to Eq. (19).opposed to, for example, kPa K or GPa K ).
As a second test of our theory, we fit the experimental MLs in Fig. 1b to the parabolic form of Eq. ( 19) as , where D 1 and D 2 are variable.Polynomial fitting in Origin produces the fits shown in Fig. 1b with the R-square values very close to 1.We note that the parabolic form (19) can be fitted to the entire range of MLs in Fig. 1a and not only to the low pressure and temperature range in Fig. 1b where our approximation was assumed to hold.
The fitted values of D 1 and D 2 give two ways to test the theory.First, Eq. ( 19) predicts D 1 = α l B l .The fitted values D 1 = α l B l are listed in Table 1 where we observe that the ratio of fitted and calculated values of α l B l is 1.3-1.8 for all systems except H 2 O where this ratio is 2.9.Second, Eq. ( 19) predicts D2 D1 = α l r 2(r−1) .Using r = 1.1 as before, we find that the ratio of α l obtained from fitted D2 D1 and α l calculated using the data of Refs.[46,80,82] is 0.4-2.6 for systems in Fig. 1b.
The agreement of fitted of α l B l and α l with their experimental values indicates that, similarly to the first test results, the theory gives sensible quantitative predictions of ML properties overall.This agreement also makes a connection to the discussion in the previous section: parameters of MLs are governed by fundamental physical constants.In this section, we have seen that experimental parameters of MLs are in the characteristic range set by α and B. Indeed, the coefficient in front of the linear term in Eq. ( 19) set by α l B l is experimentally on the order of 10 MPa K in Table 1.Similarly, the experimental coefficient in front of the quadratic term is in the range set by typical α l and B l as predicted.As discussed in the previous section, this characteristic range and the associated universality of ML parameters can be understood from Eq. (21) showing that α and B are governed by fundamental constants.

F. Relation to previous forms of the ML
The quadratic form of the ML in Eq. ( 18) interestingly compares to the commonly used empirical functions used to fit the ML such as the Simon-Glatzel (SG) equation of the form [3,8] where D 1 , D 2 and α are free fitting parameters.The fitted α is about 1.6 for noble Ar and He, 1.8 for H 2 and takes a larger value of 3 for H 2 O [79].Hence this fitted α is often intermediate between the linear and quadratic terms in our theoretical parabolic ML in Eq. ( 18).This is expected in view that both functions fit experimental data well.We also note that (a) differently to the empirical functions used, the theoretical ML derived here does not have free fitting parameters: as follows from Eq. ( 18) or (19), the coefficients in front of linear and quadratic terms are fixed by the system thermal and elastic properties and the constant term is fixed by the location of the triple point and (b) the SG equation gives dP dT > 0 only [3], whereas our theory describes anomalous melting too as discussed in Section 2 C.
A similar form of the ML follows from the soft-sphere model where interactions are purely repulsive with the inverse power-law potential (IPP) U ∝ 1 r n .The IPP model is thought to approximately describe interactions in simple systems such as Ar in a limited pressure range corresponding to the interatomic distances that are sufficiently shorter than the potential minimum where the attractive term is important and sufficiently longer than the short range where the IPP is too steep compared to the potential in real systems [3].For the IPP, the free energy trivially scales with V T 3 n [1], as do other basic properties including the equation of state [75,83].Ref. [84] explores scaling arguments applied to certain melting properties using molecular simulations.
The above scaling of the free energy and the equation of state for the IPP model results in the SG form (22) of the ML with α = 1 + 3 n [3,75]: One might wonder how the parabolic form of the ML derived here relates to the scaling form (23) in the IPP model.We recall that the parabola for the ML comes from Eq. ( 14) which is quite general: it follows from the CC equation (1) using one condition only, the equality of liquid and solid c v (or energy) across the melting line as implied by theory and experiments.On the other hand, the IPP model is a peculiar distinctive system constrained by several rigid scaling conditions.Two relevant constraints are: (a) the entropy change across the ML is constant, ∆S=const (compare this with the condition T ∆S=const used to derive boiling and subli-mation lines, see Introduction), and (b) the volumes of both phases scale as V ∝ T − 3 n along the ML [75].These constraints change the differential equation for the ML.In particular, the right-hand side of Eq. ( 6) becomes zero and derivatives dV dT in the second term simplify and become dV dT ∝ − 3 n T − 3 n −1 .The differential equation ( 6) for the IPP model then becomes and is different from Eq. ( 14) governing the ML.Eq. ( 23) for the ML follows from Eq. ( 24).The same Eq.( 23) follows if the two constraints from the IPP model, ∆S=const and V ∝ T − 3 n , are used in the CC equation (5).
It is therefore an expected result that the derived parabolic ML differs from the scaling form (23): by imposing scaling constraints, the IPP model changes the differential equation for the ML and its solution.
In passing, we note that α = 1 + 3 n in Eq. ( 23) is α = 1.25 if a typical value n = 12 is used [85] and α = 1.17 if n is taken as n = 18 (n = 18 is considered as an effective exponent of the IPP approximating the Lennard-Jones (LJ) potential [86]).These α are significantly smaller than α = 1.6 experimentally seen in Ar and He [79] where the LJ potential is considered to operate and where the IPP is assumed to approximately apply at high pressure.To be consistent with the experimental α = 1.6, n in Eq. ( 23) should be n = 5, which is significantly smaller than what's expected for the IPP in real systems.The origin of this discrepancy is unclear.This is especially so in view that the scaling property V ∝ T − 3 n of the IPP model applies to the ML of Ar well with n close to expected n = 12 [87].A potentially useful observation is that the pressure range in these experiments [79,87] was different.The earlier experiments yielding n = 12 from density scaling in Ar were up to about 1 GPa [87], whereas the P (T ) ML was fitted to Eq. ( 23) with α = 1.6 (n = 5) up to 50 GPa, or over 10,000 times the critical pressure [79].In He, the ML was fitted to Eq. ( 23) with α = 1.6 up to 24 GPa, or over 100,000 times the critical pressure (the ML could not be fitted at higher pressure up to 42 GPa) [79].This difference of pressure ranges could affect the fitted exponents, which would be qualitatively consistent with the variability of the effective n with the state point for the LJ potential [88].This then raises a more general question of the extent to which earlier models including the IPP model applies to real physical systems and MLs in particular.

FINAL REMARKS AND SUMMARY
We make two final remarks.First, recall that the properties of the sublimation and boiling lines are derived by assuming that the latent heat q in the CC equation (1) is constant [1,2].Not assuming constant q and considering a more general case was important in our theory, in Eqs. ( 5), ( 6) and all equations that followed.This suggests that our theory can be applied to those parts of the sublimation and boiling lines where q is not constant and where the textbook solution of the CC equation, P (T ) ∝ e − q T (2), deviates from experiments [5].This would also enable to discuss these two lines in more general terms and improve their understanding.
Second, the derived parabolic line for the melting line can be used to draw the Frenkel line in the supercritical state [65].Indeed, the Frenkel line starts just below the critical point and runs approximately parallel to the melting line in the double-logarithmic plot because both lines correspond to the qualitative change of particle dy-namics [89].
In summary, we proposed a general theory of the ML for the first time.This theory describes the ML in terms of thermal and elastic properties of liquid and solid phases.We showed that the approximate solutions of the theory are in agreement with experimental MLs.The agreement can be refined by fully solving Eq. ( 14) with variable a, b and d in future work.We also showed that the parameters of the ML parabola are governed by the fundamental constants.For this reason, MLs can not be too different in different systems and in this sense possess universality.
I am grateful to V. Brazhkin, J. Dyre, J. Proctor, P. Tello and G. de With for discussions and EPSRC for support.
dV dT to the standard definition of α = 1 B is the isothermal bulk modulus.Similarly, dS dT in Eq. (6) is dS dT = ∂S ∂T P ∂T P − αV dP dT .Using these dV dT and dS dT in Eq. (

FIG. 1 :
FIG.1:(a): Experimental MLs of Ar, He, H2 and H2O[79].The ML for In is from an empirical EoS with experimental parameters[80].The ML for Fe is from Ref.[81] shown in the inset due to a different temperature range and up to 60 GPa where static experiments exist.(b): the same MLs in the low-pressure and low-temperature range (circles).The lines show the fits to Eq. (19) as discussed in text.