Isomorph invariance in the liquid and plastic-crystal phases of asymmetric-dumbbell models

We present a numerical study of the asymmetric dumbbell model consisting of ``molecules'' constructed as two different-sized Lennard-Jones spheres connected by a rigid bond. In terms of the largest (A) particle radius, we report data for the structure and dynamics of the liquid phase for the bond lengths 0.05, 0.1, 0.2, and 0.5, and analogous data for the plastic-crystal phase for the bond lengths 0.05, 0.1, 0.2, and 0.3. Structure is probed by means of the AA, AB, and BB radial distribution functions. Dynamics is probed via the A and B particle mean-square displacement as functions of time and via the rotational time-autocorrelation function. Consistent with the systems' strong virial potential-energy correlations, the structure and dynamics are found to be isomorph invariant to a good approximation in reduced units, while they generally vary considerably along isotherms of the same (20%) density variation. Even the rotational time-autocorrelation function, which due to the constant bond length is not predicted to be isomorph invariant, varies more along isotherms than along isomorphs. Our findings provide the first validation of isomorph-theory predictions for plastic crystals for which isomorph invariance, in fact, is found to apply better than in the liquid phase of asymmetric-dumbbell models.

We present a numerical study of the asymmetric dumbbell model consisting of "molecules" constructed as two different-sized Lennard-Jones spheres connected by a rigid bond.In terms of the largest (A) particle radius, we report data for the structure and dynamics of the liquid phase for the bond lengths 0.05, 0.1, 0.2, and 0.5, and analogous data for the plastic-crystal phase for the bond lengths 0.05, 0.1, 0.2, and 0.3.Structure is probed by means of the AA, AB, and BB radial distribution functions.Dynamics is probed via the A and B particle mean-square displacement as functions of time and via the rotational time-autocorrelation function.Consistent with the systems' strong virial potential-energy correlations, the structure and dynamics are found to be isomorph invariant to a good approximation in reduced units, while they generally vary considerably along isotherms of the same (20%) density variation.Even the rotational time-autocorrelation function, which due to the constant bond length is not predicted to be isomorph invariant, varies more along isotherms than along isomorphs.Our findings provide the first validation of isomorph-theory predictions for plastic crystals for which isomorph invariance, in fact, is found to apply better than in the liquid phase of asymmetric-dumbbell models.

I. INTRODUCTION
This paper presents a systematic study of the asymmetric dumbbell model (ASD) consisting of two different-sized Lennard-Jones (LJ) spheres connected by a rigid bond.The ASD model is designed to be a simple molecular model [1][2][3].Because of the asymmetry, the liquid phase of the model is easily supercooled, i.e., the model does not readily crystallize [2].This makes it suitable for numerical studies as a simple, single-component glass-forming model of non-spherically symmetric constituents [2,[4][5][6][7][8].
The focus here is not on supercooling and glass formation, however, but on investigating what happens when the bond length is varied in the less viscous liquid as well as plastic-crystal phases.It was previously shown that for the bond length 0.58 (in units of the largest (A) particle's radius), the ASD model obeys the hidden-scale invariance symmetry, which manifests itself in strong virial potential-energy correlations in the thermal-equilibrium constantvolume fluctuations [2].This implies the existence of so-called isomorphs, which are lines in the thermodynamic phase diagram along which structure and dynamics in reduced units are invariant to a good approximation [2,9,10].We confirm this below and demonstrate that the existence of isomorphs is robust to bond-length changes.This is done by simulating the ASD model in the liquid phase with bond lengths 0.05, 0.1, 0.2, and 0.5.We also investigate the model's plastic-crystal phase with a focus on the existence of isomorphs; here the largest bond length simulated is 0.3 (for larger bond lengths the system was liquid at the chosen reference state point).
Our investigation reports data illuminating the degree of isomorph invariance of structure and dynamics.This is done by comparing results along isomorphs with results along isotherms of the same density variation.Specifically, the structure is investigated by monitoring the radial distribution function (RDF), while the dynamics is investigated via the mean-square displacement (MSD) as a function of time, as well as the rotational time-autocorrelation function (RAC).We find good isomorph invariance of both structure and dynamics.
The paper presents the first investigation of the consequences of hidden scale invariance for a system of molecules forming plastic crystals.In the plastic-crystal phase the molecules' centers are ordered on a crystalline lattice-in the present case a face-centered cubic lattice-while the molecular orientations vary more or less randomly because the molecules are free to rotate [11][12][13][14][15].In this phase we find excellent isomorph invariance, in fact better than in the (isotropic) liquid phase.

II. MODEL AND SIMULATION DETAILS
The asymmetrical dumbbell (ASD) is a constrained molecular model.It consists of two different-sized Lennard Jones (LJ) particles, a sphere (A) and a slightly smaller one with lower interaction energies (B).If the two spheres are connected by a rigid bond of length = 0.584, the model mimics toluene [2,16].The values of the parameters used below are as follows.For particle A the distance and energy parameters are σ AA = ε AA = 1 and the mass is m A = 1; for particle B one has σ BB = 0.788, ε BB = 0.117, and m B = 0.195; for the AB interaction one has σ AB = 0.894 and ε AB = 0.342 [16].
We define the (number) density ρ as the total number of A and B particles (N = 8000) divided by the simulation box volume V ; T is below the temperature.The system is studied by molecular-dynamics simulations in the canonical N V T ensemble using the Nosé-Hoover thermostat [17].The simulated system consists of 4000 ASD molecules in a cubic box with periodic boundaries.Simulations are performed using the open-source Roskilde University molecular dynamics software (RUMD) that runs on GPUs (graphics processing units) [18].The constant bond-length constraint was implemented using the method of Ref. [19].
The leapfrog algorithm is used with time step dt = 0.001 (in A particle LJ units).First, N V E simulations were carried out for 10 6 time steps for equilibration and 10 6 time steps for production runs with time step sizes 0.0005, 0.005, 0.001, 0.002, and 0.01 in order to investigate how constant the energy is at the reference state point defined by ρ = 1.5 and T = 1.5.The bond lengths chosen for these N V E simulations were 0.05 and 0.5.We concluded that the time step 0.001 in both cases shows a good energy conservation.For the subsequent N V T simulations generating the data, at each state point the simulations again ran for 10 6 time steps for equilibration and 10 6 time steps for the production runs.
For the post analysis of structure and dynamics, the units used are the isomorph-theory's reduced units defined as follows: energy unit: k B T , length unit: ρ −1/3 , time unit: ρ −1/3 m A /k B T .Note that these units vary with state point.Exceptions to the use of reduced units apply for density and temperature, which by definition are both unity in reduced units; these quantities, as well as bond lengths, are reported in A particle LJ units.

III. ESSENTIAL ISOMORPH THEORY
Hidden scale invariance is a symmetry of the potential energy function U (R) in which the configuration R in terms of the particle coordinates is defined by R ≡ (r 1 , . . ., r N ).If R a and R b are two same-density configurations, hidden scale invariance is defined by the following mathematical implication (in which λ is a uniform-scaling parameter) This symmetry is only rigorously obeyed if U (R) is an Euler-homogeneous function plus a constant, but it applies to a quite good approximation for many other systems, e.g., single-and multicomponent LJ systems, the Yukawa and bi-Yukawa pair potentials, the EXP and Morse pair potentials, Mie potentials involving various exponents, etc. [2,9,10,[20][21][22]. Equation (1) may also apply for molecular models like the ASD model.The derivations given below apply in the case of constant-length rigid bonds, i.e., not scaled with λ.Equation ( 1) then refers to a uniform scaling of the center-of-mass coordinates with molecular orientations and sizes kept unchanged.The excess entropy S ex is defined as the entropy minus that of an ideal gas at the same density and temperature.Any state point of the thermodynamic phase diagram is fully characterized by the two thermodynamic variables ρ and T , but other variables may of course also be used to characterize a state point, for instance density and average potential energy.We define the microscopic excess-entropy function S ex (R) by [20] S ex (R) ≡ S ex (ρ, U (R)) . ( Here the function S ex (ρ, U ) is the excess entropy of the state point (ρ, U ) in which U is the average potential energy.It follows from statistical mechanics that S ex (R) is proportional to the logarithm of the number of configurations with the same density and potential energy as R [23].Inverting Equation (2) leads to in which U (ρ, S ex ) is the average potential energy of the state point of density ρ and excess entropy S ex .
It can be shown that Equation (1) implies the function S ex (R) is scale invariant, i.e., S ex (λR) = S ex (R) [20].In this case, S ex (R) depends only on the configuration's so-called reduced coordinate vector R ≡ ρ 1/3 R (which is of course invariant upon a uniform scaling): Consequently, Equation (3) becomes All identities of the isomorph theory may be derived from Equation (5) [20].In particular, Equation (5) implies strong correlations between the constant-volume fluctuations of the virial W and the potential energy U .Thus, the equilibrium deviations from the average of these two variables, ∆W and ∆U , obey [24][25][26] ∆W ∼ = γ ∆U .
Here the so-called density-scaling exponent γ is defined and characterized by [9] γ The second (generally valid) identity allows one to calculate γ from the constant-volume thermal-equilibrium canonicalensemble fluctuations at the state point in question.
Using the identity T = (∂U/∂S ex ) ρ , a first-order Taylor expansion of Equation ( 5) leads [20] to Consider two state points (ρ 1 , T 1 ) and (ρ 2 , T 2 ) with average potential energies U 1 and U 2 and the same excess entropy S ex , and suppose that R 1 and R 2 are equilibrium configurations of the two state points with the same reduced coordinates, i.e., ρ In that case, by elimination of S ex ( R) − S ex , Equation (8) implies that with T 1 ≡ T (ρ 1 , S ex ) and T 2 ≡ T (ρ 2 , S ex ) one has This implies that in which C 12 is a constant.Equation ( 10) is the original 2009 definition of isomorphic state points [9], which stated that the canonical probabilities of configurations that scale uniformly into one another are identical along an isomorph (the value of the constant C 12 is irrelevant because probabilities are normalized).This identity of the probabilities of scaled configurations implies that S ex is constant along an isomorph.In fact, the 2014 version of the isomorph theory that introduced Equation (1) defined isomorphs as lines of constant excess entropy in the thermodynamic phase diagram and showed from this that the excess entropy is an isomorph invariant [20].
It can be shown, either from Equation (10) [9] or from Equation ( 5) [20], that the dynamics at two isomorphic state points are identical in the "same movie" sense: Filming the molecules' motion at one state point results in the same movie at a different, isomorphic state point-except for a uniform scaling of space and time.This implies several dynamic isomorph invariants and that the reduced-unit structure is isomorph invariant in reduced units.
How does one decide whether a given system is expected to have good isomorphs?According to Equation (6) this is the case when the virial and potential-energy fluctuations are highly correlated, which can be investigated by evaluating their Pearson correlation coefficient, A system is "R-simple", i.e., simple in the sense of having good isomorphs and thus an essentially one-dimensional thermodynamic phase diagram, if R > 0.9 at the state points of interest [20,24,25].This is a somewhat arbitrary criterion, though, in the sense that the property of good isomorph invariance depends on the quantity in question-in practice some reduced-unit quantities are more isomorph invariant than others.Note also that how invariant a given property is, of course, depends on how large a density range that is explored.

IV. RESULTS FOR THE LIQUID PHASE
In this section, we investigate the liquid phase for the bond lengths 0.05, 0.1, 0.2, and 0.5.First, however, Figure 1a shows a snapshot of the ASD liquid at the bond length previously studied by Ingebrigtsen et al. [2] (0.58) [2] at a liquid state point.The A particles are red, the B particles are blue.We see the typical disorder of a liquid configuration.Reducing the bond length results in a system that is increasingly like the atomic LJ liquid, but even for the smallest bond length (0.05) we still find typical molecular-liquid behavior (see below).Isomorphs, i.e., curves of constant excess entropy S ex , are traced out by numerically integrating Equation (7).In this way one avoids thermodynamic integration to determine S ex throughout the phase diagram.Specifically, the last term of Equation ( 7) specifies how the density-scaling exponent γ may be calculated from an N V T simulation at the state point in question.Integrating the first-order differential equation defined by the second equality sign in Equation ( 7) in order to determine how temperature varies with density along an isomorph, is in principle straightforward.The highly accurate fourth-order Runge-Kutta (RK4) integration method was recently implemented and used for this [27,28].We used that method here with a density change of 1%.For all four bond lengths we started the integration at the reference state point (ρ, T ) = (1.5, 1.5).Thus, all four isomorphs have this state point in common.Table I gives the pressure, density-scaling exponent, and virial potential-energy correlation coefficient for the four bond lengths at the reference state point.TABLE I. Thermodynamic parameters at the liquid reference state point, giving density, temperature, and pressure (in A particle LJ units), as well as the density-scaling exponent γ and the virial potential-energy correlation coefficient R. The isomorphs are not very different (Figure 1b), but we note that larger bond lengths result in slightly larger temperatures at a given density.This reflects larger density-scaling exponents γ for the larger bond lengths, as is evident from Figure 1c which shows γ at the reference state point for several bond lengths.A value of γ between 5 and 6 is typical for atomic LJ systems [25,26].In fact, in the liquid phase a density-scaling exponent between 5 and 6 is not only typical of atomic LJ system, but also for diatomic LJ systems that often have values somewhat higher than the atomic case [29,30].Interestingly, density-scaling exponent relatively close to 6 are often found in experiments [31][32][33].We note a slight drop of γ at the largest bond length studied, but the overall γ variation is just 15% while, as we shall see, the physics varies considerably when the bond length is changed from 0.05 to 0.5.
As mentioned, the isomorph theory only applies for systems with strong virial potential-energy correlations.To ensure that this requirement is obeyed, we evaluated the correlation coefficient R (Equation ( 11)) along the isomorphs and along T = 1.5 isotherms.The results are shown in Figure 2 as functions of the density in which the traditional R = 0.9 pragmatic limit for isomorph theory to apply is marked by horizontal dashed lines.We see that the majority of state points are above this line, suggesting that the isomorph theory applies.The strongest correlations are found for the largest bond length studied (0.5).
To investigate to which degree there is isomorph invariance of the reduced-unit structure we performed simulations of the three radial distribution functions (RDF) g(r) defined by the AA, AB, and BB particle pairs [17] along the isomorphs.The results are shown in Figure 3, which for comparison also shows data for simulations carried out at the reference-state point isotherm with the same (20%) density variation, i.e., for densities between that of the reference state point, 1.5, and 1.8.All three RDFs are nicely invariant along the isomorphs compared to their variation along the isotherms.We note a slight variance at the first peak maximum.This deviation from isomorph theory has often been observed; it is due to the fact that in the direction of increasing γ along an isomorph, the system becomes more hard-sphere like.Consequently, the left hand side of the first peak diminishes when moving along an isomorph in the direction of increasing γ.To "compensate" for this effect and maintain a constant overall coordination number (defined as the RDF integral over the first peak), the peak maximum increases.The variance of the first peak maximum has recently been put into a consistent framework by proofs that the so-called bridge function is isomorph invariant [28].Focusing henceforth on the isomorph RDFs, we note that the AA RDF does not change very much as the bond length is increased.In fact, this RDF looks much like that of the standard single-component LJ liquid, which is also very similar to that of, e.g., the hard-sphere, Yukawa, and inverse-power-law liquids [34].The observed variation of the AB and BB RDFs can be interpreted as a consequence of the near invariance of the AA particle RDF as the bond length is increased: The structure of the ASD liquid may be thought of as primarily determined by the A particles, which behave much like standard LJ particles that are not really bothered by the B particles.This is because the B particles are smaller and have significantly lower interaction energy parameters.In this picture, the B particles are to a significant degree "slaves" of the A particles, and as a consequence one expects what is observed in Figure 3: The BB RDF changes significantly with increasing bond length, and at bond length 0.5 it is almost constant.This is because the B particles place themselves in many possible positions around the A particles that, as mentioned, more or less have the RDF of a single-component LJ liquid.This effect is present at all bond lengths, but for the smaller bond lengths it does not give rise to any significant spread of the BB RDF because a given B particle is constrained to be close to the A particle of the same molecule.Indeed, the BB particle RDF is very similar to the AA RDF for the small bond length 0.05.This "B slaving A" picture is confirmed by the AB RDF, which (except for the vertical line coming from the intramolecular AB bond) gives data for the AB correlations between different ASD molecules.The AB RDF is also significantly smeared out as the bond length increases, but less so than the BB RDF because the relative order of the A particles is partly inherited by the AB RDF.
Figure 4 shows data for the (reduced) time dependence of the reduced mean-square displacement (MSD) [17] of the A and B particles in a log-log plot.The data are isomorph invariant to a good approximation, but not invariant along the isotherm except in the short-time ballistic regime where the invariance is a rigorous consequence of the use of reduced units, which implies unity thermal velocity.Focusing on the behavior along the isomorphs, we note that the A particle motion is pretty similar for all bond lengths.This is consistent with the above picture according to which the A particles to a significant extent behave as if the B particles were not present.In the ballistic regime, the B particles move faster than the A particles because of their lower mass.At long times, the A and B particles are constrained to follow each other, resulting in the same long-time MSD for all bond lengths (which, of course, also applies along the isotherms).An interesting feature appears at intermediate times for the B particle MSD, which for small bond lengths have a slight kink that at larger bond lengths develops into a hint of a plateau.We interpret this as an effect of the fact that the fast ballistic motion of the B particles around the A particles eventually "saturates".This is confirmed in Figure 5 discussed next, showing that the rotational time-autocorrelation function for short bond lengths has a negative minimum in this range of the reduced times (around 0.1).(a) The RDFs along the isomorph for the bond length 0.05; for comparison the same RDFs are shown along the referencestate-point isotherm of the same (20%) density variation.We see good, but not perfect invariance along the isomorph, unlike along the isotherm.The thick vertical line in the AB RDF comes from the fixed bond length, which in reduced units varies with density.(b-d) show similar plots for bond lengths 0.1, 0.2, and 0.5, respectively.There is good isomorph invariance of all three RDFs in comparison to their isotherm variation.We note that the first peak of the BB RDF gets lower as the bond length increases.This reflects an increased spread of the B particle positions relative to each other.At the same time, the AB RDFs are lowered and their second peak almost disappears.Figure 5 shows data for the rotational time-autocorrelation function (RAC) defined as the autocorrelation of the unit vector from particle A to particle B, n, plotted as a function of the reduced time.Because of the normalization, this quantity always starts in unity at time zero.The upper row shows data for the isotherms, the lower for the isomorphs.There is little difference and good collapse in both cases, except for the largest bond length 0.5 for which the data are significantly less invariant along the isotherm.At this bond length there is also a change of behavior in the sense that the RAC never becomes negative as it does at intermediate times for the lower bond lengths.A negative value of the RAC signals a more than 90 degree average rotation of the molecule.This finding for bond length 0.5 is consistent with the above observation that a change of the physics appears to set in at bond lengths 0.3-0.4.It is important to note that, in contrast to the RDF and the MSD, the RAC is not predicted to be isomorph invariant because the moment of inertia is not isomorph invariant (the bond length is fixed and not scaled with the density by the factor ρ −1/3 required for a constant reduced-unit moment of inertia).

V. RESULTS FOR THE PLASTIC-CRYSTAL PHASE
We proceed to report results for the plastic-crystal phase for the bond lengths 0.05, 0.1, 0.2, and 0.3 (for bond lengths larger than 0.3 the system was not crystalline at the chosen reference state point of density 2.2 and temperature 0.5).Table II gives the pressure, density-scaling exponent, and virial potential-energy correlation coefficient for the four bond lengths at the reference state point.A typical configuration is shown in Figure 6.We clearly see that the system is ordered (compare Figure 1a).A closer inspection reveals that the order is not perfect, however, because the bond directions are disordered.This is a typical signal of a plastic crystal.At a lower temperature, there must be a phase transition to a perfectly ordered phase in which the bond orientations are also ordered, but we did not investigate that transition.TABLE II.Thermodynamic parameters at the plastic-crystal reference state point, giving density, temperature, and pressure, as well as the density-scaling exponent γ and the virial potential-energy correlation coefficient R.  Before proceeding to probe the extent of isomorph invariance, we look at the properties at the reference state point.Figure 7a shows the AA particle RDFs for the four bond lengths.These show a well-defined order due to the (almost) crystalline ordering of the A particles; the AA particle RDFs are moreover very similar for all bond lengths.In contrast, the AB and BB RDFs vary significantly because of the rotation of the B particles around the A particles.In these cases, the most ordered RDFs are those of the shortest bonds.These observations show that the slaving of the B particles to the A particles (that are well ordered) applies also in the crystalline phase.Figure 8 shows the MSD of the A and B particles as functions of the reduced time.For the A particles the results are very similar for the different bond lengths, although bond length 0.3 deviates from the three smaller ones by having a larger long-time plateau.We have no good explanation for this, but deviations for the longest bond length from the three others are also noted in some of the later figures.The B particle plateaus vary considerably, with larger plateaus observed for larger bond lengths.This is because the B particles rotate around the A particles and a larger bond length gives them more freedom to do so, resulting in a larger long-time plateau.Figure 9 shows γ and R at the reference state point as a function of the bond length.We see that γ is almost independent of the bond length and has the typical LJ value between 5 and 6 [25].The correlation coefficient drops somewhat with increasing bond length, but stays comfortably above 0.9.As in the liquid phase, we generated isomorphs from the reference state point by integrating Equation ( 7) using the RK4 algorithm [27] for density steps of magnitude 0.01.Additionally, as in the liquid case, we chose to increase the density by 20% from the reference state-point density.This corresponds to realistic density changes of high-pressure experiments.The four isomorphs are shown in Figure 10.We see that they are quite similar.Figure 11 shows the reduced RDFs along the isomorphs and, for comparison, along reference-temperature isotherms of the same density variation.For the smallest bond length (0.05), we see a very good collapse of all three RDFs along the isomorph, but not along the isotherm.Moreover, the three RDFs are quite similar, which is a consequence of the fact that all B particles are constrained to be very close to an A particle.For larger bond lengths we see in all cases a very good isomorph collapse and note that, with increasing bond length, the AB and BB particle RDFs become increasingly smeared out compared to the AA particle RDF.The same was seen in the liquid phase (Figure 3), and as for the liquid phase we interpret this as deriving from rotations of the molecules.In particular, it confirms that the crystals are plastic.Interestingly, the AA particle RDF becomes more invariant along the isotherm as the bond length is increased, but for all bond lengths this quantity is less isotherm than isomorph invariant.For the MSD as a function of reduced time (Figure 12), we note that these are all isomorph invariant to a good approximation.In contrast, there is a notable variation along the reference-state-point isotherm for both A and B particles.Only at short times (in the ballistic regime) do we observe invariance along the isotherms, but as mentioned this is a consequence of the definition of reduced units.

VI. SUMMARY
We have presented a numerical study of the liquid and plastic-crystal phases of asymmetric dumbbell models of different bond lengths with the purpose of testing the isomorph-theory predictions of invariant structure and dynamics.In the liquid case, the bond lengths 0.05, 0.1, 0.2, and 0.5 were studied, while in the crystalline case the largest bond length was 0.3.At all state points the virial potential-energy correlation coefficients were above 0.84; in the vast majority of cases it was above 0.9, and for the plastic crystals it was quite close to unity.This implies that the isomorph theory is expected to apply [9,10].Indeed, we found good isomorph invariance of the reduced-unit RDF and MSD in the liquid phase and even better invariance in the plastic-crystal phase.Isomorph invariance also applies to a good approximation for the RAC, despite the fact that this quantity is not predicted to be an exact invariant because the bond length is constant in real units, not in reduced units.In contrast, these three quantities were generally not invariant along isotherms of the same (20%) density variation.The behavior of the B particles conform to a picture in which these are "slaves" of the A particles, which behave almost like single-component LJ particles, i.e., are not much affected by the B particles.This paper is the first time the isomorph theory has been shown to apply for plastic crystals of simple two-atom dumbbell molecules.In conjunction with recent papers demonstrating the applicability of isomorph theory to liquid crystals [35,36], the present findings demonstrate the generality of this theoretical framework derived from the hidden-scale-invariance condition Equation (1) and, in particular, that the isomorph theory applies independent of which phase the system is in.It would be interesting to investigate the transition from the plastic crystal to the ordinary, orientationally ordered crystal phase and to study whether the phase transition line is close to an isomorph, as has been found for the liquid-solid transition line [37].

FIG. 1 .
FIG.1.(a) Snapshot of a liquid configuration at ρ = 0.8, T = 0.8.The bond length between particle A (red) and B (blue) is in this case 0.58 (in A particle units)[2].(b) Four isomorphs traced out in the liquid regime of the ASD thermodynamic phase diagram using the fourth-order Runge-Kutta (RK4) method with density step size 0.01.Each isomorph is traced out with a different bond length (from 0.05 to 0.5), starting from the reference state point marked in red.The isomorphs are qualitatively similar, but we observe a slight increase in the temperatures of the traced state points with increasing bond length.(c) The variation of the density-scaling exponent γ (Equation (7)) at the reference state point ρ = 1.5, T = 1.5 of each isomorph for the different bond lengths.Data for a few extra bond lengths have been added in order to have a better view of the γ variation.We see that γ has a maximum around bond length 0.3-0.4,indicating a change of physics here.
FIG. 2. (a) Variation of the virial potential-energy correlation coefficient R along each isomorph plotted for the four bond lengths.(b) The variation of R plotted as a function of density along along the T = 1.5 isotherm for the four bond lengths.

FIG. 3 .
FIG.3.AA, AB and BB radial distribution functions (RDF) at the reduced pair distance ρ 1/3 r where r is the pair distance.(a) The RDFs along the isomorph for the bond length 0.05; for comparison the same RDFs are shown along the referencestate-point isotherm of the same (20%) density variation.We see good, but not perfect invariance along the isomorph, unlike along the isotherm.The thick vertical line in the AB RDF comes from the fixed bond length, which in reduced units varies with density.(b-d) show similar plots for bond lengths 0.1, 0.2, and 0.5, respectively.There is good isomorph invariance of all three RDFs in comparison to their isotherm variation.We note that the first peak of the BB RDF gets lower as the bond length increases.This reflects an increased spread of the B particle positions relative to each other.At the same time, the AB RDFs are lowered and their second peak almost disappears.

FIG. 4 .
FIG. 4. (a)The A and B reduced-unit mean-square displacement (MSD) along the isomorph and the isotherm for bond length 0.05.There is invariance along the isomorph for both the A and B reduced MSD, but not along the corresponding isotherm.(b-d) show similar plots for bond lengths 0.1, 0.2, and 0.5.There is invariance along the isomorphs for both MSDs, but not along the corresponding isochores.

FIG. 5 .
FIG. 5. a-d The rotational time-autocorrelation function (RAC) along the isomorphs and the isotherm for the four different bond lengths.As the bond length increases, the decay to zero becomes significantly slower.The largest bond length (0.5) behaves differently from the others by not going below zero, confirming the change of physics suggested by Fig. 1c.

FIG. 6 .
FIG.6.Snapshot of a plastic crystal configuration at the reference state point (ρ, T ) = (2.2,0.5).The bond length between particles A (red) and B (blue) is in this case 0.2.

3 FIG. 8 .
FIG. 8. (a) A particle MSD as a functions of time at the reference state point (ρ, T ) = (2.2,0.5) for each of the bond lengths 0.05, 0.1, 0.2, 0.3.The long-time stabilization to a plateau shows that the system is a solid.(b) B particle MSD at the same state point.

5 FIG. 9 .
FIG. 9. Variation of the density-scaling exponent γ (a) and of the virial potential-energy correlation coefficient R (b) at the reference state point plotted as functions of the bond length.Data for a few extra bond lengths have been added.

3 FIG. 10 .
FIG.10.Four isomorphs traced out in the plastic-crystal phase of the ASD thermodynamic phase diagram using the RK4 method with density step size 0.01.The isomorphs are traced for bond lengths 0.05, 0.1, 0.2, and 0.3, starting from the reference state point marked in red.

FIG. 11 .
FIG.11.AA, AB and BB RDFs as a function of the reduced pair distance.(a) The RDFs along the isomorph for the bond length 0.05 and, for comparison, along the reference-state-point isotherm of the same (20%) density variation.We see almost perfect invariance along the isomorph, unlike along the isotherm.The thick vertical line in the AB RDF comes from the fixed bond length, which in reduced units varies with density.(b-d) show similar plots for bond lengths 0.1, 0.2, and 0.3.There is a very good isomorph invariance of all three RDFs in comparison to their isotherm variation.The first peak of the BB RDF gets lower as the bond length increases, which reflects an increased spread of the B particle positions relative to each other.At the same time the AB RDFs also decrease.

FIG. 13 .
FIG. 13. (a-d)The rotational time-autocorrelation function (RAC) along the isotherm and the isomorphs for the four different bond lengths.As the bond length increases, the decay to zero becomes significantly slower.The largest bond length (0.3) behaves differently from the others by not going below zero.