To thread or not to thread? Effective potentials and threading interactions between asymmetric ring polymers

Abstract

We use computer simulations to study a system of two unlinked ring polymers, whose length and bending stiffness are systematically varied. We derive the effective potentials between the rings, calculate the areas of minimal surfaces of the same, and characterize the threading between them. When the two rings are of the same kind, threading of a one ring through the surface of the other is immanent for small ring–ring separations. Flexible rings pierce the surface of the other ring several times but only shallowly, as compared to the stiff rings which pierce less frequently but deeply. Typically, the ring that is being threaded swells and flattens up into an oblate-like conformation, while the ring that is threading the other takes a shape of an elongated prolate. The roles of the threader and the threaded ring are being dynamically exchanged. If, on the other hand, the rings are of different kinds, the symmetry is broken and the rings tend to take up roles of the threader and the threaded ring with unequal probabilities. We propose a method how to predict these probabilities based on the parameters of the individual rings. Ultimately, our work captures the interactions between ring polymers in a coarse-grained fashion, opening the way to large-scale modelling of materials such as kinetoplasts, catenanes or topological brushes.

---

1 Introduction

Since the beginnings of macromolecular science the possibility of circular polymer topology has captivated scientists attention.^1–3 While in the early days the organic rings were indicated to be doubly-folded into linear structure,⁴ later with the improvement of understanding of molecular interactions, it was recognized the rings can and do attain expanded open conformations. Such conformations have many consequences unparalleled by linear chains. For example, in the long chain limit, the conformations are characterized by the Flory exponent ν ≃ 0.588 (e.g. gyration radius R_g scaling with the polymer length as N^ν) for fixed topology even with no excluded volume,^5,6 thus forming large openings. Yet even at moderate chain lengths the ring opening impacts chain's properties, particularly in the presence of other molecules, because these can thread through the ring.

The threading occurs in diverse settings and concentration regimes. For example, in dilute regime, the threading can support topology-based supramolecular architectures, such as polycatenanes^7–9 or polyrotaxanes, with their movable threaded ring constituents, promising the development of complex synthetic molecular-sized machines.^10–12 In the biological realm, threading occurs across domains of life, typically at higher densities, from semidilute to melt conditions. The kinetoplast is the mitochondrial DNA of trypanosomes in the form of thousands of rings, concatenated, or not, depending on the species^13–19 in semi-dilute conditions. Bacterial plasmids and chromosomes are (often) circular, subject to threading constraints under a strong confinement in case of multiple chromosomes or after genome replication where the subsequent chromosome segregation is assisted by the entropic repulsion²⁰ enhanced by the circular topology.²¹ In eukaryotes the circular extrachromosomal DNA is abundant²² and the linear chromosomes contain many biologically functional domains in the form of loops that reduce the inter-domain interactions.²³ The formation of the loops and their effective repulsion arising from the nonconcatenation and threading constraints, compacts the chromosomes, segregates and disentangles the sister chromatids²⁴ in mitosis and meiosis. The average large scale conformation of long chromosomes is characterized by the same exponents and territorial arrangement as that of long nonconcatenated and unknotted rings in melt^25,26 which in both systems arise from topological (uncrossability) constraints.

Besides their vast biological implications, threadings play a decisive role for determining mechanical properties of materials composed of polymers with ring topology, particularly above the overlap concentration. While the theories on the equilibrium static properties of the ring melt, neglecting ring threadings, are accurate, their generalizations on dynamics are falling short.^27–32 It has been conjectured^33–36 that the threadings become crucial for long ring lengths (correlate with slower dynamics already at shorter lengths^37,38) and might allow for the formation of highly threaded conformations restricting relative ring motion to the extent of topological glass. Some models of the topological glass are based on an assumption that the dynamics of a ring is impacted differently by the nature of the threading: if a ring is actively threading another one, in comparison to the one that is passively threaded.^39,40 As the topological glass has not yet been observed, the role of the nature of the threading has not yet been resolved in this case. However, for example mixtures of rings and linear chains^41–45 exhibit viscosity exceeding that of the components, due to rings passively threaded by the linear chains. The passive ring threading by linear chains is even more evident in nonlinear rheology, in shear or extensional flows.^43,46–48 In other systems, the relevance of (active or passive) threadings, or lack thereof, is manifested too. Tadpole-shaped polymers exhibit a strong slow-down as a function of the tail length due to head-tail threadings. Polymer rings diffuse faster when supercoiled due to reduction of threadings and reduction in entanglements.⁴⁹ The significance of threadings is highlighted even more for systems out of equilibrium. Rings in extensional flow form supramolecular daisy-chains by threadings dramatically increasing the extensional viscosity.^42,50,51 The threading is also enhanced in a system of rings with driven segments leading to the active topological glass.^52,53

While most of the above examples occur at semi-dilute to melt conditions, therefore involving many-body effects, already the simplest case of two threading rings in dilute conditions is not fully understood or described with an effective theory. Here we aim at making this step. Threadings are related to the collective conformational entropy, which is inherently tied to the effective interaction potential. Although the results from the dilute regime cannot be directly transferred to the interesting melt conditions, exactly because of the many-body effects causing strong distortion of the conformations, it has been proved that the effective potential description is accurate to about five times the overlap concentration.^54,55 Therefore understanding the effective potentials and threadings in the simplest two-ring conditions fills at least the bottom part of the interesting concentration range. Moreover, the rings in the melt conditions are largely distorted only on the scales above the entanglement length, while below the ring segments form loops which mutually thread exactly on these scales.^27,56 Therefore, our understanding of threading in dilute conditions is likely relevant also for the loopy ring segments in the melt below the entanglement scale.

Besides the specific interest in threadings, because of to the ubiquitous presence of interacting rings in nature and materials,^57–59 ranging in complexity,^60,61 a focused effort has already been put to elucidate the effective interactions in the simplest case of two unlinked flexible rings of equal^20,54,62,63 and unequal lengths.⁶⁴ The effective potential is repulsive not only due to excluded volume effects, but with a significant topological contribution⁶⁵ due to the non-concatenation constraint. To accommodate these constraints at small separations, one of the two rings swells and the potential has a plateau or a slight attraction depending on the ring topology. If the rings are of equal length, their swelling roles exchanges over time, while for unequal lengths⁶⁴ or topologies,⁶³ dominantly the longer or simpler ring swells respectively. The swelling and the slight attraction are caused by the entropy gain of the encompassed ring on the expense of a relatively smaller swelling-related entropy loss of the larger ring. The swelling of the ring is clearly connected to its (passive) threading that, as mentioned above, can impact the dynamics. The threading was thoroughly investigated,⁶⁴ fitting the picture that the longer ring was passively threaded by the shorter one, but only for a fully flexible polymer model.

The free energy balance, however, must depend not only on the relative sizes, but also on the stiffness of the two rings, in particular, for non-asymptotic ring lengths. The entropy loss of the swelling ring might be compensated, to some extent, by the decrease of the bending energy. In such a case, would it be possible to induce the exchange of the swelling roles, i.e. a shorter ring would swell to encompass a longer one? And related to that, is it possible to change an actively threading ring into a threaded one? The answer is not straightforward as with increasing stiffness, the rings adopt a more planar shape, leading to anisotropic potential and stiffness-induced mutual parallel orientation^66–68 decreasing the threading. At higher densities the orientation can form cluster glasses⁶⁹ due to columnar phases occurring both, in monodisperse and polydisperse solutions. In the latter, the smaller rings are found inside the columns of larger ones,⁶⁸ but mixtures with different stiffnesses have not been investigated. Varying the stiffness also affects properties of rings in films,^70,71 melts⁵⁶ or catenanes.⁷² To understand such complex, many-body effects, we have to at first understand the simpler two-body problem.

As illustrated above, the swelling, threading, orientation and effective potential of rings are a nontrivial function of the relative ring length and stiffness. In realistic setting, be it biological or synthetic, the interacting and threading loops have frequently different (non-asymptotic) sizes and stiffnesses simultaneously, e.g. the chromatin loops vary in sizes, their stiffness depends on the sequence and the histone association extent, simple experimental systems are polydisperse etc. Therefore, it is very important to understand how does the combination of length and stiffness affect the conformations and threadings, and particularly, to generalize these findings for broad applicability in the various circumstances. Here we explore the effective interactions and concentrate on the threading properties of two unknotted rings when their length and stiffness is varied. We will focus on relatively short rings (≤100 beads) with stiffness ranging from flexible polymers like polyethyleneglycol (with characteristic ratio C_∞ ∼ 6), through semi-flexible like polystyrene (C_∞ ∼ 9) up to rather stiff like poly(dodecyl methacrylate) (C_∞ ∼ 14).

After presenting the Model and the methods in the next section, we split the Results into subsections focusing on the specific aspects of the ring–ring interaction: (3.1) the effective potentials; (3.2) the minimal surfaces and ring shapes; (3.3) probability and the role of threading and lastly, (3.4) where we present the generalization of our results in our model of the threading roles.

2 Model & methods

2.1 Microscopic model

Our system contains two unknotted non-concatenated ring polymers, which are modeled by means of a standard bead-spring coarse-grained polymer model.⁷³ A ring is composed of N ∈ {25, 50, 100} point particles, connected by finitely extensible non-linear elastic (FENE) bonds governed by the potentialwhere x is the instantaneous inter-particle distance, K_FENE = 30k_BT/σ² and R₀ = 1.5σ, where we set the unit length as σ = 1 and it represents the range of the steric inter-monomer interaction. Every pair of monomers interact via purely repulsive non-bonded Weeks–Chandler–Andersen potential defined aswhere ε = k_BT and H(·) is the Heavside step-function. The combined potential well of the two above potentials gives rise to uncrossable bonds of mean lengths 〈b〉 ≈ 0.96σ. Finally, we emulate the intramolecular stiffness with harmonic cosine potential

U_bend(ϕ) = K_bend(1 − cos(ϕ))², (3)

where K_bend ∈ {0, 10, 20, 30}k_BT is the bending spring constant and ϕ is the instantaneous angle between consecutive bond vectors. In Table 1, we relate the bending spring constants to usual descriptors of macromolecular stiffness – Flory's characteristic ratio, C_∞, and the persistence length, l_per.⁷⁴ We estimated C_∞ as the long-chain limit of characteristic ratio for a corresponding linear ideal chain, and subsequently using the C_∞ value to calculate l_per of a corresponding freely rotating chain analytically, as we show in Fig. S1 in the ESI.† For a detailed description of this procedure, we refer the reader to the Section IA in the ESI† and to previous studies.^67,69

Table 1 Mapping between the input force constant, K_bend, in bending potential eqn (3) and emerging descriptors of polymer flexibility. See Section IA for more details

βK bend	C ∞	l per/b
0.0	1.0	0.0
10.0	9.4	4.8
20.0	12.9	6.6
30.0	15.4	7.9

---

Our two rings are placed in a cubic simulation cell of size L = 150σ ≳ 10R_g, where R_g is the radius of gyration of the larger of the two rings, hence assuring that a ring does not interact with itself through the periodic boundary condition. These conditions correspond to extremely dilute solution, where the (good) solvent is treated only as uniform continuum.

Ultimately, our problem is akin to exploration of a 4-dimensional space, where we independently vary the length, N_A,N_B and stiffness, C_∞(A),C_∞(B) of both rings. We will focus on three specific slices through the parameter space:

(1) fully symmetric case: rings of equal lengths N_A = N_B ∈ {25, 50, 100} and equal flexibilities C_∞(A) = C_∞(B) ∈ {1.0, 9.4, 12.9, 15.4}

(2) asymmetric case: rings of equal lengths N_A = N_B = 100, one flexible C_∞(A) = 1.0 and other rigid C_∞(B) ∈ {9.4, 12.9, 15.4}

(3) fully asymmetric case: long flexible ring N_A = 100, C_∞(A) = 1.0 and a shorter ring N_B ∈ {25, 50} of arbitrary flexibility C_∞(B) ∈ {1.0, 9.4, 12.9, 15.4}

Finally, in Fig. 1 we show a representative snapshots of the system.

Fig. 1 Representative snapshots of the system at ring–ring separation r ≈ 0: (a) fully symmetric: N_A = N_B = 100,C_∞(A) = C_∞(B) = 15.4, (b) asymmetric: N_A = N_B = 100,C_∞(A) = 1.0,C_∞(B) = 15.4, (c) fully asymmetric: N_A = 100,N_B = 25,C_∞(A) = 1.0,C_∞(B) = 15.4. Ring A is always in blue, while ring B in red.

2.2 Shape parameters

To characterize the shape and size of ring polymers, we opt for the collective variables used in our previous studies.^63,75,76 For an instantaneous conformation of each of the rings, let us consider the gyration tensor, Ŝ, with matrix elements:where α,β ∈ {x, y, z} are the Cartesian coordinates of particles indexed by i and j. Diagonalization of the gyration tensor yields three eigenvectors corresponding to the principal axes of the instantaneous conformation of the macromolecule, and the associated eigenvalues [small klambda, Greek, circumflex]₁ ≥ [small klambda, Greek, circumflex]₂ ≥ [small klambda, Greek, circumflex]₃ ≥ 0, which relate to the radius of gyration as R̂_g² = Tr(Ŝ) = [small klambda, Greek, circumflex]₁ + [small klambda, Greek, circumflex]₂ + [small klambda, Greek, circumflex]₃, and to the ensemble average R_g = 〈R̂_g〉. In addition to the two ring systems above, we will use R_g,0 to label the mean radius of gyration of a single ring in the infinite dilution.

To describe the shape of the rings, we will use instantaneous prolateness (and its ensemble average, = 〈〉) defined as:

which attains zero for a spherically symmetric object. Deforming the sphere by pulling its poles from each other will create an elongated cigar-like object – a prolate, for which 0 < ≤ 2, where the upper limit is reached if the object becomes a line. Conversely, deforming the sphere by pressing the poles towards each other will create a disk-like object – an oblate, for which −1/4 ≤ < 0, where the lower bound is reached if the object is a circle.

2.3 Minimal surfaces & threading analysis

To obtain the minimal surface,⁷⁷ depicted in Fig. 2a, of each of the rings in a given configuration, we use the Surface Evolver package^78,79 following the procedure introduced in ref. 37 and used also in previous publications.^52,80 First we span an initial triangular mesh over the ring, using the monomers and a set of additional ghost particles as the vertices of the mesh. The monomers of the rings stay immobilized, while ghost particles are iteratively evolved under the imposed virtual surface tension. The whole mesh is progressively restructured by vertex averaging and weeding out the small triangles, driving the final surface close to the Delaunay triangulation as detailed in Section IB in the ESI.† We tested that for the rings of our size, the algorithm always converges to the same local minima, which is surmised to be sufficiently close to the global minimum as explained in ref. 37.

Fig. 2 (a) Snapshot for fully symmetric case: N_A = N_B = 100,C_∞(A) = C_∞(B) = 15.4 at ring–ring separation r ≈ 0. The yellow ghost particles are delimiting the minimal surface of the blue (passive) ring, which is being threaded by the red (active) ring – the threader. (b) Scheme depicting the threading depths {L_t^(·)} and the separation length L_sep emerging on the (active) threader – the ring which is piercing the minimal surface of the other (passive) ring.

We say that ring A is threading ring B, if any of the bond vectors of ring A is intersecting the internal area of any of the triangles of the minimal surface mesh of ring B. In such a threading event, we would coin the ring A, the threader, as the active partner, whereas the ring B would be denoted as the passive partner. To detect the intersections we use Plücker coordinates and the side operators as elaborated in ref. 80 and Section IC in the ESI.† It has been observed that the minimal surface with a disk topology can be also threading itself, however in the current study we focus on the inter-ring threading rather than intra-ring ones.

Finally, to quantify the threading effects, we define the threading depths, {L_t}, and the separation length, L_sep, both being the characteristics of the active partner – the threader. We find the former by localizing all the intersections of the threader and the pierced surface, and then counting the number of monomers between the intersectional points ordered along the polymer contour yielding the set {L_t⁽¹⁾,L_t⁽²⁾,…,L_t^(2k−1),L_t^(2k)}. The latter is defined as , and it corresponds to the total amount of the material of the threader on either side of the minimal surface of the passive partner. Both quantities are schematically depicted in Fig. 2b.

2.4 Simulation method

To sample the configurations of the introduced model, we use LAMMPS⁸¹ implementation of Langevin dynamics. The governing equations of motion for each of the particles are

mẍ_i(t) = F_i−γmẋ_i(t) +Y_i(t), (6)

where γ is the friction coefficient and Y is a random force obeying 〈Y_i^α(t)〉 = 0 and 〈Y_i^α(t)Y_i^β(t′)〉 = 2γm_ik_BTδ_ijδ_αβδ(t − t′), where α,β ∈ {x, y, z} and δ is Kronecker delta. Finally, the force F_i is the deterministic force originating from the gradient of potentials eqn (1–3) and the external biasing potentialwith stiffness K_bias = 2.5k_BT/σ², which is introduced between the centers of mass of the two rings. For each system, we carried out a set of independent simulations for r_j ∈ {0.0, 0.5, 1.0, 1.5,…, 30.0}σ, hence allowing us to sample conformations of the rings effectively in a broad range of ring–ring separations r. For the calculation of the effective potentials, the bias was subsequently removed using the self-consistent histogram method^82–85 with a procedure described in Section ID in the ESI.† Finally, the mass of each monomer is m = 1, which sets the time unit as τ = (mσ²/k_BT)^1/2 = 1/γ = 1. We use a time step Δt = 0.01τ and the simulations were 5 × 10⁸ time steps long, which translates to 5 × 10⁶τ yielding several thousand of uncorrelated values of R_g in a typical run, as determined using the autocorrelation analysis from ref. 86.

3 Results

3.1 Effective potentials

We commence with a discussion of the effective potentials. We use only a single effective coordinate – the distance between the centers of mass of the rings, r, hence obtaining the isotropic effective potential defined as

V^eff(r) = −k_BT ln(g(r)), (8)

where g(r) is the pair correlation function between the centers of mass of the two rings in the limit of large dilution.^87,88 In Fig. 3a we show the effective potentials between two rings of length N_A = N_B = 100 for different values of the bending stiffness. The shape of the potential for flexible rings, C_∞ = 1.0, has been discussed in ref. 54 and 66 in a great detail, herein we just summarize the key observations for a reference. The presence of the other ring in the vicinity of a ring reduces the available free space, resulting a loss of conformational entropy, which increases down to the separations of r ≈ R_g,0. While flexible rings at large separations possess characteristic self-avoiding conformations, rings at r ⪅ R_g,0 typically exhibit threading-induced opening and inflation. The associated loss of conformational entropy can be rationalized in terms of free energy of confinement of ring in a cavity of size of the other ring, hence resulting in a plateau in the effective potential as explained in ref. 54. In Fig. 3a we observe that for stiff rings, the trend is qualitatively similar to the flexible ones, however the amplitude of the potentials decreases with increasing stiffness, as predicted in ref. 66. Compared to their flexible counterparts, stiff rings naturally attain more open conformations, larger in size which facilitates the threading. Additionally, distribution of monomers over larger volume decreases self-density at the center of mass of stiffer rings, further decreasing free energy penalty for placing such rings at low separations.

Fig. 3 The effective isotropic potential between the centers of mass of the two ring polymers. For the fully symmetric case, the distance is normalized by the radius of gyration of a single ring at infinite dilution. For the asymmetric cases, a similar normalization is applied, but using the average of the infinite-dilution radii of the rings from eqn (9). The analogous plots for different ring sizes, N_A,N_B are provided in Fig. S2 in the ESI.†

For the asymmetric cases, C_∞(A) = 1.0,C_∞(B) ≠ C_∞(A), let us first define the average radius of gyration

R_AB = 0.5(R_g,0(A) +R_g,0(B)), (9)

as a length scale to rationalize the interaction range. In Fig. 3b we see that the effective potentials between a stiff and a flexible ring of N_A = N_B = 100 qualitatively correspond to potentials between identical stiff rings with N_A = N_B = 100 from Fig. 3a, however the onset of the plateau is shifted to slightly larger separations. Notably, the amplitude of the potentials decreases with increasing stiffness of the ring B. Both of the observations can be understood with the help of Fig. 1b. At r ≈ 0, the stiff ring (red) attains large open conformations, while the more compact flexible ring (blue) with smaller R_g occupies the interior of the larger ring. With a stiffer ring B, the self-density of ring B in the center of mass of ring A is lower, and at the same time the confinement of ring A is weaker, because it is less likely for B sub-chains to bend and make an excursion into the pervaded volume of ring A. Albeit working in synergy, these two effects pose different free energy penalties, as we will explain in the following subsections.

Finally, for the fully asymmetric case, with one flexible ring of N_A = 100 and N_B = 25 in Fig. 3c, the trend of the previous cases is reversed – the amplitude of the effective potentials increases with the stiffness of ring B. As Fig. 1c shows, the threading roles at r ≈ 0 are now flipped and it is the large flexible ring A (blue), that is the passive partner, while the ring B (red) is the active threader. In such an arrangement, stiffer ring B poses a larger obstacle which has to be accommodated by ring A, resulting in a free energy penalty paid at the expense of conformations of ring A. Nevertheless, stiffer ring B is simultaneously more confined by ring A also contributing to the higher effective potential. An important foreshadowing from the above cases is that increasing the asymmetry between the microscopic parameters of the rings, can lead to the emergence of the distinct threading roles with one ring being almost exclusively the passive partner accommodating the active partner – the threader as will be discussed in Section 3.4.

3.2 Minimal surface, shape and size of the rings

The effective potentials reflect on the conformational changes of the rings, which we further explore by analyzing the areas of minimal surfaces of both rings, shown in Fig. 4. First, for the symmetric case, Fig. 4a shows that as the rings approach each other, they open up, manifested by the increase of the mean area of the ring surfaces. Surface inflation is notable for r ⪅ R_g,0 suggesting that the onset of threading actually emerges at distances larger than the onset of plateau of the effective potentials, which is visible only at r ≈ 0.5R_g,0. Stiffer rings possess more extended conformations and their surfaces have larger areas, compared to the flexible rings. However, in the relative terms, bringing the rings from r ≈ 3R_g,0 to center-of-mass incidence at r ≈ 0 results only in ∼5% area increase for the stiff rings with C_∞ ∼ 15.4, but in ∼13% increase for the flexible rings, signaling that the threading-induced conformational changes are more severe for the flexible rings, as also captured by higher amplitudes of the effective potentials.

Fig. 4 The mean area of the minimal surface of the rings as a function of ring–ring separation. For the fully symmetric case, we show the area averaged over both rings and the distance is normalized by the radius of gyration of a single ring at infinite dilution. For the asymmetric cases, we show areas of both rings respectively and a similar normalization of distance is applied, but using the average of the infinite-dilution radii of the rings from eqn (9). The area is normalized by the mean bond length squared. The color identifies the system, whereas line style (solid and dashed) in the asymmetric cases differentiates the rings A and B. The analogous plots for different ring sizes, N_A,N_B are provided in Fig. S4 in the ESI.†

The inflation upon ring nearing is even more visible in the upper panel of Fig. 5, where we correlate the instantaneous values of radius of gyration of the two rings for ∼1000 configurations. The purple point cloud corresponds to rings at separations r ≈ 2R_g,0, where they affect each other only weakly, hence the marginal distributions of R̂_g(A) and R̂_g(B) respectively are almost identical to the distribution of R̂_g at the infinite dilution. Compared to the reference distribution for r ≈ 2R_g,0, the green point cloud for r ≈ 0 is shifted to the higher values of radius of gyration, testifying to the swelling. We notice that the cloud at r ≈ 2R_g,0 has roughly isotropic shape, whereas the cloud for r ≈ 0 is clearly elongated in one direction (perpendicular to the dashed line). This is because in a typical configuration, the passive partner inflates more than the active partner, when compared to their conformations at large separations, as can be seen on the inset snapshots in Fig. 5. Nevertheless, we can appreciate that the whole cloud stays symmetric along the main diagonal (dashed line) because in the fully symmetric case, the rings are of the same kind, hence the average properties of the system are invariant with respect to the inversion of the rings.

Fig. 5 Instantaneous values of radius of gyration (above) and prolateness from eqn (5) (below) of ring A plotted against the corresponding values for ring B from the same configuration, shown for the fully symmetric case N_A = N_B = 100 and C_∞(A) = C_∞(B) = 1.0. The data are presented for two different mean ring–ring separations, r ≈ 0 (green) and r ≈ 2R_g,0 (purple) for ∼ 1000 configurations each. The special points marked by symbols ⊗ and ⊙ feature the instantaneous values for the selected representative snapshots in the insets – sample configurations at r ≈ 2R_g,0 and r ≈ 0 respectively. The analogous plots for the asymmetric and fully asymmetric case are in the ESI,† Sections IV and V.

An analogous analysis can be conveyed for the instantaneous values of prolateness in the lower panel of Fig. 5. Rings at r ≈ 2R_g,0 exhibit prolateness mainly in the interval 0 ⪅ ⪅ 0.25, meaning that their conformations are almost spherical, with a slight elongation towards a cigar-like object. For r ≈ 0 we witness the emergence of (green) wings in the probability distributions, spanning the regions where one ring resides in disk-like (oblate) interval −0.25 ⪅ (i) ⪅ 0, while the other in 0.25 ⪅ (j) ⪅ 0.75 holding a cigar-like shape. To conclude, in addition to the passive partner strongly inflating upon the threading, it also flattens up into a disk-like object, whereas the active partner inflates only weakly and it resembles an elongated cigar-like object.

In Sections V and IV in the ESI,† we present similar two-dimensional histograms also for the both asymmetric cases, nevertheless, due to the asymmetry in N and C_∞, it is difficult to ascribe universal meaning to the observed trends. For this reason, we will again focus on the areas of the surface as an indicator of the conformational changes. The areas in the asymmetric case, Fig. 4b, show distinct types of behaviour for the flexible ring A and stiff ring B. As shown in Fig. 1b, at small separations, it is ring A (blue) which is the active threader, piercing the surface of the ring B (red), whose minimal surface is growing in the area, when compared to the rings being far apart. On the other hand, the area of the surface of ring A is almost constant irrespectively of the ring–ring separation and stiffness of the ring B. The latter means that in the asymmetric case, it is mainly the stiff ring B, the passive partner, opening up, while the contribution of confinement free energy of the active partner is the smaller of the two. Similarly, for the fully asymmetric case in Fig. 1c, the large flexible ring A (blue), now undertaking the role of the passive partner is the one whose surface is getting inflated as the ring approach. The small stiff ring B (red), being the active threader, has almost the same area in the whole range of separations.

3.3 Probability of threading & threading roles

Finally, we quantify the threading between the rings. In Fig. 6 we show the probability that ring A is threading ring B, P(A → B), and vice versa, P(B → A), as a function of ring–ring separation. First, we note that for the symmetric case in Fig. 6a, a non-zero probability of threading emerges at separations r ⪅ 1.5R_g,0 and steadily increases up to r ⪅ 0.5R_g, where it reaches its maximum. One would expect the probability of threading at r ≈ 0 to be close to 1/2 since in the fully symmetric case, the rings are equivalent and are expected to ergodically change their threading roles of an active and passive threader. Such a change of roles is a rare event, which takes place on a time scales significantly larger than the Rouse time of a single ring. The dynamics of this transition will be the main topic of our future work, while here we focus on the static properties. Moreover, the maximal probability is actually slightly higher than 1/2, hence P(A → B) +P(B → A) > 1 at r ≈ 0, because of the cases with mutual threading, where A → B and B → A at the same time. The insets in Fig. 6 show that the probability of mutual threading at r ≈ 0 is of the order of percents, accounting for ⪅10% of the configurations even in the case of C_∞(A) = C_∞(B) = 9.4, which has the highest fraction of conformations with mutual threading. On the other hand, the cases with no threading at all at r ≈ 0 are rare, hence the most probable threading macrostates for the fully symmetric case comprise of one ring being exclusively an active partner, while the other being exclusively passive. The latter helps to contextualize the emergence of the plateau in the effective potentials in Fig. 3. It is connected to the onset of states, where the threading is immanent and the idea of confinement free energy from ref. 66 becomes valid, since one of the rings is always threading the other, hence being confined by it. Finally, the data for stiff rings in Fig. 6a show a universal trend with slightly larger range of threading interactions compared to the fully flexible rings which are systematically shifted to the lower probabilities.

Fig. 6 Probability of ring i threading ring j as a function of ring–ring separation. For the fully symmetric case, we show the probability averaged over both permutations (i → j, j → i) and the distance is normalized by the radius of gyration of a single ring at infinite dilution. For the asymmetric cases, we show the probabilities for both possible permutation respectively and a similar normalization of distance is applied, but using the average of the infinite-dilution radii of the rings from eqn (9). The color identifies the system, whereas line style (solid and dashed) in the asymmetric cases differentiates the permutation. The inset shows the probability of mutual threading at r ≈ 0 ring–ring separation. The analogous plots for different ring sizes, N_A,N_B are provided in Fig. S13 in the ESI.†

The departure of the flexible rings from the trend of the stiff ones is based in the different nature of piercing of the surface. In Fig. 7, the separation length of the active threader at r ≈ 0 is slightly lower than N/2, since the minimal surfaces of the passive ring do not always slice the active threader in half, and L_sep is defined as a minimum of the lengths at the both sides of the surface (Fig. 2). At larger separations, L_sep for flexible rings decays more rapidly than for the stiff rings, similarly to the probability of threading in Fig. 6a. The reason is that flexible rings can easily fold and slip out from the passive ring, hence L_sep attains low values more frequently with respect to the stiff rings, bringing the mean L_sep down. For instance, flexible rings can easily attain L_sep ∼ 1 when a short fold just touches the surface of the other ring, while such conformation of a subchain is penalized by bending for the stiff rings. A similar phenomenon is captured also by the inset of Fig. 7, which shows that for the fully flexible rings, the most probable threading depth is of the order of one monomer. A flexible threader pierces the minimal surface of its passive partner typically several times, back and forth with shallow threading of lengths L_t ∼ 1. On the other hand, the most probable threading depth for the stiff rings is around ≈N/2, with the peak height at L_t ≈ 0 in the inset of Fig. 7 decreasing with increasing stiffness. Accordingly, stiff rings typically exhibit a small number of deep threadings in contrast to many shallow ones of the flexible rings. We note, however, that the latter observation might be affected by finite size effect, since our stiff rings are further away from the asymptotic regime, when compared to the fully flexible counterparts.

Fig. 7 Mean separation length, as defined in Fig. 2, as a function of ring–ring separation, shown for the fully symmetric case. The distance is normalized by the radius of gyration of a single ring at infinite dilution. The inset shows the probability distributions of threading depths, as defined in Fig. 2, at r ≈ 0. The analogous plots for different ring sizes, N_A,N_B are provided in Fig. S14 (ESI†)

For the asymmetric case in Fig. 6b the clear definition of the threading roles is finally clearly visible. The probability that the flexible ring A is threading the stiff ring B (solid lines) closely approaches the unity at r ≈ 0, whereas the probability that the stiff ring B is threading the flexible ring A (dotted lines) is close to zero in the whole range of separations. For the fully asymmetric case in Fig. 6c, however, we see that the roles of A and B (solid and dashed lines) are flipped – it is the short stiff ring B that is threading the flexible ring A, while ring A is the passive threader. The difference between the asymmetric case and the fully asymmetric case is in the length of ring B. Without the loss of generality, for N_A = 100,C_∞(A) = 1.0 and C_∞(B) = 15.4, ring B of length N_B = 100 is the passive threader for the active flexible ring A, while ring B of length of N_B = 25 is the active partner for ring A. If we imagine changing N as continuous process, there should be a critical size , where the threading role at r ≈ 0 changes, and P(A → B) changes from ≈0 to ≈1 and respectively P(B → A) from ≈1 to ≈0. One would expect that for values of N_B around the critical value, , the threading roles are not clearly defined, and ring A (or ring B) can be either the active or the passive partner both with probabilities ≈1/2. How can one identify the critical will be the concern of the following subsection.

3.4 To thread or to be threaded?

In addition to the sizes and flexibilities listed in Section 2.1, we carried out a new set of simulations of the asymmetric case for N_B ∈ {30, 40, 45, 55, 60, 65, 70, 75, 80, 85} and C_∞(B) ∈ {2.6, 4.0, 6.6, 11.4, 16.9} exclusively for ring–ring separation r ≈ 0, as listed in Section IA in the ESI.† In Fig. 8a and Fig. 8b we show the threading probabilities at r ≈ 0 for different sizes N_B. As expected, for all of the flexibilities, we see that for small N_B it is the short stiff ring B that is threading the long flexible ring A, and the roles are indeed switched for sufficiently large N_B. We see that the critical length , at which P(B → A) ≈ P(A → B) ≈ 1/2 increases with decreasing stiffness of the ring B. Furthermore, the transition of the dominant macrostate from B → A to A → B is rather sharp for stiffer rings. In contrast, the rings with softer bending exhibit a larger window in which both of the above macrostates take place with probabilities of the same order of magnitude (P(B → A)/P(A → B) ∼ 1). Nevertheless, the data feature a degree of universality, as evinced in Fig. 9, where we show the relative probabilities of realization of the two above threading macrostates. The ratio of the probabilities is plotted as a function of the ratio of the radii of gyration of the involved rings in the infinite dilution limit. For all of the values of bending stiffness, if R_g,0(A) ≫ R_g,0(B), then B is threading A and vice versa, with the crossover regime of P(A → B) ≈ P(B → A) approximately coinciding with R_g,0(A) ≈ R_g,0(B). We realize however, that for stiffer rings (C_∞ ≥ 11.4), the exchange of the threading roles is slightly shifted to R_g,0(A)/R_g,0(B) ≈ 0.9 ⪅ 1.

Fig. 8 Probability of ring i threading ring j at ring–ring separation r ≈ 0 as a function of the length of ring B. We show the probabilities for both permutations respectively for the values of bending stiffness from Table 1 and also for an additional set of values from Section 1A in the ESI.† The color identifies the system, whereas line style (solid and dashed) differentiates the permutation.

Fig. 9 Ratio of threading probabilities for both permutations B → A and A → B at r ≈ 0 ring–ring separation plotted as a function of the ratio of the gyration radii of the respective rings B and A at infinite dilution. The figure essentially collapses the data from Fig. 8 by dividing the dashed line by the solid line for each color.

Evidently, one can manipulate the threading roles of the rings by modifying their stiffness or length. In principle, it should be possible to design a system with a specific dominant threading roles, once we know where is the crossover length, for the studied system. However, carrying out the full scale two ring simulations for many N_B values as in Fig. 8 and subsequent surface minimization and threading detection is time-demanding and cumbersome. For that reason we propose two arguments, how to numerically obtain an estimate of using only single ring properties N_A, N_B and C_∞(A), C_∞(B), without the need of carrying out two ring simulations.

The first method is empirical and based on Fig. 9, where we noticed that simple mean radius of gyration of individual rings in the infinite dilution is almost a quantitative indicator for the assignment of the threading role. For this reason, we constructed a new system with only one ring, and we carried out simulations to measure its radius of gyration for different flexibilities C_∞ and lengths N as shown in Fig. 10a. Next, we interpolated the data to find the value of N_B which corresponds to R_g,0(B) equal to the R_g,0(A) with N_A = 100 and C_∞(A) = 1.0. This N is our first estimate of .

Fig. 10 Auxiliary single-ring calculations for estimation of as described in Section 3.4. (a) Radius of gyration of a single ring B in its infinite dilution as a function of ring length. The solid black line denotes the R_g,0 of the ring A in the infinite dilution, hence the points where the coloured lines intersect the black line are the points where we interpolate N_B = . (b and c) Probability of ring i threading ring j estimated using the penalty argument from eqn (13) and (14). The color identifies the system, whereas line style (solid and dashed) differentiates the permutation in analogy to Fig. 8.

The second method expands on the first one as we move from comparing the mean radii of gyration to comparing their full distributions. The assumption is the following: we know from Section 3.3, that at r ≈ 0, the states with mutual threading and the states with no threading at all comprise only a few percent of all of the states. Therefore, there are only two threading macrostates with non-negligible probabilities: A → B and B → A and whichever one is associated with lower free energy is the dominant one. In our simplified picture, in the case of A → B, ring A is located inside the ring B. As a result, ring B experiences the presence of an obstacle of a mean size R_g,0(A) inside, and hence loses all of the conformations with R̂_g(B) < R_g,0(A), which would be otherwise accessible at the infinite dilution, with no presence of the threader. Similarly, the ring A, the active threader, experiences a confinement in the cavity of mean size R_g,0(B) imposed by the ring B, the passive partner. Therefore, ring A cannot attain conformations where R̂_g(A) > R_g,0(B). Let us define the penalty function of bringing the two rings from the infinite dilution and consequently imposing the threading A → B as:

where P is the probability density of radius of gyration of a single ring at infinite dilution. In the same way, we can estimate the penalty function for the threading imposed in the opposite way:as schematically depicted in Fig. S15 in the ESI.† One can indeed confirm that if R̂_g(B) < R_g,0(A) ∀ R̂_g(B) and R̂_g(A) > R_g,0(B) ∀ R̂_g(A) then [scr P, script letter P](A → B) = 1 and [scr P, script letter P](B → A) = 0, hence there is no penalty if ring B penetrates ring A but the penalty for ring A threading ring B is 1, implying that ring B always threads ring A and ring A would never thread ring B. Moreover, for the identical rings we have P_A ≡ P_B, hence [scr P, script letter P](A → B) = [scr P, script letter P](B → A) preserving the symmetry. Now let Ω > 0 be the total number of microstates at r ≈ 0, Ω(A → B) the number of microstates exhibiting the threading A → B and vice versa Ω(B → A) for B → A such that Ω = Ω(A → B) + Ω(B → A). We are now making a drastic assumption in asserting that the penalty functions [scr P, script letter P](i → j) alone can be used as estimators of the probabilities of threading and therefore of the crossover point in which the order of the probabilities reverses itself. To this goal, we first define the sum [scr Q, script letter Q] of the two threading penalties

[scr Q, script letter Q] = [scr P, script letter P](A → B) + [scr P, script letter P](B → A). (12)

Subsequently, we assert that the probability P(A → B) ≡ Ω(A → B)/Ω can be expressed asand the complementary probability P(B → A) ≡ Ω(B → A)/Ω as

Eqn (13) and (14) above express the approximation that the probability that ring i threads ring j at r ≈ 0 is given as the ratio of the penalty on j threading i over the sum of the two penalties. This approximation is indeed a drastic one and it also constitutes a dramatic simplification of the problem. The calculation of the two penalties is carried using the single ring quantities in the infinite dilution, hence this approximate method provides an estimate of a complex many-body threading probabilities using only single-ring properties. It is worth noting, however, that the method is exact at three limits, namely for the cases [scr P, script letter P](i → j) → 0 or 1 but also, by construction, for the symmetric case A = B, where the ratio on both sides of eqn (13) and (14) is equal to 1/2. Microscopically, whichever of the two threading macrostates possess more microstates is the dominant one. Based on our simplified approach, this is the state that has the lower penalty, as the latter is defined in eqn (10) and (11) above.

In Fig. 10b and c we plot, the threading probabilities obtained with eqn (13) and (14), which should be compared with the corresponding curves based on simulations of interacting rings in Fig. 8. The results on the probability based on the penalty estimate are in semi-quantitative agreement with the real ones with the general characteristic that there are steeper than the latter, i.e., they approach their asymptotic values of 0 and 1 within a narrower N_B transition domain. In parallel to Fig. 8, the points where the dashed and solid lines cross, determine the critical point, where threading A → B and B → A are equally probable, hence yielding our second estimate of .

In Fig. 11 we compare the measured in the two ring simulation with our two estimates of based on single ring calculations. For rather flexible rings (C_∞ ⪅ 6.6), both the argument using the mean R_g and the one using the penalty function comparison agree well with the two ring simulations. For the stiffer rings, our single ring ideas overestimate the actual by length of ⪅10 monomers. The reasons for the overall slight deviation of the size-based estimates of the quantity from the measured ones as C_∞(B) grows can be traced back on the results shown in Fig. 9. Indeed, although in general the value of the size ratio of the two rings and the value of the probability ratio of the threadings meet at the point (1,1), for stiff B-rings a deviation can be noticed: there is a narrow region 0.9 ⪅ R_g,0(B)/R_g,0(A) ≤ 1 for which P(A → B) > P(B → A), i.e., the larger, flexible ring A threads the smaller, very stiff ring B more frequently than the reverse. The reason for this lies in the fact that ring B is sufficiently stiff, so that its configuration is almost planar and the opening undergoes only weak fluctuations in its size and shape, offering thereby an entropically favoured pathway to be threaded by the flexible A ring, despite the latter being larger in overall size. In other words, it is not only the invariant trace of the gyration tensor of each ring that plays a role in determining the relative magnitude of the threading probabilities but additional invariants, such as the prolateness or the asphericity also play an additional role. Nevertheless, our estimates demonstrate that the dominant contribution is indeed given by the relative sizes of the two rings in isolation, providing therefore an accurate estimate for the roles of passive and active threader that the rings will undertake at separations r ≈ 0.

Fig. 11 Critical length of ring B as a function of its flexibility estimated from the full two ring simulations from Fig. 8 and using the single ring calculations from Fig. 10.

4 Conclusions

We have analyzed the conformations and interactions of a broad class of symmetric and asymmetric ring polymers of moderate sizes and of varying discrepancy in molecular weight and bending rigidity. We found that the emergence of a plateau in the effective interaction at small separations between the centers of mass is a common feature in an otherwise diverse landscape or resulting effective interactions, having its physical origin at the threading between the molecules, which is inadvertent for close approaches. To accommodate this, one of the rings opens up, assumes an oblate shape and becomes threaded by the other, which becomes prolate for symmetric or slightly asymmetric cases. For strong rigidity asymmetries, the flat and open shape of the rigid ring in comparison with the fractal shape of the flexible one, results mostly in the former assuming the role of the threaded one. Our conclusions are supported by analysis of the minimal surfaces, the shape parameters and threading probability profiles. The key result of our work is the generalization of the observed threading trends, which allow us to predict the probabilities of the threading roles at zero separation, based on properties of the rings at infinite separation, transforming thereby a 2N-body-problem into two N-body problems and thus reducing its complexity.

Possible extensions of the current work are manyfold. On the one hand, there is the realm of the concentrated-ring problems, in which it will be very interesting to see how much of the effective potential notion will continue to be of practical importance as the density grows, i.e., the role that ring deformations will play in bringing forward strong many-body effective interactions, depending on the details of the ring architecture. While it is known that the effective potential description is accurate for up to about five times the overlap concentration, whether this fully applies to the threading properties as well, deserves future verification. However, our effective theory that uses only single-ring properties suggests this is the case. Although many interesting applications are at the melt conditions, we stress that in melts the ring conformations are strongly distorted only above the entanglement length scale, while below the ring segments form loops that allow for mutual threading. It is very interesting to compare our results to the threading statistics in melts at these scales. Further interestingly, we notice that for the effective potentials between long rings with bending asymmetry, the value of the inter-species potential (flexible-stiff) is lower than the average of the intra-species potentials (flexible-flexible and stiff-stiff). Using the Flory–Huggins-like mixing argument, we would expect that flexible and stiff rings should be miscible, in contrast to their linear counterparts which are known to demix into two phases. Along similar lines, the extension of the present coarse-graining approach to the currently topical issue of poly[n]catenanes and their mixtures with linear polymers is also an interesting line for future investigations. An additional interesting future application of our work could be the attempt to predict the elastic properties of polymer networks made of ring polymers. Similarly to the case of networks of linear chains, where the significant contribution to the shear modulus comes from the topological constraints (pairwise links between chains),^89,90 we speculate that the threadings of rings can play the role of the constraints. A starting point for a general theory would be the threading propensity for two loops at different distances as a function of their size and stiffness – data that we provide here, albeit only in a dilute regime.

On the two-body level, the exploration of the effects of solvent quality, torsional degrees of freedom causing supercoiling⁹¹ as well as the quantitative description of dynamics of the exchange of threading roles are immediately relevant questions. Work on some of the above problems is currently under way. Another possible direction would be to explore the threading of two rings under nonequilibrium conditions, such as shear flow. There the ring swells due to hydrodynamic interactions coupled with the ring topology.^92,93 The interplay of this effect with the ring size changes due to threading would be very interesting to analyze. Our present work serves as the equilibrium reference point for such considerations.

Author contributions

RS, JS, CNL conceptualized the project. RS ran the simulations, conducted the investigation, analyzed the data with methodological and software (minimal surface analysis) contributions from JS. RS wrote the original draft. JS and CNL supervised the research and reviewed the manuscript.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

RS acknowledges the financial support by the Doctoral College Advanced Functional Materials-Hierarchical Design of Hybrid Systems DOC 85 doc.funds funded by the Austrian Science Fund (FWF). The computational results presented have been achieved using the Vienna Scientific Cluster (VSC). The authors acknowledge Terpsichori (Chara) Alexiou for the scientific discussions.

References

1. W. H. Carothers, J. Am. Chem. Soc., 1929, 51, 2548–2559.
2. E. Mohr, J. Prakt. Chem., 1918, 98, 315–353.
3. H. Staudinger, Ber. Dtsch. Chem. Ges., 1929, 62, 2893–2909.
4. J. R. Katz, Angew. Chem., 1929, 42, 826–831.
5. N. T. Moore, R. C. Lua and A. Y. Grosberg, Proc. Natl. Acad. Sci. U. S. A., 2004, 101, 13431–13435.
6. J. M. Deutsch, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1999, 59, R2539–R2541.
7. G. Liu, P. M. Rauscher, B. W. Rawe, M. M. Tranquilli and S. J. Rowan, Chem. Soc. Rev., 2022, 51, 4928–4948.
8. Q. Wu, P. M. Rauscher, X. Lang, R. J. Wojtecki, J. J. de Pablo, M. J. A. Hore and S. J. Rowan, Science, 2017, 358, 1434–1439.
9. P. M. Rauscher, S. J. Rowan and J. J. de Pablo, ACS Macro Lett., 2018, 7, 938–943.
10. K. Mayumi, K. Ito and K. Kato, Polyrotaxane and Slide-Ring Materials, The Royal Society of Chemistry, 2016, pp. P001–206.
11. A. Harada, Acc. Chem. Res., 2001, 34, 456–464.
12. L. F. Hart, J. E. Hertzog, P. M. Rauscher, B. W. Rawe, M. M. Tranquilli and S. J. Rowan, Nat. Rev. Mater., 2021, 6, 508–530.
13. J. Lukeš, D. L. Guilbride, J. Votýpka, A. Zíková, R. Benne and P. T. Englund, Eukaryotic Cell, 2002, 1, 495–502.
14. A. R. Klotz, B. W. Soh and P. S. Doyle, Proc. Natl. Acad. Sci. U. S. A., 2020, 117, 121–127.
15. J. Chen, P. T. Englund and N. R. Cozzarelli, EMBO J., 1995, 14, 6339–6347.
16. A. H. Fairlamb, P. O. Weislogel, J. H. Hoeijmakers and P. Borst, J. Cell Biol., 1978, 76, 293–309.
17. R. E. Jensen and P. T. Englund, Annu. Rev. Microbiol., 2012, 66, 473–491.
18. J. M. Polson, E. J. Garcia and A. R. Klotz, Soft Matter, 2021, 17, 10505–10515.
19. I. Yadav, D. Al Sulaiman, B. W. Soh and P. S. Doyle, ACS Macro Lett., 2021, 10, 1429–1435.
20. M. Bohn and D. W. Heermann, J. Chem. Phys., 2010, 132, 044904.
21. S. Jun and A. Wright, Nat. Rev. Microbiol., 2010, 8, 600–607.
22. S. Cohen and D. Segal, Cytogenet. Genome Res., 2009, 124, 327–338.
23. Q. Szabo, D. Jost, J.-M. Chang, D. I. Cattoni, G. L. Papadopoulos, B. Bonev, T. Sexton, J. Gurgo, C. Jacquier, M. Nollmann, F. Bantignies and G. Cavalli, Sci. Adv., 2018, 4, eaar8082.
24. A. Goloborodko, M. V. Imakaev, J. F. Marko and L. Mirny, eLife, 2016, 5, e14864.
25. A. Rosa and R. Everaers, PLoS Comput. Biol., 2008, 4, 1–10.
26. J. D. Halverson, J. Smrek, K. Kremer and A. Y. Grosberg, Rep. Prog. Phys., 2014, 77, 022601.
27. T. Ge, S. Panyukov and M. Rubinstein, Macromolecules, 2016, 49, 708–722.
28. A. Y. Grosberg, Soft Matter, 2014, 10, 560–565.
29. J. Smrek and A. Y. Grosberg, J. Phys.: Condens. Matter, 2015, 27, 064117.
30. J. D. Halverson, W. B. Lee, G. S. Grest, A. Y. Grosberg and K. Kremer, J. Chem. Phys., 2011, 134, 204905.
31. M. Müller, J. P. Wittmer and M. E. Cates, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1996, 53, 5063–5074.
32. M. Müller, J. P. Wittmer and M. E. Cates, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2000, 61, 4078–4089.
33. S. P. Obukhov, Talk at KITP Santa Barbara, accessed Aug. 2022, 1997, https://online.kitp.ucsb.edu/online/jamming2/obukhov/.
34. S. P. Obukhov, M. Rubinstein and T. Duke, Phys. Rev. Lett., 1994, 73, 1263–1266.
35. D. Michieletto and M. S. Turner, Proc. Natl. Acad. Sci. U. S. A., 2016, 113, 5195–5200.
36. D. Michieletto, N. Nahali and A. Rosa, Phys. Rev. Lett., 2017, 119, 197801.
37. J. Smrek and A. Y. Grosberg, ACS Macro Lett., 2016, 5, 750–754.
38. D. Michieletto and T. Sakaue, ACS Macro Lett., 2021, 10, 129–134.
39. W.-C. Lo and M. S. Turner, EPL, 2013, 102, 58005.
40. E. Lee and Y. Jung, Polymers, 2019, 11, 516.
41. J. D. Halverson, G. S. Grest, A. Y. Grosberg and K. Kremer, Phys. Rev. Lett., 2012, 108, 038301.
42. T. C. OConnor, T. Ge and G. S. Grest, J. Rheol., 2022, 66, 49–65.
43. K. R. Peddireddy, M. Lee, C. M. Schroeder and R. M. Robertson-Anderson, Phys. Rev. Res., 2020, 2, 023213.
44. D. Kong, S. Banik, M. J. San Francisco, M. Lee, R. M. Robertson Anderson, C. M. Schroeder and G. B. McKenna, Macromolecules, 2022, 55, 1205–1217.
45. D. Parisi, J. Ahn, T. Chang, D. Vlassopoulos and M. Rubinstein, Macromolecules, 2020, 53, 1685–1693.
46. D. Parisi, M. Kaliva, S. Costanzo, Q. Huang, P. J. Lutz, J. Ahn, T. Chang, M. Rubinstein and D. Vlassopoulos, J. Rheol., 2021, 65, 695–711.
47. Y. Zhou, C. D. Young, M. Lee, S. Banik, D. Kong, G. B. McKenna, R. M. Robertson-Anderson, C. E. Sing and C. M. Schroeder, J. Rheol., 2021, 65, 729–744.
48. C. D. Young, Y. Zhou, C. M. Schroeder and C. E. Sing, J. Rheol., 2021, 65, 757–777.
49. J. Smrek, J. Garamella, R. Robertson-Anderson and D. Michieletto, Sci. Adv., 2021, 7, eabf9260.
50. Q. Huang, J. Ahn, D. Parisi, T. Chang, O. Hassager, S. Panyukov, M. Rubinstein and D. Vlassopoulos, Phys. Rev. Lett., 2019, 122, 208001.
51. T. C. O'Connor, T. Ge, M. Rubinstein and G. S. Grest, Phys. Rev. Lett., 2020, 124, 027801.
52. J. Smrek, I. Chubak, C. N. Likos and K. Kremer, Nat. Commun., 2020, 11, 26.
53. I. Chubak, C. N. Likos, K. Kremer and J. Smrek, Phys. Rev. Res., 2020, 2, 043249.
54. A. Narros, A. J. Moreno and C. N. Likos, Soft Matter, 2010, 6, 2435–2441.
55. A. Narros, C. N. Likos, A. J. Moreno and B. Capone, Soft Matter, 2014, 10, 9601–9614.
56. M. A. Ubertini, J. Smrek and A. Rosa, Macromolecules, 2022 DOI:10.1021/acs.macromol.2c01264.
57. L. Tubiana, F. Ferrari and E. Orlandini, Phys. Rev. Lett., 2022, 129, 227801.
58. C. Chen, M. K. Singh, K. Wunderlich, S. Harvey, C. J. Whitfield, Z. Zhou, M. Wagner, K. Landfester, I. Lieberwirth, G. Fytas, K. Kremer, D. Mukherji, D. Y. W. Ng and T. Weil, Nat. Commun., 2021, 12, 3959.
59. B. A. Krajina, A. Zhu, S. C. Heilshorn and A. J. Spakowitz, Phys. Rev. Lett., 2018, 121, 148001.
60. F. Guo, K. Li, J. Wu, L. He and L. Zhang, Polymers, 2020, 12, 2659.
61. E. Orlandini and C. Micheletti, J. Phys.: Condens. Matter, 2021, 34, 013002.
62. N. Hirayama, K. Tsurusaki and T. Deguchi, J. Phys. A: Math. Theor., 2009, 42, 105001.
63. A. Narros, A. J. Moreno and C. N. Likos, Macromolecules, 2013, 46, 3654–3668.
64. I. Chubak, C. N. Likos and J. Smrek, Mol. Phys., 2021, 119, e1883140.
65. M. D. Frank-Kamenetskii, A. V. Lukashin and A. V. Vologodskii, Nature, 1975, 258, 398–402.
66. M. Bernabei, P. Bacova, A. J. Moreno, A. Narros and C. N. Likos, Soft Matter, 2013, 9, 1287–1300.
67. P. Poier, C. N. Likos, A. J. Moreno and R. Blaak, Macromolecules, 2015, 48, 4983–4997.
68. P. Poier, P. Bačová, A. J. Moreno, C. N. Likos and R. Blaak, Soft Matter, 2016, 12, 4805–4820.
69. M. Z. Slimani, P. Bacova, M. Bernabei, A. Narros, C. N. Likos and A. J. Moreno, ACS Macro Lett., 2014, 3, 611–616.
70. J. H. Choi, T. Kwon and B. J. Sung, Macromolecules, 2021, 54, 11008–11018.
71. J. Kim, J. M. Kim and C. Baig, Soft Matter, 2021, 17, 10703–10715.
72. P. Chiarantoni and C. Micheletti, Macromolecules, 2022, 55, 4523–4532.
73. K. Kremer and G. S. Grest, J. Chem. Phys., 1990, 92, 5057–5086.
74. M. Rubinstein and R. H. Colby, Polymer Physics, Oxford University Press, 2003.
75. M. Liebetreu and C. N. Likos, Soft Matter, 2020, 16, 8710–8719.
76. M. Liebetreu and C. N. Likos, ACS Appl. Polym. Mater., 2020, 2, 3505–3517.
77. M. Lang, Macromolecules, 2013, 46, 1158–1166.
78. K. A. Brakke, The Surface Evolver Version 2.70 August 25, 2013, https://facstaff.susqu.edu/brakke/evolver/evolver.html/, 1989-2013, [Online; accessed 01-August-2022].
79. K. A. Brakke, Exp. Math., 1992, 1, 141–165.
80. J. Smrek, K. Kremer and A. Rosa, ACS Macro Lett., 2019, 8, 155–160.
81. S. Plimpton, J. Comput. Phys., 1995, 117, 1–19.
82. A. M. Ferrenberg and R. H. Swendsen, Phys. Rev. Lett., 1989, 63, 1195–1198.
83. S. Kumar, J. M. Rosenberg, D. Bouzida, R. H. Swendsen and P. A. Kollman, J. Comput. Chem., 1992, 13, 1011–1021.
84. D. Frenkel and B. Smit, Understanding Molecular Simulation: From Algorithms to Applications, Academic Press, 2002.
85. M. Andrec, The Weighted Histogram Analysis Method (WHAM), 2010, https://core.ac.uk/download/pdf/22873967.pdf.
86. W. Janke, Quantum simulations of complex many-body systems: from theory to algorithms, 2002, vol. 10, pp. 423–445.
87. R. Blaak, B. Capone, C. Likos and L. Rovigatti, Accurate Coarse-Grained Potentials for Soft Matter Systems, 2015, pp. 209–258.
88. C. N. Likos, Phys. Rep., 2001, 348, 267–439.
89. W. W. Graessley and D. S. Pearson, J. Chem. Phys., 1977, 66, 3363–3370.
90. R. Everaers and K. Kremer, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1996, 53, R37–R40.
91. K. R. Peddireddy, M. Lee, Y. Zhou, S. Adalbert, S. Anderson, C. M. Schroeder and R. M. Robertson-Anderson, Soft Matter, 2020, 16, 152–161.
92. M. Q. Tu, M. Lee, R. M. Robertson-Anderson and C. M. Schroeder, Macromolecules, 2020, 53, 9406–9419.
93. M. Liebetreu and C. N. Likos, Commun. Mater., 2020, 1, 4.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d2sm01177h

---