2.1.

Coarse-Grained Molecular Dynamics Model

A previously developed coarse-grained molecular dynamics model is used in this work.^25,26 Each polymer is composed of unified beads or atoms that are each equivalent to four monomeric units. Polymers are coarse grained at this level since the statistical segment length of many common polymers, such as poly(styrene) and poly(methyl methacrylate), is less than four monomers. An illustration of this coarse-graining can be seen in Fig. 1(a). The beads in a polymer are connected to each other using a harmonic bond potential of the form

Fig. 1

(a) An illustration of the coarse-graining used in this work. Blue atomistic beads represent carbon, white represent hydrogen, and red represent oxygen. The blue and red transparent bubbles represent one coarse-grained polystyrene and poly(methyl methacrylate) bead, respectively. The coarse-grained bead is placed at the center of mass of the four monomeric units. Parts (b) through (d) show the force that the three potentials in this work apply to beads. (b) A harmonic bond potential is applied to every pair of bonded beads. (c) A harmonic angle potential is applied to every three consecutively bonded beads. (d) A nonbonded potential is applied to every pair of beads, excluding beads in the same bond or angle potential. In this case, the nonbonded forces between the dark blue bead and the other neighboring six beads shown are illustrated.

In addition to the bond potential, every three consecutively bonded beads have a harmonic angle potential acting on them of the form

Except for bead pairs that already interact through a bond and/or an angle potential, every pair of beads experiences a nonbonded potential (i.e. a bead does not experience a nonbonded potential with a nearest neighbor or next-nearest neighbor bead, if they exist, on its polymer chain). The nonbonded potential is of the form

In this paper, during thermodynamic integration, ${ϵ}_{AB}$ is varied from a starting value of ${ϵ}_{AB}=0.5\mathrm{kcal}/\mathrm{mol}$, where $\chi =0$, to an ending value of ${ϵ}_{AB}=0.35\mathrm{kcal}/\mathrm{mol}$, where $\chi =0.55$. The polymers used have a degree of polymerization ($N$) of 64 monomers (16 beads) making $\chi N$ have a value of 35. The value of $\chi$ used here is larger than that of most BCPs currently used in the literature.²⁷ The free energy of defects has been shown to have a strong dependence on $\chi$;²⁴ therefore, the free energy measured in this paper will likely be higher than most known BCPs. However, this higher value of $\chi$ should not affect how the free energy changes with the size of the defect and the design of the underlayer, which are the trends observed in this work.

Simulations are built using code written in MATLAB.²⁸ Molecular dynamics simulations are run using HOOMD-Blue on a GPU cluster.^29,30 Simulation results were viewed using code written in MATLAB.

Simulations in this work are composed of two pieces, the underlayer and the film. The underlayer is a brush underlayer, with each brush being a seven bead ($N=28$) chain that has one end fixed in space. Brushes are placed with a brush density of $0.44{\text{brushes}/\mathrm{nm}}^2$ and a square grid (leading to a grid spacing of 1.5 nm). A hard surface is simulated below the brushes by placing a second layer of fixed beads.²⁶ These fixed beads are located 0.83 nm beneath the brushes on a square grid that is offset from the brushes so that each bead is in the center of four brushes. This layer of fixed beads helps prevent any chains from inverting and facing downward into the substrate. The brush is composed of the same bead types as are used in the BCP film; however, they are labeled differently, so they do not experience the same external potentials that will be discussed later. To simulate chemoepitaxy, the underlayer types are assigned based on position. The underlayers being simulated are composed of an alternating pattern of a thin, highly preferential pinning stripe, and a wider neutral region. The pinning stripes are always $0.5·L_0$ wide, whereas the background region is $(n-0.5)·L_0$ wide, where $n$ is the degree of density multiplication. In this work, the background region is always composed of 50% A beads and 50% B beads. The pinning stripe is always 100% of either A or B beads depending on which bead type needs to be the pinning type for the simulation. Illustrations of three of these underlayers are shown in Fig. 2. Figure 2(a) shows an $8\times$ density multiplying underlayer with its only pinning stripe located at position 0 (the center of the defect). Figure 2(b) shows the same density multiplication but with the pinning stripe in position 2 (two lamellae from the center of the defect). Figure 2(c) shows the same pinning stripe position (PSP) as Fig. 2(b) but on a $4\times$ density multiplying underlayer that has an additional pinning stripe located at position $-6$.

Fig. 2

Illustrations of three different patterned underlayers used in this work. The top portion of each illustration is the BCP film with a dislocation of order 1. The bottom portion represents the brush underlayer. The background region is represented by the gray region and is neutral to both the white and black blocks of the BCP. The pinning stripe is the white region and is highly preferential to the white block of the BCP. (a) An $8\times$ density multiplying underlayer with pinning stripe position 0. (b) An $8\times$ density multiplying underlayer with pinning stripe position 2. (c) An $4\times$ density multiplying underlayer with pinning stripe position 2.

BCP films are generated by creating chains using a random walk. The chains are then randomly placed with their center of mass inside a simulation box with dimensions of the desired film until the appropriate density is reached. For the polymer used here, the density is $1.35{\text{beads}/\mathrm{nm}}^3$. Next, the BCP chains are periodically wrapped in the $x$- and $y$-dimensions. Finally, the film is placed on top of the desired underlayer. These films are initially well mixed. The majority of films were built with dimensions of $8·L_0\times 6·L_0\times 0.75·L_0$, where $L_0$ is the pitch of the BCP, which is 11.8631 nm for this BCP.

Simulations undergo two minimization steps after being built. The first minimization gently pushes apart beads that were placed too close together initially. To do this, ${\sigma }_{ij}$ is initially set to a very low value, which prevents beads from feeling too high of a force from the highly repulsive regime of the nonbonded potential. A series of short minimizations (50 steps each) are run, increasing ${\sigma }_{ij}$ in each run until the desired value is reached. Minimizations are run using the HOOMD FIRE minimizer (parameters: $dt=5\times {10}^{-6}$, $\mathrm{ftol}=1\times {10}^{-2}$, $\mathrm{Etol}=1\times {10}^{-7}$, $\mathrm{finc}=1.99$, $\mathrm{fdec}=0.8$, $\text{alpha}\_\text{start}=0.01$, and $\mathrm{falpha}=0.9$).

The second minimization step is a more generic minimizer run at the final parameter values. This minimization runs for 20,000 steps using the FIRE minimizer with parameters $\mathrm{dt}=5\times {10}^{-4}$, $\mathrm{ftol}=1\times {10}^{-2}$, and $\mathrm{Etol}=1\times {10}^{\mathrm{e}-7}$. After completing both minimizations, a brief simulation is run with $\chi =0$ for 200,000 timesteps (10 ns) to collapse the film to the brush underlayer and allow the film to better equilibrate, eliminating any unrealistic density fluctuations generated during the initial build. This simulation is run with the HOOMD standard NVT integrator, which uses a Nose–Hoover thermostat. For these simulations, the temperature set point is $T=500\mathrm{K}$, the controller coupling constant is $\mathrm{tau}=0.2$ timesteps, and the integration timestep is 0.05 ps.

2.2.

Thermodynamic Integration

Thermodynamic integration is a method that can be used to calculate the free energy difference between two different states in a molecular dynamics simulation. In thermodynamic integration, the molecular dynamics simulation is brought from one state to another in a reversible manner. In this paper, thermodynamic integration is used to calculate the difference in the free energy of a defective BCP film versus a defect-free lamellar state. The thermodynamic path used is summarized in Fig. 3. A simulation is started in the mixed state at $\chi =0$. This simulation is then phase separated using only an external potential by ramping the strength of the external potential from $A=0$ to $A=1$. Next $\chi$ is ramped from a value of $\chi =0$ to $\chi =0.55$. Finally, the external potential is turned off (by ramping $A=1$ to $A=0$). This is performed once using a defect-free external potential and a second time using a defect external potential.

Fig. 3

Images of the thermodynamic integration pathway used. Simulations start in a mixed state and then phase separated using an external potential. In the next branch, $\chi$ is turned on, which primarily just sharpens the interfaces. Finally, the external potential is turned off, which allows the film to relax and reach its equilibrium state.

The external potential used to generate defect-free lamellae is given by the following equation:

In the first branch, the external potential is turned on by ramping $A$ from $A=0$ to $A=1$. It was found that $A=1$ is a sufficient strength for the external potential to fully phase separate the BCP film. The external potential was turned on in 50 steps, with a stepsize of $\mathrm{\Delta }A=0.02$. Each step in $A$ was run for 1,000,000 timesteps or 50 ns. It was found that equilibrium was typically reached within the first few nanoseconds. The external potential in the simulation was then logged every 100 timesteps. The external potential energy was averaged over the final 45 ns of each step. The integral in Eq. (6) was then evaluated using the trapezoidal rule.

The second branch ramps on $\chi$ while keeping the external potential at full strength. The value of ${ϵ}_{AB}$ is decreased in this branch in increments of $\mathrm{\Delta }{ϵ}_{AB}=0.005\mathrm{kcal}/\mathrm{mol}$. Each step in ${ϵ}_{AB}$ is run for 80,000 timesteps or 4 ns. Every 100 timesteps the nonbonded potential between $A$ and $B$ beads is measured. This value is averaged over the final 2.5 ns to give $⟨V_{\mathrm{non}-\mathrm{bonded}}^{AB}({ϵ}_{AB})⟩$, which is then used in the following integral to get the free energy of this branch:

The final branch ramps down the external potential strength. The value of $A$ is decreased from $A=1$ to $A=0$ in 30 steps. Each step in $A$ is run for 80,000 timesteps or 4 ns. Every 100 timesteps the nonbonded potential between $A$ and $B$ beads is measured. This value is averaged over the final 2.5 ns to give $⟨V_{\mathrm{ext}}(A)⟩$, which is then used in the following integral to get the free energy of this branch:

The difference in free energy of one path versus another is given by the following equation:

Thermodynamic integration is a useful tool for measuring free energy in molecular dynamics models. However, there are a few significant drawbacks. First, as is the case in all free energy models, the size of the simulation volume can affect the free energy of the system.²³ This effect is lessened by using larger simulation volumes. Second, thermodynamic integration gives increases in accuracy by increasing the number of integration steps and by increasing the number of timesteps to average over at each integration step. Both of these methods of increasing the accuracy of free energy measurements involves a significant increase in simulation time. In addition, in this work there is the possibility that a defect will anneal out of a film (leaving the film in a defect-free state) when the external potential is turned off. While the defect annealing is a rare event, if the simulation is run longer (to increase the accuracy of the measurement), then the likelihood of this rare event occurring increases. When the defect does anneal, it violates reversibility, invalidating the free energy measurement.

3.1.

Pinning Stripe Location

Figure 5 shows the free energy difference for an $8\times$ density multiplying underlayer for various PSPs. In these cases, there is only one pinning stripe in the system. The position of the pinning stripes is their distance from the center of the defect (divided by $0.5·L_0$), as shown in Fig. 4. It should be noted that the odd numbered defects have pinning stripes centered on integer values ($n$), while even numbered defects have pinning stripes on the half values ($n+0.5$). The individual branches of the thermodynamic integration for the second branch (ramping up $\chi$) and third branch (turning off the external potential) are shown in Figs. 5(a) and 5(b), respectively, whereas the overall free energy difference is shown in Fig. 5(c). The first branch of the integration is not plotted since $\mathrm{\Delta }{F}_1$ is not a function of pinning stripe location. This is because the entire first step is done with $\chi =0$, so there is no differentiation between underlayer beads. The values of $\mathrm{\Delta }{F}_1$ are 895, 898, 957, 1018, and $1168\mathrm{kcal}/\mathrm{mol}$ for defects of order 1, 2, 3, 4, and 5, respectively. All points shown in the figure and the data point mentioned above are the average of three replicates.

Fig. 5

Plots showing (a) the second branch of the free energy calculation, (b) the third branch of the free energy calculation, and (c) the total free energy difference between the defect and a defect-free state versus the position of the pinning stripe for an $8\times$ density multiplying underlayer.

Looking at the overall free energy in Fig. 5(c), it can be seen that the dislocations have similar trends with the exception of the dislocation with a DO of 1. This is due to many $\mathrm{DO}=1$ simulations annealing out during the third branch, that is to say, the defects were annihilated, leaving behind straight lamellae. While defect annihilation (as well as defect growth) is possible for all DOs, it was only observed in the $\mathrm{DO}=1$ case due to the close proximity of the pair of dislocations of which $\mathrm{DO}=1$ is composed. Due to this close proximity, one dislocation will “feel” when the shape of the other dislocation is perturbed by a nearby pinning stripe. For a larger DO, there are multiple lamellae separating the dislocations that damp out the perturbation from the pinning stripe. However, in the $\mathrm{DO}=1$ case, there is little material to damp it out. While a dislocation pair can withstand some perturbation in one dislocation, when both are altered, it becomes far more likely that the dislocation will change in order.

The effect of $\mathrm{DO}=1$ annealing out in the third branch when a pinning stripe is close can be seen in the third branch of the thermodynamic integration [Fig. 5(b)]. During the third branch, there should be almost no rearrangement of beads in the system if the defect that was built using the external potential is close to the true equilibrium defect shape. If there is no rearrangement, the third branch should be constant with respect to the PSP. This is approximately the case for all defects other than $\mathrm{DO}=1$, as will be discussed further later. As mentioned before, when the pinning stripe approaches the center of $\mathrm{DO}=1$, the defect will frequently anneal out while ramping down the external potential. This spontaneous transition violates reversibility, which invalidates the free energy measurement from thermodynamic integration in these cases.

On the other hand, for defects other than $\mathrm{DO}=1$, the third branch is mostly constant with respect to PSP. There are some variations, in particular, a lowering of $\mathrm{\Delta }{F}_3$ when the pinning stripe approaches the terminating block of the defect. In these cases, the defect will necessarily adapt its shape slightly to accommodate the pinning stripe. However, these variations are small in magnitude and gradually happen in a reversible manner, allowing thermodynamic integration to still be valid. Since the $\mathrm{\Delta }{F}_3$ is roughly constant with regard to PSP, and what few variations are present appear small, it can be assumed that the second branch of the integration, $\mathrm{\Delta }{F}_2$, is the primary contributor to the total free energy with respect to the PSP. Therefore, while the total free energy of $\mathrm{DO}=1$ is likely invalid, the second branch of the integration can be looked at as a proxy for the total free energy.

A few observations can be made by analyzing the second branch of the integration [Fig. 5(a)]. First, as the pinning stripe moves further from the defect, toward higher PSP values, the free energy begins to level off at its minimum value. This is due to the pinning stripe interacting less and less with the defect and, therefore, having a negligible effect on the free energy. The PSP where the free energy levels off increases with an increasing DO but appears to be within $0.5·L_0$ of the terminating block. While $\mathrm{DO}=1$ seems to level off around $\mathrm{PSP}=3$, each increase in DO seems to increase this transition by approximately 0.5 PSP. This suggests that pinning stripes beyond a distance of $L_0$ outside of the terminating block (marked by the gray dots in Fig. 4) have little effect on the free energy of the defect.

As the pinning stripe approaches the defect, the free energy greatly increases due to the increased interactions of the pinning stripe with the defect. The magnitude of this increase is on the order of 50 to $100\mathrm{kcal}/\mathrm{mol}$, implying that these defects will almost never occur naturally above the pinning stripe at equilibrium.

For some of the larger defects, such as $\mathrm{DO}=3$ or 5, there appears to be a maximum in the free energy when the pinning stripe is not at the center of the defect ($\mathrm{PSP}=0$). The position of these maximums are roughly $\mathrm{PSP}=1$ and 2 for $\mathrm{DO}=3$ and 5, respectively, implying that the maximum likely shifts further right as the DO is increased. The location of the maximum corresponds to the lamellae interior to the terminating block of the defect, which is at $\mathrm{PSP}=2$ and 3 for $\mathrm{DO}=3$ and 5, respectively. Since the pinning stripe being directly under or on the interior of the terminating block consistently gives the maximum in free energy, it suggests that a pinning stripe in these locations destabilizes the defect. This pinning stripe location also likely has the highest kinetics for defect annihilation since a pinning stripe in this location will increase the likelihood of bridge formation, which is a necessary step in defect annihilation.³¹

The location where the free energy levels off can be a good approximation for what density multiplication is large enough to allow the defect to form. For example, since $\mathrm{DO}=1$ levels off around $\mathrm{PSP}=3$, it can be assumed that this dislocation is more likely to form if the density multiplication is $3\times$ or greater since the defect can form entirely between two pinning stripes. However, if the density multiplication is less than $3\times$, then there is guaranteed to be a pinning stripe in either $\mathrm{PSP}=0$, 1, or 2, which are in the increased free energy region. Similarly, for $\mathrm{DO}=3$ the free energy levels off around $\mathrm{PSP}=4$, and for $\mathrm{DO}=5$ it levels off around $\mathrm{PSP}=5$. This implies that decreasing the density multiplication should decrease the number of the larger defects significantly. However, since the maximum in the free energy is not at $\mathrm{PSP}=0$ for these larger defects, there exists the possibility of a dislocation pair that straddles a pinning stripe. For example, with $\mathrm{DO}=5$, if the pinning stripe is at $\mathrm{PSP}=0$, the dislocation pair can fit nicely on a $3\times$ underlayer with one pinning stripe in the center of the defect and another two on the outside of the dislocation pair (at positions $-6$ and $+6$). This particular example defect has been naturally observed in simulations using this model and has been shown to be very stable since the defect persisted without changing the DO for very long simulation times.

Looking at the overall free energy, $\mathrm{\Delta }F$, Fig. 5(c), where the pinning stripe is far from the defect (approaching $\mathrm{PSP}=8$), it can be seen that the free energy is lowest for $\mathrm{DO}=1$, increases for $\mathrm{DO}=2$, increases further for $\mathrm{DO}=3$, stays approximately the same for $\mathrm{DO}=4$, and then decreases for $\mathrm{DO}=5$. It is hypothesized that this trend is due to two competing factors. First, increasing the DO increases the amount of curvature in the BCP system, which enthalpically would increase the free energy. Second, the two halves of the dislocation pair likely destabilize each other when they are too close. This hypothesis would lead to the observed results where initially the free energy increases from defect free to $\mathrm{DO}=1$ to $\mathrm{DO}=3$ due to the increasing interfacial area. When $\mathrm{DO}=4$ is reached, the increase in interfacial area is counteracted by the stabilizing effect of the dislocation pairs being further separated. Finally, with $\mathrm{DO}=5$, the defect begins to get more stable due to the further separation of the dislocation pairs. This has large practical consequences on defect annihilation kinetics. If it is assumed that the pathway to annealing $\mathrm{DO}=5$ is for the defect to first shrink to $\mathrm{DO}=4$, then to $\mathrm{DO}=3$, and so on, as has been observed in simulation,^31,32 then $\mathrm{DO}=4$ could be considered a transition state in the process of annealing $\mathrm{DO}=5$ to a defect-free state. While $\mathrm{DO}=5$ has a much higher free energy than defect free ($\approx 300\mathrm{kcal}/\mathrm{mol}$) and should therefore be virtually nonexistent at equilibrium, the fact that the transition state $\mathrm{DO}=4$ has a far higher free energy ($\approx 340\mathrm{kcal}/\mathrm{mol}$) suggests that there will be very slow kinetics in annealing out $\mathrm{DO}=5$. These data help support the theory that defects in BCP films are kinetically trapped.

It has been shown that thermodynamic integration simulations can suffer from a limited simulation volume.²³ This was a concern, particularly for the larger DOs since it is more likely that the defect could see itself across the periodic boundary. However, it should be noted that the presence of a pinning stripe near the periodic boundary of the system lessens this concern since the pinning stripe decreases the amount a defect could influence itself across the boundary. To test the simulation volume effect, a limited number of larger simulations were run. These larger simulations were $12·L_0$ wide, while the typical simulations were only $8·L_0$ wide. These larger simulations took roughly 1.5 times as long as the smaller simulations. A single replicate was run for DOs 3, 4, and 5 where the size effect should be most apparent. This is also where the free energy begins to level off and then decrease, as detailed in the previous paragraph, so these larger simulations should help verify that trend. An underlayer with $12\times$ density multiplication was used, with the pinning stripe located away from the defect on the edge of the simulation. The results of the larger volume simulations showed the same general trends as were observed in the smaller simulations, with $\mathrm{DO}=3$ and $\mathrm{DO}=4$ having roughly equivalent free energies ($\mathrm{\Delta }F=256$ and $\mathrm{\Delta }F=250\mathrm{kcal}/\mathrm{mol}$, respectively) and $\mathrm{DO}=5$ having significantly lower free energy ($\mathrm{\Delta }F=170\mathrm{kcal}/\mathrm{mol}$). However, while the trend is the same, the absolute values of the free energy differences were approximately $100\mathrm{kcal}/\mathrm{mol}$ lower than before. This suggests that the smaller simulations are indeed experiencing this size effect; however, it appears they do capture trends accurately.

3.2.

Density Multiplication

Underlayers with $1\times$, $2\times$, $4\times$, and $8\times$ density multiplication were simulated. In addition, an unpatterned underlayer was simulated to approximate an infinite density multiplication. These results are shown in Fig. 6 where the free energy differences between a defective film and a defect-free film are shown for various DOs as a function of density multiplications. The points shown for each density multiplication are the free energy difference corresponding to the lowest average free energy over all PSPs for that density multiplication. It can be seen that, as the density multiplication increases, the difference in free energy decreases. This is due to higher density multiplications having less pinning stripes per area and, therefore, a lowered driving force for pattern correction when a defect is present. Intuition would suggest that the free energy difference should asymptotically approach the unpatterned (infinite density multiplication) case as the density multiplication increases, which appears to be the case in these simulations.

Fig. 6

Free energy differences for thermodynamic integration for various density multiplications. Each point represents the minimum free energy difference of all pinning stripe positions for that density multiplication. The horizontal dashed lines near the bottom represent the free energy difference of the defect measured on an unpatterned underlayer.

Assuming a Boltzmann’s distribution, the defect density can be estimated using the free energy difference and Eq. (1). From Fig. 6, these free energy differences can be used to determine what the return is for putting in the increased lithographic effort to make a more highly defined underlayer pattern. For instance, for $\mathrm{DO}=5$, the difference in free energy between $1\times$ and $2\times$ density multiplication is approximately $100\mathrm{kcal}/\mathrm{mol}$. Since in these simulations $k_B·T\approx 1\mathrm{kcal}/\mathrm{mol}$, it is $e^{100}$ times more likely to find this defect on a $2\times$ density multiplication than $1\times$. This type of comparison can be useful when looking at higher density multiplications to determine when increasing density multiplication will no longer have an effect on the equilibrium defect density. For all defects shown here, the $8\times$ underlayers had nearly the same defect-free energy difference as the unpatterned cases, meaning that at equilibrium the defect density for any density multiplication greater than $8\times$ will all be nearly the same. With the exception of $\mathrm{DO}=2$, decreasing the density multiplication from $8\times$ to $4\times$ increased the free energy difference by $20\mathrm{kcal}/\mathrm{mol}$ or more. This suggests that these defects will occur at least $5\times {10}^8$ times more on an $8\times$ underlayer than a $4\times$ underlayer. On the other hand, for $\mathrm{DO}=2$, the defects are still small enough that a decrease in density multiplication from $8\times$ to $4\times$ has no effect on the free energy, and even higher smaller density multiplications would be needed to affect the free energy.

Figure 7 takes a closer look at the effect of density multiplication in regard to PSP for $\mathrm{DO}=3$ by looking at the second branch of thermodynamic integration. The previous section has already discussed the 8× underlayers and has explained why having a pinning stripe in positions 1 or 2 ($\mathrm{PSP}=1$ or $\mathrm{PSP}=2$) has the highest free energy difference due to its proximity to the terminating block. As the density multiplication decreases, the number of pinning stripes present in the underlayer increases. For $4\times$ underlayers, there are two pinning stripes throughout our simulation, yet the $4\times$ underlayers mostly follow the same trend as $8\times$ underlayers. This implies that the extra-pinning stripe is still too far away from the defect’s terminating blocks to cause it any substantial instability.

Fig. 7

Free energy differences for the second branch in the thermodynamic integration for $\mathrm{DO}=3$. The unpatterned underlayer free energy is represented by a horizontal dashed line at the bottom of the graph.

For $2\times$ underlayers, four pinning stripes are present in the simulations underlayer. It can be seen that the free energy significantly increases as the PSP shifts to the right. This trend can be explained by looking at the locations of all the pinning stripes in the system, but especially the two closest to the defects terminating blocks. For $2\times$ density multiplication, the pair of pinning stripes are always $\mathrm{PSP}=x$ and $\mathrm{PSP}=x\pm 4$ (see Fig. 4, Defect Order 3 for reference). In the case of $\mathrm{PSP}=0$ and $\pm 4$, the closest distance from one pinning stripe to a terminating block is relatively far ($L_0$). On the other hand, $\mathrm{PSP}=1$ and $-3$ has both pinning stripes $0.5·L_0$ away from the terminating block (one on the inside and one on the outside). While in the $8\times$ simulations, a pinning stripe in position 3 does not significantly increase the free energy, a pinning stripe in position 1 does, which explains why this point is higher in free energy than $\mathrm{PSP}=0$ and $\pm 4$. A similar analysis can be done for $\mathrm{PSP}=2$ and $-2$. Since both of these pinning stripes are beneath a terminating block, and from the $8\times$ simulations we see that this is near where the maximum in the free energy occurs, it follows that this underlayer would have the highest free energy of the $2\times$ underlayers.

A similar analysis can be performed for the $1\times$ underlayer. In this case, the nearest pinning stripes are at $\mathrm{PSP}=x\pm 2$. Using this, when $\mathrm{PSP}=0$ and $\pm 2$, there are two pinning stripes directly underneath the terminating block as well as one pinning stripe in the center of the defect, which greatly contributes to increasing the free energy. When $\mathrm{PSP}=-3$, $-1$, 1, and 3 while there are two pinning stripes in highly unfavorable spots ($\mathrm{PSP}=\pm 1$), the two other pinning stripes are getting to the region where they will have less effect on the free energy, causing this position to have an overall lower free energy.