Potential of Mean Force between Bare or Grafted Silica/Polystyrene Surfaces from Self-Consistent Field Theory

We investigate single and opposing silica plates, either bare of grafted, in contact with vacuum or melt phases, using self-consistent field theory. Solid–polymer and solid–solid nonbonded interactions are described by means of a Hamaker potential, in conjunction with a ramp potential. The cohesive nonbonded interactions are described by the Sanchez-Lacombe or the Helfand free energy densities. We first build our thermodynamic reference by examining single surfaces, either bare or grafted, under various wetting conditions in terms of the corresponding contact angles, the macroscopic wetting functions (i.e., the work of cohesion, adhesion, spreading and immersion), the interfacial free energies and brush thickness. Subsequently, we derive the potential of mean force (PMF) of two approaching bare plates with melt between them, each time varying the wetting conditions. We then determine the PMF between two grafted silica plates separated by a molten polystyrene film. Allowing the grafting density and the molecular weight of grafted chains to vary between the two plates, we test how asymmetries existing in a real system could affect steric stabilization induced by the grafted chains. Additionally, we derive the PMF between two grafted surfaces in vacuum and determine how the equilibrium distance between the two grafted plates is influenced by their grafting density and the molecular weight of grafted chains. Finally, we provide design rules for the steric stabilization of opposing grafted surfaces (or fine nanoparticles) by taking account of the grafting density, the chain length of the grafted and matrix chains, and the asymmetry among the opposing surfaces.


Introduction
Grafting polymer chains on solid surfaces is a standard procedure for the steric stabilization of nanocomposite systems [1][2][3]. Various methods for the experimental synthesis of such systems are reported in the literature [4][5][6]. Understanding the behavior of grafted polymer brushes requires a thorough investigation of the thermodynamics of the system under different conditions. In the case of two or more nanoparticles embedded inside a polymer melt, the challenge is to keep them separated by overcoming their tendency to form aggregates. When nanoparticles are bare, the attractive Van der Waals forces [7] drive them to come closer to each other. One of the possible ways to get around this behavior is to graft polymer chains on the surface of the nanoparticles. Achieving a proper dispersion of nanoparticles inside the polymer melt is associated with a considerable enhancement of its properties [8][9][10][11][12].
Major computational research has been conducted on systems comprising a single grafted nanoparticle embedded in a solvent or a homopolymer matrix, using theoretical formulations [13,14] or atomistic simulations [15][16][17][18]. Moreover, considerable work has addressed the behavior of grafted and matrix chains in systems comprising multiple grafted solid surfaces [19][20][21][22][23][24][25]. Munao et al. [19] demonstrated the effect of a third nanoparticle, when inserted in a system of two interacting grafted silica nanoparticles, while Martin et al. [23] investigated the effect of polydispersity of grafted chains on the structural properties of the nanocomposite system.
In this work, we undertake a self-consistent field theory [26][27][28] (SCFT)-based analysis to examine the structural properties and thermodynamics of polystyrene (PS) melt (matrix) confined between two silica plates, either grafted or bare, as well as of the same plates in the absence of melt. In fact, when no matrix chains are present between the two grafted surfaces, or they are present but their molecular weight is significantly larger than that of the grafted chains, the latter behave as if they were in contact with a very bad solvent. Following Flory's theory [29], when matrix chains are present and their molecular weight is similar to the one of grafted chains, then the system is analogous to one wherein the grafted chains are embedded inside a theta solvent. On the contrary, when matrix chains are much shorter than the grafted chains, the latter are starting to swell towards the matrix as if they were in contact with a good solvent.
Specifically, we derive the potential of mean force between the plates by varying the distance between them. Studying the thermodynamics in such a planar geometry is quite important in the field of biomembranes and other biological applications [30,31]. One could consider this planar geometry study as the equivalent of investigating the potential of mean force between spherical (fine) particles, whose radius is large enough, in comparison to chain dimensions, for their curvature to be negligible. Materials consisting entirely of matrix-free grafted plates or hairy nanoparticles (also referred to as "particle solids") exhibit interesting mechanical and optical properties, while they behave as tough glasses when assembling in specific configurations [32][33][34][35][36]. Barnett and Kumar [37,38] have published several works, where they report the use of such materials in the design of membranes for separation processes. In their recent work, Biltchak et al. [39] studied the effect of the addition of matrix chains to a neat grafted nanoparticle-based membrane on its selectivity in separations of gases of different molecular size. Mydia et al. [40] developed a two-layer theoretical model to describe the configurations of the grafted chains in the vicinity of the grafted nanoparticles and at intermediate distances between them and they compared their model with atomistic molecular dynamics simulations.
The problem of polymer grafted chains in planar geometries has been addressed in the past by several works [41,42]. By removing the incompressibility assumption and imposing Dirichlet boundary conditions at the solid surfaces, we make a step forward towards the investigation of systems with realistic interfacial free energies. Furthermore, we explicitly describe the solid-polymer and solid-solid interactions via the Hamaker potential, and we explore its influence on the potential of mean force. Without these considerations, we would not be able to apply our methodology in systems comprising exclusively grafted chains, i.e., systems in the absence of polymer matrix chains, since the incompressibility condition requires that the total segment density profile be uniform across the entire domain of interest [43]. Furthermore, such approaches-when applied to three dimensions-allow for the investigation of systems of complex geometry, where the use of Fourier based methods is not recommended, since no symmetry appears and no periodic boundary conditions can be implemented.
In the majority of computational works, it is assumed that the plates have the same grafting density and molecular weight as the grafted chains [41,42]. In reality, it is rather hard for experimentalists to prepare such a perfectly symmetric system. Herein, we increase the degrees of freedom of the system by allowing the grafted chains on each plate to have different molecular weights, whereas each plate may have its own grafting density as well. The goal is to reveal the influence of these kinds of asymmetries on the potential of mean force between the two grafted plates when varying the distance between them, and to propose a scaling law that accounts for all of these design aspects.
As a first step, we build our thermodynamic reference system by deriving the free energy of single bare or grafted plates, either isolated or in contact with polystyrene melt. Next, we demonstrate the potential of mean force of a system which contains matrix chains exclusively. The melt is assumed to be at equilibrium with a bulk melt phase at all times and it is gradually squeezed by reducing the distance between the two silica surfaces. Finally, we demonstrate results concerning the structure and potential of mean force in a system of two grafted silica plates over a large parameter space, involving the molecular weight of matrix chains, the molecular weight of grafted chains on each plate, and the grafting density of each plate. Understanding the behavior of matrix and grafted chains in this planar geometry is a stepping-stone towards the corresponding spherical case, i.e., two or more polystyrene-grafted silica nanoparticles embedded in a polystyrene matrix.
In the context of our SCFT calculations, we employ the Sanchez-Lacombe (SL) equation of state in conjunction with square gradient theory (SGT), which has been shown to reproduce correctly the density profiles and surface tension of a number of different industrially important thermoplastic polymers, when compared with atomistic molecular dynamics simulations or experiment [44]. In this work, we introduce an additional potential term which allows for a refinement of the polymer-solid interactions; by adjusting the strength of these interactions, we examine interphases with low (LW), high (HW) and perfect (PW) wetting. For comparison purposes, we evaluate several of our findings with Helfand's (HFD) free energy density [45] as well, using a constant compressibility [15]; this free energy density is most commonly used in SCFT studies whenever the incompressibility condition is relaxed. It is noted at this point, that in cases where matrix chains are not present in the system, only the SL EoS is applicable, since HFD fails to describe correctly the density profiles of grafted chains in contact with vacuum and the corresponding free energy of the system.
It is mentioned here that, since the employed model is one-dimensional, a smearing of the grafting points parallel to the solid surfaces is required. As discussed in our previous work [46], this smearing may become an inaccurate approximation when addressing grafted surfaces with low σ g R g 2 values (σ g being the surface grafting density and R g 2 the mean squared radius of gyration of grafted chains) and small nanoparticles (large curvatures). Nonetheless, the effect of smearing is weaker in the present work, since it deals with planar surfaces (zero curvature) where it is less probable for the grafted chains to adopt a mushroom-like configuration. The manuscript is structured as follows. Section 2 discusses the systems under study, presents the formulation for grafted chains in the presence and in the absence of a melt phase, whereas Section 3 provides details regarding the calculations. The first part of the Results section (Section 4.1) discusses single surfaces, either bare (polymer/vacuum and polymer/solid interphases), or grafted and in the presence or absence of melt chains. The structural and thermodynamic properties derived from these systems are used as reference quantities for the analysis of two opposing surfaces. Subsequently, we present results regarding the potential of mean force (PMF) of approaching silica surfaces, either bare (Section 4.2), or grafted, in the presence (Section 4.3) or absence (Section 4.4) of melt chains. Figure 1 illustrates a mesoscopic bead-spring representation of the systems under study. Two opposing silica plates at distance h ss are grafted with polystyrene chains, while matrix chains may or may not exist in the space between the grafted plates. Hereafter, the matrix chains are denoted by m, the chains grafted on the left silica plate are denoted by g -, and the chains grafted on the right silica plate are denoted by g + . When polymer melt exists between the two plates, it comprises a total number of n m monodispersed chains of length equal to N m monomers/segments. The melt, when present, is at equilibrium with a bulk polymer phase of chain length N m of temperature T and pressure P. The left/right surfaces are grafted with n g − /n g + chains of length N g − /N g + segments, whereas the corresponding grafting density equals σ g ∓ = n g ∓ /S solid with S solid denoting the area of each surface. The grafting density and molecular weight of grafted chains are allowed to vary between the two plates. It is obvious from the density profiles in Figure 1b,d-corresponding to systems with Figure 1a and without Figure 1c matrix chains-that, in the former case, the extension of each brush towards the opposing plate is favored. Figure 1. (a) Bead-spring representation of two opposing grafted silica walls embedded in a melt (GMG system) comprising matrix chains of length N m . The silica wall on the left (right) is grafted with σ g − S solid (σ g + S solid ) grafted chains of length N g − (N g + ). (c) The same system as (a), but in absence of matrix chains (GVG system). Additionally, the corresponding reduced density profiles from the opposing grafted silica walls in the (b) presence and (d) absence of the matrix chains are shown.

Model Systems
In the rest of this manuscript, the situation at the boundaries will be referred to either as V (vacuum), S (bare solid surface) or G (grafted solid surface), and the intermediate region between the plates will be referred to as V (vacuum) or M (melt). For example, systems comprising vacuum/melt, bare solid/melt, grafted solid/melt and grafted solid/vacuum phases will be referred to as VM, SM, GM and GV, respectively. In the same spirit, more complex systems: solid/melt/solid, grafted/melt/grafted ( Figure 1a) and grafted/vacuum/grafted (Figure 1b) will be referred to as SMS, GMG and GVG, respectively.

Self-Consistent Field Convergence Scheme
The polymer chains are described with the Gaussian thread model. Consequently, their propagation inside a three-dimensional space is described by the Edwards diffusion equation [47] subject to the field w ifc : where q c is the restricted partition function, R G,c is the radius of gyration of chains, N c is the chain length measured in monomer units, and c denotes the kind of the chains (m: matrix, g -/g + : chains grafted to the left/right side). Unless otherwise stated, the initial configuration of the field is set equal to zero across the domain. For more details concerning the iterative solution of Equation (1) in the presence of the field, the reader is referred to Section 2.1.1 of Reference [46]. In the Section 2.1.2 of the same work [46], we explain how Equation (1) reduces to its one-dimensional analogue, and in its supporting information we elaborate on the numerical solution of Equation (1) via a finite difference scheme and on the stability criteria for the proper equilibration of the field w ifc .
Dirichlet boundary conditions are imposed on the edges of the domain when they correspond to solid or gas phases and Neumann (with zero flux) when they correspond to bulk polymer phases (in cases of single interfaces). The initial conditions of matrix chains are given by q m (r, 0) = 1, while for grafted chains the following Equation (2) holds in an one-dimensional domain.
with S h g being the surface area of the plane over which the smearing of grafting points takes place (equal to S solid in planar geometries), and q m (h g , N g ) is the chain propagator for matrix chains with length N g at distance h g from a solid plate. The delta function δ h − h g is analytically evaluated as the inverse of the width ∆h of the interval assigned to the points of the mesh which are closest to the two solid surfaces and is zero everywhere else. The reduced density, ϕ c = ρ c /ρ seg,bulk , can be calculated from the corresponding restricted partition function using the convolution integral in the following Equation (3).
with ρ seg,bulk = ρ mass,bulk N A /m monomer being the segmental number density in the bulk.
With ϕ c known for each chain type, the field is estimated based on the free energy density, f [ρ, ∇ρ], from a suitable EoS: For more details about the corresponding EoS, the reader is referred to Sections 2.4.1-2.4.3.

Systems with Melt: Grand Canonical Ensemble
Following the formulation developed in our previous work [46], when melt/matrix chains exist between the two plates, then we describe the equilibrium of the system in the context of the grand canonical ensemble. Hence, the thermodynamics of the polymergrafted planar surfaces in contact with a matrix phase is described by a grand potential of the form: Note that the grand potential has been defined with respect to a reference system at the same temperature, which comprises the two solids at infinite distance from each other so that they do not interact, a bulk phase of matrix chains of length N m with chemical potential µ m N m , occupying the same polymer-accessible volume as the interfacial system, and a set of n g − and n g + unperturbed end-pinned grafted chains of length N g − and N g + , respectively, exposed to bulk polymer. The segment density of the bulk polymer at temperature T and chemical potential µ m N m is ρ seg,bulk .∆Ω coh in Equation (5) describes the nonbonded cohesive interactions among the polymer segments as governed by the free energy density of the equation of state: ∆Ω field describes the interaction between the chemical potential field and the density field [26], ∆Ω m encompasses the free energy associated with the translational and conformational entropy (with respect to the entropy of a bulk melt) of matrix chains and it is given by the following Equation (8).
∆A g is associated with the conformational entropy of n g − and n g + grafted chains subject to the field w ifc , and it is given by Equation (9).
Note that the last term in the right-hand side of Equation (9) is added in order to normalize the free energy term, ∆A g , with respect to the distance of the grafting point from the solid surface and therefore render it independent of the spatial discretization of the system [46]. U s is the interaction energy between the polymer chain segments and the left and right solid surfaces, U ss describes the direct interactions between the opposing semiinfinite solid surfaces. Expressions for u s (r) and U ss are given in Sections 2.5.1 and 2.5.2 below. For the full SCFT formulation, the reader is referred to Reference [46] and its supporting information.

Systems without Melt: Canonical Ensemble
In order to address the system where no matrix chains exist, we need to modify the above formulation and describe the thermodynamics of the system by means of the canonical ensemble. In more detail, the thermodynamics of a melt-free system with polymer chains grafted to a surface is described by the following expression for the Helmholtz free energy, (11) with A bulk being the Helmholtz energy of a reference system at the same temperature, consisting of the two solids at infinite distance from each other and a set of n g − and n g + unperturbed end-pinned grafted chains exposed to bulk melt with length N g − and N g + , respectively [46]. The first term on the right-hand side of Equation (11) describes the cohesive interactions among the polymer segments as governed by the free energy density of the equation of state, the second term describes the interaction between the chemical potential field and the density field [26]: (13) while the rest of the terms remain the same as those in Equation (5).
Since the article addresses planar geometries, it would be meaningful to report the free energies per unit area. Hereafter, the free energy terms in Equations (5)-(13) will be referred to as: γ sys α = E sys α /S solid (14) where E is Ω, A, or U, "α" denotes the kind of contribution (coh, field, m, g −/+ and s −/+ ) whereas omitting "α" corresponds to the total free energy. The "sys" symbol denotes the kind of system (see Section 2.1) and hence implicitly the statistical ensemble; i.e., a system with/without melt is described in the grand canonical/canonical ensemble. A few examples: γ VM , the total free energy of a melt in contact with vacuum (surface tension); γ SM , the total free energy of a melt in contact with a silica plate (minus adhesion tension); γ GVG coh = ∆A coh /S solid , the contribution from cohesive interactions in the polymer to the Helmholtz energy of a system consisting of two grafted silica surfaces with vacuum in between them; γ GMG coh = ∆Ω coh /S solid , the cohesive term for a melt capped between two grafted silica surfaces; etc.

Nonbonded Free Energy Density
The following subsections illustrate the contributions of the square gradient theory and the free energy density functionals to the cohesive (Equations (6) and (12)) and field (Equations (7) and (13)) terms. These equations represent the nonbonded interactions between polymer segments in the context of the field theoretical formulation applied herein.

Square Gradient Theory
Given the steepness of the field and the segment density profiles, we expand the free energy density functional governing the nonbonded interactions by considering it a functional of the density gradient too, as shown in the following Equation (15) [48,49]. The contribution of the gradient term to the free energy of the system is weighted by the influence parameter, κ. The value of this parameter is reported in Table 1. It has been obtained after a thorough comparison of SCFT with atomistic molecular dynamics simulations conducted by the authors in the past for a wide range of polymer melts [44].

Sanchez-Lacombe
The Sanchez-Lacombe (SL) equation of state [50] and the corresponding free energy density are given by Equations (17) and (18), respectively. The symbols ρ = ρ/ρ * , P = P/P * , and T = T/T * represent the reduced density, pressure and temperature, respectively, where ρ * , P * , and T * correspond to the characteristic SL density, pressure and temperature. For polystyrene, the values of those characteristic parameters are reported in Table 1.
Using Equation (17), we can derive the mass density across the bulk fluid phase, ρ mass,bulk , for each chain length from the vapor-liquid equilibrium of a SL fluid; more details in the Supplemental Material Section S1 in Reference [44].

Helfand Free Energy Density
A simpler EoS which can be used in the context of a compressible model for the description of the nonbonded polymer interactions is the Helfand (HFD) equation appearing in the following Equation (19) [45].
where κ T is the isothermal compressibility of the polymer and ρ seg,bulk is the polymer segment density in the bulk melt region. For a detailed comparison between the SL and HFD equations of state, the reader is referred to previous works of the authors [44].

Polymer-Solid Interactions
In the context of Hamaker theory [7], the attractive and repulsive interactions between a sphere of radius a M and a semi-infinite solid surface are given by Equations (20) and (21).
being the distance between the surface of the sphere and the solid surface, h being the distance of the center of the sphere from the surface, and σ SM = 0.5(σ S + σ M ) with σ S and σ M being the collision diameter of Lennard-Jones particles constituting the solid and melt, respectively. Even though the Hamaker potential allows for a proper description of the long-range interactions, it cannot by itself restore the proper interfacial thermodynamics of the system, due to spurious interactions with the self-consistent field, which can result in incorrect adhesion tension and density profiles. This issue prevails regardless of the EoS and whether SGT is enabled, as we demonstrate later on. Therefore, to circumvent the aforementioned limitations, we introduced an additional potential term which is given by Equation (22).
This term has a local character and by adjusting v ramp , it allows to tune the polymersolid interactions in the vicinity of the interface and restore the proper interfacial properties of the system.
Putting it all together, in this work the interactions between the polymer chains and the solid surface are described by the Hamaker potential in conjunction to the ramp potential: Hence allowing to preserve the long-range interactions and at the same time the short-range structure and thermodynamics.

Solid-Solid Interactions
In the context of Hamaker theory, the interactions between two semi-infinite solid surfaces separated by a medium are described by the following equation: The Hamaker constant for these interactions is usually approximated as . This approximation is, however, inappropriate for describing the surface-surface interactions in our case, since the polymer/solid interactions are computed explicitly for each surface and it would lead to double counting of these interactions. In other words, the aforementioned approximation accounts-apart from the solid-solid interactions-for the polymer/solid interactions in an effective manner [51]. This can be understood more clearly by the following example: In the trivial case where the solid and medium phases are exactly the same (in terms of structure and chemical constitution, A S = A M ), the total free energy balance (sum of the polymer-solid and solid-solid interactions) assumes a constant value that is independent of h ss . In that case, using the above approximation, A SMS = 0, and therefore, U ss = 0. However, the medium-solid interactions arising from the integration of Equation (10) vary with h ss , indicating that the free energy balance is not satisfied.
Since in our implementation the polymer-solid interactions are described explicitly, the proper treatment for A SMS is to set it A SMS = A S . This way, all kinds of interactions (polymer-solid and solid-solid) are calculated in an explicit manner and there are no double counting issues.

Calculation Details
The calculations have been performed with our in-house RuSseL code for a onedimensional formulation. The same code has been extended to address arbitrary 3D geometries as well [52]. The values of the parameters that we used in our SCFT calculations are presented in the following Table 1. The temperature, T, of the system was always equal to 500 K and when matrix chains exist in the system, the calculations were performed under the grand canonical ensemble, which we find most suitable in such interfacial solidpolymer systems. On the other hand, in absence of matrix chains, we employ the canonical ensemble to describe the equilibrium properties of the system.
The position of the grafting point is set equal to the coordinate of the first discretization node of the one-dimensional spatial mesh. Ideally, the grafting point should be attached exactly on the solid surface, but numerical issues would occur, since Dirichlet boundary conditions are imposed on the solid; hence the propagator of matrix chains appearing in the denominator of Equation (2) becomes zero on the solid surface [41,58]. In our work, we use an alternative approach to address the boundary conditions at the interfaces, rather than applying reflecting boundary conditions [41]. Mathematically preventing all chains from touching the solid surfaces via Dirichlet boundary conditions, while preserving the severity of delta function initial conditions, may seem to be an overkill in such one-dimensional and symmetric geometry systems, but it is necessary if complex systems of arbitrary geometry, dimensionality or chemistry are to be addressed in the future.
The intensity of PS-SiO 2 interactions in the context of the Hamaker potential are dictated by the corresponding σ PS , A PS and σ SiO 2 , A SiO 2 parameters, whose values are reported in Table 1. The effective radius of a PS segment was estimated as a PS = 3 3/ 4πρ seg,bulk .
The range (σ ramp ) of the ramp potential was set to a distance where intensity of the Hamaker potential becomes equal to -0.005 k B T (see critical adsorption distance in Reference [46]), and the well depth to v ramp = 0 J, -2.481 × 10 −20 J and -3.975 × 10 −20 J for cases of low (LW), high (HW) and perfect (PW) wettability conditions [54], respectively. v ramp were optimized for each W A value via the secant optimization method. Retaining the full repulsive part of this potential would result in a very steep configuration of the self-consistent field near the solid surface. For that reason, we replace the repulsive part of the solid potential with a repulsive wall positioned at a distance equal to, h HS~0 .4 nm, from the solid silica plate, where it assumes a value of u s (h HS ) = 5 k B T (exactly). Unless otherwise stated, the evaluations will be performed in the absence of the ramp potential, whereas the solid-solid interactions are active in all cases of two opposing solid surfaces.

Single Surfaces
The present section discusses the thermodynamics and structure (in terms of the brush dimensions) of vacuum/melt (VM), solid/melt (SM), grafted/melt (GM) and grafted/ vacuum (GV) interphases. The quantities extracted from these systems will be used as reference for the potential of mean force calculations involving two approaching silica surfaces, either grafted or bare. Figure 2 presents the free energy of VM (surface tension) and SM (minus adhesion tension) interphases as a function of chain length. In Figure 2a, the free energies have been evaluated with HFD. In Figure 2b the SL-SGT model is considered with the original Hamaker potential corresponding to low wetting conditions (LW), and with the addition of the ramp potential that has been adjusted to reproduce the work of adhesion of PS melts in contact with treated (high wetting, HW) and untreated (perfect wetting, PW) silica [54].

Vacuum/Melt (VM) and Solid/Melt (SM) Interphases
In interphases with low wetting, the free energies are qualitatively similar for VM and SM; they differ by about γ SM s across the chain molar mass range explored herein. This is attributed to the fact that the density profiles for these cases are very similar (e.g., compare Figure 2c with Figure 2d, and Figure 2e with the "LW" curve in Figure 2f); thus, any differences in the free energies can be attributed to the polystyrene/silica Hamaker interactions.
Increasing the strength of the ramp potential enhances the polymer-solid interactions as indicated by the more negative free energies in Figure 2b, and the more pronounced peaks of the density profiles in Figure 2f. The free energy from SL in Figure 2b appears to be an increasing function of chain length, and this behavior is anticipated, since the cohesion of the polymer increases with increasing chain length [59][60][61]. In contrast, HFD (Figure 2a) exhibits the opposite trend and this is attributed to the fact that all evaluations have been performed using a constant isothermal compressibility. In Section S1 of the Supplementary Materials we demonstrate that tuning the HFD compressibilities (with and without employing SGT) based on the predictions of SL, or even fitting them directly to the experimental surface tensions, allows to restore the proper chain length dependence in the thermodynamic behavior of the systems in terms of the interfacial free energies. Furthermore, as shown in Section S2 of the Supplementary Materials, for the longer chain lengths considered here, one must apply special care to converge the SCF efficiently; the density profiles tend to shift further towards the bulk region with increasing N m and, as a result, the SCF is prone to entrapment in metastable states for a large number of iterations.
The thermodynamics of these films can be understood better in terms of the four macroscopic wetting functions, i.e., the work of cohesion (W C ), the work of adhesion (W A ), the work of spreading (W S ) and the work of immersion (W I ): where σ VM ≡ γ VM is the surface tension, (σ SVσ SM ) ≡ -γ SM is the adhesion tension [62], and θ the contact angle of the corresponding solid-fluid-vapor interface. Table 2 illustrates the wetting functions and contact angles of the PS/SiO 2 interphases studied here for N m = 384; they remain practically the same with larger chain lengths as shown in Figure S4. W A corresponds to the reversible work required to separate two phases in contact. It is noted that in the absence of the ramp potential, W A is significantly lower than in the HW and PW cases, where v ramp has been fitted to experimentally reported values of W A . W S quantifies the spontaneity of the wetting process: positive values indicate spontaneous spreading across the interphase (perfect wetting), while negative values indicate finite contact angles (partial or no wetting). In the LW and HW interfaces, W S remains negative across the full molar mass regime investigated here, indicating that the corresponding solid-fluid-vapor interphase will form finite contact angles. The PW interface, on the other hand, exhibits positive W S ; thus, PS will spread spontaneously on the silica surface.

EoS
Wetting   The behavior of grafted chains can be classified into three distinct regimes depending on the combinations of σ g and N g : Mushroom regime (low σ g N g ). The density of the profiles is less than the bulk density and chains assume randomly coiled conformations [63]. Increasing σ g and N g has a minor effect on the thickness of the profiles but rather makes them more pronounced.

2.
Dense brush regime (high σ g N g ). The brushes become stretched [46,63] and feature extended regions with bulk density. The thickness of the profiles depends strongly on both σ g and N g , reaching the limiting scaling behavior~σ g 1 N g 1 ; i.e., the dimensions of the brushes become proportional to their mass.

3.
Crowding regime (very high σ g ). In this regime, the crowding of the chains becomes so intense that the density of the grafted chains surpasses slightly the bulk one, as observed in our previous work [46]. This happens because the entropic penalty due to stretching overcomes the enthalpic penalty due to deviations from the bulk density.
There are, however, noticeable differences between the two systems. In GV, the thickness of the density profiles becomes commensurate with the number of the grafted PS segments which equals the product σ g N g . Indeed, the profiles in Figure 3a collapse together for constant σ g N g values, i.e., for constant amount of material. In GV, the tails of the profiles feature a sigmoid region on the order of 1 nm at the polymer/vacuum interface [44,61], whereas in GM they are much more expanded across the bulk region [15,46]. Finally, in both GV and GM, the profiles become slightly more pronounced with increasing σ g in the vicinity of the grafting points (~0.4 nm).
The dimensions of the grafted chains can be quantified in terms of their root mean squared brush thickness, which is calculated via the following Equation (29).
Additionally, we employ another measure for the extension of grafted chains, namely the height, h 99% , which is defined here as the distance between the surface of the solid (∂R solid ) and a surface (∂R h 99% ) that encloses 99% of grafted chain segments. It is determined via the following Equation (30).
with R 99% denoting the region between ∂R solid and ∂R h 99% . Figure 4a,b illustrates evaluations of h g 2 1/2 and h 99% against the scaling law~σ g N g in GV (left column) and GM with the matrix chains either being equal in length to the grafted ones (N m = N g ; central column), or varying between 24 and 1536 segments (right column). Figure 4c depicts the ratio h 99% / h g 2 1/2 which can be thought of as a measure of the shape of the profile. A striking difference between the two systems is that in the first one, the measures of the brush thickness collapse to a single master curve across the full regime. This is because in GV, the shape of the collapsed films is a function of the mass of the film (~σ g N g ) and does not depend on the ratio (σ g /N g ); see Figure 3a.
Across the mushroom regime, the brush thickness is practically independent of σ g . The thickness in vacuum is independent of N g as well, indicating collapse of the sparsely grafted chains on the surface. In the case where melt chains are present, the thickness exhibits a random walk-like N g -dependence, scaling approximately as~N g 0.5 . The shape of the brushes as quantified by the ratio h 99% / h g 2 1/2 , is quite sensitive to the grafting density in the presence of matrix chains; in GV, the ratio decreases with respect to the predictions from Alexander's model [64,65] for incompressible brushes, whilst the opposite behavior is exhibited in the GM system.
Across the dense brush regime the brush dimensions depend strongly on both σ g and N g . For very large σ g and N g the brush thickness scales as~σ g N g , no matter the solvent conditions (vacuum or melt). This means that the dimensions of the brushes become proportional to their mass and the ratio h 99% / h g 2 1/2 approaches the limiting value of √ 3 predicted by Alexander's model [64,65]. Regarding the effect of the matrix chains on the scaling of the brushes, it appears that, as long as N g ≤ N m , the brush dimensions are practically independent of N m . For N m N g , on the other hand, the brushes expand with decreasing N m due to the fact that the matrix chains can easily penetrate the brushes, thus the latter swell towards the bulk region. A similar behavior has been recently observed by Biltchak et al. [39]. Therefore, modulating N m allows for the tuning of the solvent conditions, from theta solvent (N m = N g ) up to good solvent (N m < N g ) conditions. Regarding the ratio h 99% / h g 2 1/2 , it features a complicated behavior with varying N m , where it decreases/increases for low/high σ g .
It seems that, regardless of the choice of the free energy density equation (e.g., compare lines (HFD) with markers (SL-SGT) in Figure 4) or the strength of the polymer-solid interactions (see Figure S5), the structural features and the scaling behaviors of GM systems are quantitatively very similar. In addition, the shortcoming of using HFD with constant κ T does not have a practical effect on the structural properties of grafted chains and on the potential of mean force of the system, as will be shown below. Figure 5 illustrates the total grand potential per unit area of grafted/vacuum (γ GV ) and grafted/melt (γ GM ) interphases as a function of N g and σ g , as well as the partial contributions from the polymer/solid interactions (γ s ), and the entropy of matrix (γ m ) and grafted chains (γ g ). Some key remarks regarding the evaluation of each term with SL-SGT can be summarized as follows: • γ s is purely of enthalpic origin and thus it is a functional of the total density profiles. γ s becomes more attractive with increasing σ g , since the profiles become more pro-nounced in the vicinity of the solid (e.g., see Figure 3). In GV, γ GV s decreases with increasing N g , since increasing N g increases the amount of material near the solid. In GM,γ GM s is independent of N g , since the total density profiles are also invariant to N g , e.g., compare profiles in Figure 2c-f for different chain lengths. • γ m describes the entropic contribution of the matrix chains. In GM, it decreases precipitously with increasing N g /N m ratio, since the grafted chains occupy more space in the vicinity of the Interphase leaving the matrix chains with fewer available conformations. Note that, by keeping N g fixed, γ m scales with N m about as γ m~Nm −1 ; see Equation (8). In GV systems, this term is of course zero since there are no matrix chains at all. • γ g quantifies the entropic contribution of the grafted chains and it is an indicator of the stretching of the brush [46]. It increases monotonically with increasing N g and σ g , since the grafted chains expand and stretch further towards the bulk region. In addition, γ GM g increases with decreasing N m , since the grafted chain-melt interactions are enhanced and, as a result, the brushes swell as shown in Figure 4b,c (rightmost column). • γ is the total free energy of the interfacial systems per unit area. It increases with increasing molecular weight of grafted chains and appears to be dominated by the conformational entropy term of grafted chains, γ g . Overall, the evaluations of the free energy terms with HFD are in good qualitative agreement with SL-SGT. Furthermore, it appears that evaluating γ GM s with HFD is more negative than when the SL-SGT EoS is used, because the density profiles lie closer to the silica surface in the former model (e.g., compare Figure 2d,f). In their recent work concerning neat grafted nanoparticles, Mydia et al. [40] report that for constant grafting density, the stretching energy does not increase monotonically with the chain length. This behavior is attributed to the curved space around the nanoparticles since, at some point, the grafted chains do not experience the presence of each other and become unperturbed. For planar surfaces, however, the threshold chain length becomes infinitely large, since no curvature is involved and the chains will always experience the presence of each other, considering that the dimensionless quantity σ g R g 2 is above a threshold value as well [46].

PS Melt Capped between Two Bare Silica Surfaces (SMS)
It is worth studying the potential of mean force (PMF) between two approaching bare silica surfaces. This means that only matrix chains are present in the system. This situation corresponds to the limiting case of very low grafting densities, where allophobic dewetting occurs and the enthalpic interactions between the two solid surfaces prevail. The matrix chains are gradually restricted in terms of available conformations as the distance between the two plates decreases. For these calculations, the PMF SMS is expressed with respect to the free energies of the single SM interphases for the same chain length: with γ SM being the free energy of a SM interphase in presence of matrix chains of length N m , depicted in Figure 2a,b. Note that, with PMF SMS known, the disjoining pressure can be calculated as: Due to numerical stability issues, the initial configuration of the field can affect the outcome of the converged solution. To investigate this effect, we performed these calculations using two different compression methods:

1.
In the first method, the calculations were performed in a decremental fashion, in which the initial configuration of the field was set to the field corresponding to the converged calculation for a slightly larger domain, w ifc,init (h ss ) = w ifc,final (h ss + ∆h) 2.
In the calculations corresponding to the second method, the initial configuration of the field was set to zero across the domain, w ifc,init = 0. Using the first method, it is easier to derive a solution that corresponds to a stable configuration. The second method, on the other hand, can provide a measure of the stability of the films in terms of their tendency to collapse and their sensitivity to fluctuations about equilibrium (e.g., their response during the formation of a cavity). Figure 6 illustrates evaluations of PMF SMS with HFD (a), and with SL-SGT under low (b, LW), high (c, HW) and perfect (d, PW) wetting situations, using the first method for decremental compression. Results concerning the second method are presented in Supplementary Materials Section S5, along with a relevant discussion.
In the case of HFD, regardless of the compression approach, PMF SMS increases with decreasing plate-plate distance, suggesting the manifestation of a repulsive force that resists the attractive interactions between the surfaces. These repulsive forces are dominated by the loss of polymer-solid interactions; with decreasing h ss , the mass of the film decreases and there are fewer interaction sites. The sign of these forces depends on an interplay between the strength of the polymer-solid and the solid-solid interactions. If the latter become much stronger than the former, the solid-solid forces will dominate and PMF SMS will become attractive. The steadily increasing forces in this case can be, however, a misleading result, because HFD does not account for the gaseous phases which may arise during the process. According to evaluations with SL-SGT under low wetting (LW) conditions (Figure 6b), the functional dependence of PMF SMS is quite similar to that of HFD. However, below a critical plate-plate distance (h crit ss ), PMF SMS decreases abruptly, indicating a phase transition. At these distances it is impossible for SCF to maintain a metastable film; hence, a cavity is formed and the calculation converges to the more stable solution, ϕ = 0. The films remain stable above h crit ss on the order of 3.5 nm, regardless of N m (Figure 6b). Upon departure of the melt from the gap between the plates, the only contribution to the free energy is due to plate-plate interactions which are described here by means of the Hamaker potential; therefore, leading to the eventual contact of the adjacent solid surfaces. Note that, for the LW surfaces, PMF SMS (ϕ = 0) = γ SMS ss − γ SMS inf , which is depicted by dashed lines in Figure 6b is negative, indicating that the polymer-solid interactions are really weak for the LW films and these films are actually metastable with respect to cavitation across the entire range of thicknesses.
A similar picture has been reported by past simulations from a variable-density lattice based SCF model [66]. In that model, the interactions in flat geometries become insignificant for h ss slightly larger or equal than 4 R g 2 0.5 , whereas the maximum recorded force per radius before cavitation in a crossed cylinder geometry was on the order of~0.01 mN/m when considering high energy surfaces and 0.1 mN/m for low energy surfaces. With enhanced polymer/solid interactions, the stability of the capped polymer films is reinforced considerably. According to the more reliable solution scheme where the domain length is adjusted decrementally, the HW and PW films remain stable throughout the full h ss range examined here, always converging to the stable, polymer-filled equilibrium solution of the problem (see Figure 6c,d). In contrast to LW films, the PMF SMS of the HW and PW ones increases steadily for h ss less than 2.5 nm, whereas in the limit of low h ss the density decreases significantly due to the entropy of confinement. A significant free energy barrier of approximately −γ SMS inf , on the order of 20 mN/m in the HW case and 80 mN/m in the PW case, has to be overcome for the polymer to be expelled completely and the solids to come in direct contact at PMF SMS (ϕ = 0) = γ SMS ss − γ SMS inf . There is, however, a striking difference between the HW and PW films. The HW films can be considered as metastable with respect to cavitation, since, after crossing a barrier of~−γ SMS inf with decreasing h ss , they can be trapped in the global minimum of attractive Hamaker interactions between the bare solids (see minima of the dashed lines in 6c), as in the case of the LW film. The PW film, on the other hand, is indeed stable down to thicknesses of ca. 1.5 nm, since the minimum of the Hamaker potential lies way higher (~+60 mJ/m 2 ) than the plateau PMF SMS value at large h ss .

Interacting Grafted Surfaces in Melt (GMG)
The current section presents evaluations of the potential of mean force of approaching grafted surfaces. To facilitate comparisons across the wide parameter space considered in this work, the PMF GMG will be expressed in terms of the reduced surface-surface distance which is defined by the following Equation (33).
with h g − 2 0.5 and h g + 2 0.5 being the root mean squared brush thickness of the single grafted surfaces (infinite surface-surface distance) at the same temperature and in the presence of matrix chains of length N m (brush thickness from Figure 4). A similar normalization can be obtained by dividing h ss with R g 2 0.5 , since R g 2 0.5 ∼ N g and N g ∼ h g 2 0.5 across the dense brush regime. However, normalizing the plate-plate separation distances as shown in Equation (33) is a more natural approach of making such comparisons, since it allows evaluating the tendencies of the brushes to interpenetrate. In addition, h g 2 0.5 takes account of the chain perturbations when varying the melt conditions and leads to comparisons that are less sensitive to the particular equation of state used in the nonbonded free energy density model.
As in the SMS systems, PMF GMG will be expressed relative to the free energy of the isolated G − M and GM + systems in the presence of matrix chains of length N m (see Figure 5).
In other words, PMF GMG is expressed with respect to the free energy of the system at infinite plate-plate separation. Regarding the segment/substrate interactions-unless otherwise stated-LW conditions are employed throughout this section.

Symmetric Surfaces
We present the PMF GMG for the simplest case where the opposing surfaces are grafted symmetrically with respect to σ g and N g for varying N m /N g . Figure 7 illustrates PMF GMG as a function of the plate-plate distance, h ss , whereas Figure 8 depicts the corresponding density distribution across the examined parameter space for N m = N g . Figure 7. Potential of mean force against the reduced surface-surface distance in a symmetric system of approaching grafted silica surfaces in a melt from the SL EoS. Colors correspond to evaluations for N m /N g = 1/2 (red), 1 (blue) and 2 (green), whereas the labels in brackets denote σ g (nm −2 ) and N g . Lines are guides to the eye.  Figure 7, as a function of the plate-plate distance (abscissa) and the distance from the left wall (ordinate) in reduced units. Red/green/blue colors correspond to regions with high density (ϕ c = 1) of c = g − /m/g + chains, respectively. Grey denotes regions which lie outside the modeled domain or have not been evaluated at all. Labels in brackets denote σ g (nm −2 ) and N g .
According to Figure 7, the brushes start to experience the presence of each other at distances in the order of 4-5 h ss , while for larger h ss , PMF GMG 0. At lower separation distances the brushes interact strongly with each other and PMF GMG increases. In systems with low σ g (mushroom regime), the brushes are relatively soft and thus PMF increases at lower h ss . In dense systems, on the other hand, the brushes are more compact and as a result the PMF increases abruptly at larger h ss , on the order of 3-3.5 (see bottom right panel of Figure 7). This regime coheres with the predictions of Alexander's model for incompressible brushes [64,65], where, regardless of N g and σ g , the predicted separation distance of two Alexander brushes in contact is h ss,min = 2 √ 3 ∼ 3.4; the calculation is illustrated in Appendix A.
Regarding the effect of varying the chain length of matrix chains on PMF GMG , for N m /N g ≤ 1, PMF GMG becomes strictly repulsive, with the exception of the case for densely grafted and long grafted chains, σ g = 0.4 nm 2 , N g = 768. Increasing N m /N g leads to the manifestation of attractive interactions as indicated by the formation of a minimum in PMF GMG (autophobic dewetting). These interactions become slightly stronger with increasing N g , and even stronger with increasing σ g . The enhancement of the attractive interactions with increasing grafting density has been also observed in several theoretical [22,25,67] and experimental studies [68][69][70][71][72]. Figure 9 illustrates the contribution of individual free energy terms to PMF GMG for a case with σ g = 0.2 nm -2 and N g = 192. It appears that the attractive part of PMF GMG is dominated by the entropic contribution of the grafted chains, γ GMG g , shown in Figure 9c. In addition, the terms responsible for the cohesive and the chemical potential-density field interactions change due to minor variations in the mean density profile and exhibit an opposite trend to that of γ GMG g , albeit weaker. That PMF GMG becomes more repulsive when N m < N g is attributed to the matrix chains being able to penetrate the space occupied by the brush, compelling the grafted chains of the two surfaces to expand in space until they interact with each other, thus keeping the two surfaces separated. In other words, the solvent conditions improve with decreasing N m /N g ; hence, the brushes prefer to interact with the solvent molecules than with each other. In general, the evolution of the term associated with the polymer-solid interactions, γ GMG s , with decreasing plate-plate distance, depends on an interplay between two processes: (i) The thinner the polymer film becomes, the fewer polymer segments can interact with the solid, thus leading to increased γ GMG s (less attractive). (ii) The brushes can become slightly denser with decreasing h ss , hence leading to decreased γ GMG s (more attractive).
Indeed, for N m > N g , γ GMG s is dominated by the first process, since it increases monotonically with decreasing h ss . However, for N m < N g it features an interesting behavior where it initially increases with decreasing h ss , while for h ss < 3.8 it decreases, indicating that the second process dominates.

Asymmetric Surfaces
The present section investigates the effect of asymmetry of the opposing grafted surfaces in terms of the relative chain lengths and grafting densities. In particular, Figure 10 illustrates evaluations of PMF GMG for constant chain length and grafting density at the lower face (N g − = 96, σ g − = 0.2 nm −2 ) and varying N g + /N g − , N m /N g − and σ g + /σ g − . Similarly, Figure S7 and Figure S8 depict the same evaluations but for four times larger N g − and for two times larger σ g − , respectively. Figure 10. Potential of mean force against the reduced surface-surface distance for approaching grafted silica surfaces with σ g − = 0.2 nm −2 and N g − = 96. Colors correspond to evaluations for N m /N g − = 1/2 (red), 1 (blue), 2 (green) and 4 (purple), and the labels in brackets denote the ratios σ g + /σ g − and N g + /N g − . Lines and markers correspond to evaluations with the HFD and SL EoS, respectively. Inset: each bead denotes a segment of 24 PS monomers.
Irrespectively of the degree of asymmetry, adjacent brushes experience the presence of each other at distances commensurate with 4-5 h ss , similarly to the symmetric case in Figure 7. By keeping σ g − , N g − constant, the ratio N m /N g − fixed at various values (same colors), and varying σ g + and N g + , some general trends are emerging: the attractive interactions between the surfaces become stronger with increasing σ g + /σ g − ratios (from top to bottom) and with decreasing N g + /N g − (from right to left).
That PMF GMG becomes more attractive with increasing σ g + /σ g − is to be expected: upon increasing σ g + the mean grafting density increases as well, thus PMF GMG becomes more attractive, as in the case for the symmetric surfaces in Figure 7. However, that the interactions become more attractive with decreasing N g + /N g − has to be reconciled with the findings reported in Figure 7 for symmetrically grafted surfaces. To interpret this effect, one should take into account that during these evaluations the ratio N m /N g − was fixed, whereas varying σ g + and N g + can have direct implications for the effective ratio of N m with respect to the average size of grafted chains.
Let N g be the average chain length of grafted chains, estimated as the weighted average of grafted chain length with respect to the grafting densities: N m /N g is a measure of the length of matrix chains in relation to grafted chains. Based on Equation (35) it is evident that, upon decreasing N g + /N g − at constant N g − and σ g + /σ g − , N g decreases as well and thus the effective ratio N m /N g increases. Increasing the size of matrix chains was shown to enhance the attractive interactions in Figure 7, and the results shown in Figure 10 are consistent with this trend. An alternative way to interpret and isolate the effect of asymmetry is to vary the ratios σ g + /σ g − and N g + /N g − , but fix the effective ratio N m /N g based on Equation (35). Figure 11 depicts such evaluations, wherein the top/left panel depicts the reference symmetric case (σ g = 0.2 nm -2 , N g = 96), while the asymmetry with respect to σ g (N g ) increases from topto-bottom (left-to-right). As can be seen, varying σ g + /σ g − or N g + /N g − individually has a minor effect on PMF GMG . This finding is important, since it shows that minor deviations in σ g and N g do not affect PMF GMG provided N m /N g is fixed. On the other hand, PMF GMG can become very attractive in extreme cases of asymmetry where both the asymmetry in σ g and N g increases (e.g., see bottom right panel in Figure 11), and this is mainly attributed to the increased average grafting density. Figure 11. Potential of mean force against the reduced surface-surface distance for approaching grafted silica surfaces. Colors correspond to evaluations for N m /N g = 1/2 (red), 1 (blue), 2 (green) and 4 (purple), whereas the small legends in brackets below the insets denote the ratios σ g + /σ g − and N g + /N g − , respectively. The top left panel corresponds to the reference symmetric case with σ ref g = σ g − = σ g + = 0.2 nm -2 and N ref g = N g − = N g + = 96. The actual N g ± and σ g ± for each case can be retrieved as follows: N g ± = N ref g n ±1/2 and σ g ± = σ ref g m ±1/2 , with [n, m] = [N g + /N g + ,σ g + /σ g − ] being the numbers at the legends under the insets. It should be noted that inset schematics belonging to panels other than the corner ones are only approximate; in these cases, σ g +/− and N g +/− are scaled by a factor of ± √ 2. Figure 12 presents the PMF GMG for symmetric systems as a function of the energy of the solid surface. It seems that, regardless of σ g , N g and the ratio N m /N g , the strength of the polymer-solid interactions (wetting situation) has a minor effect on PMF GMG . The PMFs become less pronounced with enhanced wetting situation, albeit the effect is minor. Figure 12. Potential of mean force against the reduced surface-surface distance of approaching symmetric grafted silica surfaces in a melt from the SL EoS for σ g set to (a) 0.1, (b) 0.2 and (c) 0.4 nm -2 and N g = 192. Colors correspond to evaluations for N m /N g = 1/2 (red), 1 (blue) and 2 (green), whereas different styles denote interfaces with low (solid lines), high (dashed lines) and perfect (dots) wetting.

Interacting Grafted Surfaces in Vacuum (GVG)
The present section discusses evaluations of the PMF GVG of approaching grafted surfaces separated by vacuum (GVG). Figure 13 presents evaluations of PMF GVG with N g ranging from 24 to 768 and σ g from 0.1 to 0.4 nm -2 , whereas Figure 14 illustrates the corresponding density distributions for these cases. LW conditions have been used for surface-segment interactions.   Figure 13, as functions of the plate-plate distance (abscissa) and the distance from the left wall (ordinate) in reduced units. Red/blue colors correspond to regions with high density (ϕ c = 1) of c = g -, g + chains. White corresponds to vacuum. Grey denotes regions which lie outside the modeled domain or have not been evaluated at all. Labels in brackets denote σ g (nm −2 ) and N g .
The evolution of PMF GMG with decreasing plate-plate distance can be classified in three distinct regimes:

1.
For large separation distances the adjacent polymer brushes interact weakly with each other and the dominant contribution to PMF GVG arises due to polymer/solid and solid/solid Hamaker interactions (e.g., compare with the dotted lines in Figure 6b-d).

2.
Below a critical plate-plate distance, PMF GVG decreases abruptly, indicating the manifestation of a phase transition where the adjacent brushes interpenetrate each other and form a single film in the central region of the system. In addition, low density regions are formed in the vicinity of the solid surface, indicating that the brushes have been stretched significantly towards the bulk region.

3.
Decreasing the plate-plate separation further makes the brushes more compact. The low-density regions next to the solid become suppressed, until the free energy becomes commensurate to minus the mean surface free energy of the isolated brushes, −0.5 γ G − V + γ VG + . Further squeezing of the brushes leads to increased reduced densities above unity, as indicated by the vertical dashed lines in Figure 13.
The chain configurations across these regimes are illustrated schematically in the inset of Figure 15f. A more detailed picture can be unveiled by inspecting the evolution of individual contributions to the energy terms shown in Figure 15. The thermodynamics of the merger seems to be dominated by cohesive interactions. According to Figure 15a, below some critical distance the abrupt drop of the cohesive term (γ GVG coh ) indicates the enthalpic gain upon film merging. At the same time, the more positive polymer-solid interactions in Figure 15b (γ GVG s ) indicate the enthalpic penalty due to the departure of a large portion of the brushes from the solid surface. The term associated with the entropy of the grafted chains (γ GVG g ) in Figure 15c is of particular interest. At first glance, it does not reflect the entropic penalty due the stretching of the grafted chains. However, this is attributable to the fact that γ GVG g does not reflect the total conformational contribution to the grand potential, since it is evaluated in the presence of the field [46]. The conformational component of the grafted chains can be retrieved as follows [46]: with γ GVG g,field being the field experienced by the grafted chains: Indeed, as indicated in Figure 15e, the conformational free energy of the grafted chains increases abruptly below a critical distance (entropic penalty due to stretching) and then decreases with decreasing h ss , as the film becomes more compact and the grafted chains become less stretched.
In summary, the manifestation of the phase transition depends on an interplay among three dominant factors: • An enthalpic gain due to the lower surface area of the merged brushes that increases with increasing surface tension.
• An enthalpic loss due to the detachment of the grafted film from the surface that depends on the polymer/solid interactions. • A conformational penalty due to chain stretching.
As far as the equilibrium plate-to-plate distance after the merger is concerned, for low σ g N g products this is on the order of 1.5 h ss , or 0.75(h g + + h g − ). In this case the plates come considerably closer than the sum of the individual root mean squared brush thicknesses, since the brushes lie in the mushroom regime and can interpenetrate each other easily.
In the limit of large σ g N g , however, the brushes are more compact, and therefore one could make meaningful predictions using Alexander's model for incompressible brushes (see Appendix A) [64,65]. Indeed, the denser brushes investigated here become compact (ρ > ρ seg,bulk ) at separation distances on the order of h ss,min ∼ 2 √ 3 ∼ 3.4, e.g., compare with the vertical line in the bottom-right panel of Figure 13.
Overall, the effect of asymmetry on the equilibrium distance is expected to be minor. Indeed, according to Figure 13, brushes with equal N g σ g products become compact at similar distances. For example, compare the case (σ g ,N g ) = (0.4 nm -2 , 48) with (0.1 nm -2 , 192) and the case (0.4 nm -2 , 192) with (0.1 nm -2 , 768).

Discussion
All of the above results concerning the potential of mean force between two planar surfaces (or equivalently of particles with radius large in relation to the chain dimensions) imply some design rules that one can be guided by when addressing such nanostructured systems.
As mentioned before, addressing the system of bare surfaces can give us an insight on what happens in the limit of very low grafting density. At moderate distances the behavior of the PMF depends on an interplay between the strength of the solid-polymer, polymerpolymer, and solid-solid interactions, as quantified by the wetting/spreading phenomena taking place at a single solid/polymer interface. For low wetting conditions (θ > 90 • ), the PMF becomes weakly repulsive from distances of the order of 3.5-10 nm, but the polymer film is metastable with respect to the chains evacuating the gap and the solid surfaces snapping into direct contact with each other. The spontaneous manifestation of cavities can lead to eventual collapse of the surfaces, leading to agglomeration. This is a manifestation of the allophobic dewetting phenomenon observed in the low wetting situation examined here (see SL-SGT (LW) case). On the contrary, for high wetting (θ < 90 • ) and perfect wetting (spreading) conditions, the PMF is practically zero at large distances and starts rising steeply below ca. 2 nm. A free energy barrier on the order of 20-80 mJ/m 2 has to be overcome for the solid surfaces to come into direct contact. The system is stabilized in relation to the solid surfaces sticking to each other, and becomes trapped in a potential well with a depth in the order of ca. 20 mJ/m 2 with respect to a melt-free system, e.g., with respect to the dotted lines shown in Figure 6. Matrix chains adhering to the solid surfaces resist compression and screen the solid-solid attractive interactions, as if they were grafted. The potential of mean force between bare solid surfaces does not depend strongly on the length of matrix chains; therefore, varying the molecular weight of matrix chains does not have a significant effect on the stability of the system.
In the case of low wetting of the solid by the polymer, grafted chains are necessary to stabilize the dispersion. In particular, when the grafted chain lengths and grafting densities are sufficiently high (e.g., as shown in Figure 12), the grafted chains effectively screen the solid-solid interactions, preventing the two plates from approaching each other at a distance where they would experience the full depth of the plate-plate potential. In many cases, PMF can become strictly repulsive, but a prerequisite for this is that the grafted chains be longer than the matrix chains. The longer the grafted chains in comparison to the matrix chains, the steeper the repulsive potential of mean force that develops and the longer the range over which it manifests itself. Short matrix chains are able to penetrate the brush of long grafted chains and swell it, increasing the range of the repulsive interaction. In this case, the most important design rule that has to be met, as already found experimentally [68][69][70][71][72], is that the grafted chains need to be larger the matrix chains, N m /N g ≤ 1. On the other hand, when matrix chains start becoming larger than the grafted ones, then immediately an attractive well is exhibited in the potential of mean force of the system (autophobic dewetting); solvent conditions become worse for the grafted chains. Furthermore, this behavior is intensified at high surface grafting densities. In these cases, the density of grafted chain segments near the interface is so high, that matrix chains are not able to penetrate into the region occupied by grafted chains, even if they have lower molecular weight than grafted chains. Thus, the use of excessive grafting densities is to be avoided for the purposes of steric stabilization, even if N g > N m .
For the asymmetric cases, where the grafting density or the grafted chain molecular weight of grafted chains differs on each silica plate, it seems that the introduction of asymmetries does not give rise to a minimum in the potential of mean force, as long as deviations from the symmetric system are small, and the effective ratio N m /N g from Equation (33) is fixed. Individually, adjusting the asymmetry on the grafting density or on the molecular weight of grafted chains does not alter significantly the potential of mean force, with the former having a stronger influence. This implies that, when experimentalists are trying to stabilize dispersion of the two slabs (or large grafted particles) in a melt, there is some room for deviation from symmetry, especially as regards the molecular weight of grafted chains. On the other hand, when there are larger discrepancies between both the grafting densities and the molecular weight ratios, then the system will eventually exhibit an attractive well. Again, this phenomenon is more pronounced when the molecular weight of matrix chains increases.
The characteristics of the well depth of the PMF GMG reflect the "softness" of the brushes as well as the associated tendency to penetrate into each other. [22] Wrapping together all the parameters which influence the attractive well of the PMF GMG between the two plates, one can generate empirical design rules regarding the prediction of stable configurations of opposing plates (membranes or fine particles) as a function of the mean grafting density, the chain length of grafted chains and the chain length of matrix chains. According to Hasegawa et al. [25], PMF is expected to become repulsive for σ g ≤ b k −2 N k, g −1/2 , with N k,c = N c C ∞ l c−c 2 /b k 2 being the number of Kuhn segments that comprise a type c chain. In addition, as it has been demonstrated in Section 4.3, PMF GMG becomes more attractive with increasing chain length of the matrix chains with respect to the effective length of grafted chains, N m /N g (or N k,m /N k, g in Kuhn units). Note that the effective ratio, N m /N g , takes account of the asymmetry in both σ g ± and N g ± . Putting all these together, it would be instructive to present the depth of the attractive well as a function of the dimensionless quantity, σ g b k 2 N k, g 1/2 N k,m /N k, g , or σ g b k 2 N k, g −1/2 N k,m , for simplicity. Such comparisons are shown in the master plot of Figure 16 against all the data gathered here for both symmetric and the asymmetric surfaces. According to Figure 16 the attractive interactions become negligible for σ g b k 2 N k, g −1/2 N k,m < 5, therefore for such combinations the surfaces are expected to stabilize. For larger values, on the other hand, in most cases the plates stick to each other (aggregation). Nevertheless, we must take into account that in the limit of very small average grafting density, the system will be led towards the case of bare solid plates in contact with vacuum, thus the melt is expected to evacuate the gap and the plates will spontaneously come in contact to each other. According to our calculations, the PS melt remains stable even for low grafting densities on the order of 0.1 nm −2 (∼ 0.33 b k −2 ), thus for these systems the region of stability could be traced along the range 0.33 ≤ σ g b k 2 ≤ 5 N k, g 1/2 N k,m −1 . Our findings conform with experimental studies concerning systems of the same [69][70][71] or similar [68,70,72] chemical constitution and with theoretical works [22,25,67], whilst accounting for the effect of asymmetry as well.

Conclusions
The present work investigates the thermodynamics of single and opposing silica surfaces, either bare or grafted, under various wetting conditions and solvent situations. The calculations have been realized with the self-consistent field theory. For comparison purposes, the nonbonded dispersive interactions are based on the Helfand (HFD), or the Sanchez-Lacombe (SL) free energy density. The polymer-solid and solid-solid interactions were described with the Hamaker potential in conjunction with a ramp potential. By adjusting the depth of the ramp potential, we examined cases with low wetting (θ = 140-160 • , depending on the EoS), high wetting (θ~67 • ), and perfect wetting (spreading), whereas for each case we derived the corresponding macroscopic wetting functions.
We investigated situations of grafted surfaces in contact with vacuum or a melt phase, in terms of their structure and interfacial free energies. In vacuum conditions, the overall shape of the profiles and the brush thickness are dictated by the product, σ g N g ; the brushes appear collapsed and feature a sigmoid region of ca. 1 nm at the vacuum interface. In melt conditions, the brushes expand significantly towards the bulk phase. Nonetheless, in the high σ g N g regime, the scaling of the brush thickness conforms to the prediction of Alexander's model for incompressible brushes, independent of whether they are in contact with vacuum or melt. The EoS and the wetting conditions appear to have a minor effect on the overall structure of the single brushes.
We performed a systematic and comprehensive investigation of PMF over a broad parameter space in terms of the following aspects: • Grafting densities; very low (σ g = 0, bare surfaces), moderate and high. • Asymmetry regarding the grafting densities and grafted chain lengths. • Solvent conditions; good (N m /N g < 1), theta (N m /N g = 1), bad (N m /N g < 1 ), and very bad (vacuum) solvents. • Wetting degree; low, high and perfectly wetted interfaces.
In doing so, we isolated the key parameters that are important to the resulting PMF and the tendencies of the opposing surfaces (or the corresponding fine particles) to agglomerate.
Initially, we evaluated the potential of mean force (PMF) between two bare silica plates with melt between them, each time varying the wetting situation. HFD was found inadequate to describe these systems because it cannot describe the occurring cavitation, whereas SL can. The LW films collapse for surface-surface distances lower than ca. 3.5 nm. In contrast, the HW and PW films are retained during the whole compression experiment. A difference between HW and PW is that the latter are stable with respect to cavitation down to thicknesses of ca. 1.5 nm, since the minimum due to Hamaker potential at low distances is located higher than the plateau free energy at large distances. In these situations, the matrix chains adhere strongly to the surfaces-as if they were grafted-and are able to effectively screen the solid-solid attractive interactions.
Regarding the PMF between grafted surfaces with melt in between them, a general conclusion is that with decreasing molecular weight of matrix chains, the solvent conditions for the grafted chains are enhanced, meaning that the grafted chains prefer to interact with the chains of the melt than with themselves. This behavior is manifested through the repulsive potential of mean force between the two grafted solid plates. On the other hand, when the melt chains are of much larger molecular weight than the grafted chains or when vacuum exists between the two grafted plates, then the grafted chains prefer to interact with themselves, adopting the corresponding configurations. The effect of asymmetry is minor when adjusting the lengths and grafted chains individually, and given that the effective ratio N m /N g is fixed, indicating that there is some room for error during the experimental processes. The role of wetting and the employed EoS appear to have a minor effect on PMF. Based on our calculations we find a region of stability (steric stabilization) traced along the following range:0.33 ≤ σ g b k 2 ≤ 5 N k, g 1/2 N k,m −1 . This predictive behavior is in agreement with relevant studies reported in literature [22,25,[67][68][69][70][71][72].
Regarding the grafted silica plates without melt chains in the space between them, it would be interesting to confirm experimentally the predictions derived here on how their equilibrium distance is affected by the grafting density and the matrix to grafted chains molecular weight ratio. For the lowest grafting densities investigated here the equilibrium distances are on the order of 1.5 the mean brush thickness, whereas with increasing N g and σ g the brushes become more compact and the equilibrium distance increases to about 2 √ 3 mean brush thickness.
Future prospects of the study include the investigation of the potential of mean force between opposing surfaces (or even pairs of nanoparticles) embedded in melt or vacuum phases (particle solids) using the three-dimensional version of RuSseL [52]. These calculations will allow us to investigate situations with homogeneous, and inhomogeneous distributions of grafted points across the surface of the NP (plane).

Conflicts of Interest:
The authors declare no conflict of interest.