Advances in Tribology

Volume 2013 (2013), Article ID 794151, 27 pages

http://dx.doi.org/10.1155/2013/794151

Review Article

A Description of Multiscale Modeling for the Head-Disk Interface Focusing on Bottom-Level Lubricant and Carbon Overcoat Models

Myung S. Jhon,1,2 Pil Seung Chung,¹ Robert L. Smith,¹ and Lorenz T. Biegler¹

¹Department of Chemical Engineering and Data Storage Systems Center, Carnegie Mellon University, Pittsburgh, PA 15213, USA

²School of Advanced Materials Science and Engineering, Sungkyunkwan University, Suwon 440-330, Republic of Korea

Received 1 October 2012; Accepted 19 December 2012

Academic Editor: Bruno Marchon

Copyright © 2013 Myung S. Jhon et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

The challenges in designing future head disk interface (HDI) demand efficient theoretical modeling tools with flexibility in investigating various combinations of perfluoropolyether (PFPE) and carbon overcoat (COC) materials. For broad range of time and length scales, we developed multiscale/multiphysical modeling approach, which can bring paradigm-shifting improvements in advanced HDI design. In this paper, we introduce our multiscale modeling methodology with an effective strategic framework for the HDI system. Our multiscale methodology in this paper adopts a bottom to top approach beginning with the high-resolution modeling, which describes the intramolecular/intermolecular PFPE-COC degrees of freedom governing the functional oligomeric molecular conformations on the carbon surfaces. By introducing methodology for integrating atomistic/molecular/mesoscale levels via coarse-graining procedures, we investigated static and dynamic properties of PFPE-COC combinations with various molecular architectures. By bridging the atomistic and molecular scales, we are able to systematically incorporate first-principle physics into molecular models, thereby demonstrating a pathway for designing materials based on molecular architecture. We also discussed future materials (e.g., graphene for COC, star-like PFPEs) and systems (e.g., heat-assisted magnetic recording (HAMR)) with higher scale modeling methodology, which enables the incorporation of molecular/mesoscale information into the continuum scale models.

1. Introduction

The continuous increase in the areal recording density specification beyond 1 Tb/in² has led to ever decreasing head media spacing (HMS) requirements at the head disk interface (HDI). The key material components of the HDI are the carbon overcoat (COC) and lubricant layers, which protect the magnetic media from corrosion and tribological damage. Perfluoropolyethers (PFPE) with both functional and nonfunctional groups are standard HDI lubricants due to their low vapor pressure and low surface tension as well as good chemical and thermal stability. To make a more reliable product, improved lubricant and COC materials must have self-healing capability and lubricant-COC adhesion in addition to molecularly thin spreading layer thickness.

The challenges involved in designing improved HDI materials require efficient theoretical modeling tools which allow flexibility in investigating various pairs of PFPE-COC materials. Due to the broad range of time and length scales of interest in the HDI components, a multiscale/multi-physical modeling approach can be utilized to produce paradigm-shifting improvements in advanced HDI design [1–13]. We believe that the multiscale modeling methodology as described in this paper is an effective strategy for this system, although it is far from being complete. The essence of multiscale modeling involves integrating low- and high-resolution models of a system so that information is passed seamlessly between adjacent scales. This can be accomplished via either a top to bottom, bottom to top, or middle-out approach supplemented by optimization as illustrated in Figure 1. In the top to bottom approach, low-resolution descriptions are established, and inputs from lower-scale descriptions are incorporated to pass information to the upper-scale model. An advantage of this approach is that it is often more practical to construct a model closer to the device/macroscopic scale first so that the device performance determines the lower scale information that is needed. However, it may be more difficult to make fundamentally significant changes in material design in the top to bottom procedure with the present modeling paradigm. In the inverse approach (bottom to top), the detailed, high-resolution descriptions are constructed first, and via coarse-graining procedures, the outputs from the model are inputs for the development of lower resolution descriptions. This multiscale integration style is particularly useful in impacting device performance via nanoscale materials science and engineering, since the molecular architecture is at the root of material properties. A unique middle-out approach employs an intermediate scale model and outputs information to models at adjacent scales above and below the primary model allowing for a flexible integration strategy.

Figure 1: An illustration of the several multiscale modeling methodologies.

For the HDI system, constructing a multiscale description involves establishing an integrated model for both lubricant and overcoat materials. In this paper, we adopt a bottom to top approach to outline our methodology, although top to bottom or middle-out approaches with optimization are also suitable for HDI design. Due to the polymeric structure of the PFPE system, we begin at the high-resolution, fully atomistic scale with the construction of models which describe the intramolecular PFPE degrees of freedom governing the chain conformations. We transition to intermolecular interactions between PFPEs as well as between PFPE-COC combinations to predict lubricant self-healing (mobility) and adhesive behavior. We then discuss methods for integrating the atomistic models with molecular/mesoscale models via coarse-graining procedures. The coarse-grained descriptions allow us to approach the time and length scales most relevant for HDI operation. By bridging the two scales, we are able to systematically incorporate first-principle physics into molecular models, thereby demonstrating a pathway for designing materials based on molecular architecture. This approach allows us to describe the conformations, dynamics, morphology, rheology, and thermal performance of HDI materials.

We begin the paper with an overview of intramolecular and intermolecular studies of PFPE and COC materials at the fully atomistic scale where ab initio methods have been utilized. The atomistic review is followed by a survey of works investigating the HDI at the molecular/mesoscale via molecular dynamics simulations. This is followed by an overview of novel integration methods and a discussion of the perceived impact of multiscale modeling on the challenges facing advanced HDI design.

2. Atomistic Scale Models

The atomistic scale models include the most physically detailed description necessary for linking material performance with molecular architecture. Investigations at the atomistic scale have included the intramolecular and intermolecular PFPE interactions as well as PFPE-COC interactions. The intramolecular study is motivated by the need to describe the conformational behavior of PFPEs. In the study reviewed here, Z-dol and Z-tetraol serve as examples of this method, although similar treatment can be applied to other lubricants. This is followed by a PFPE-PFPE interaction energy study, since the strength of the intermolecular interactions influences the dynamic behavior of the molecules and thereby, the mobility and self-healing capability of the PFPE in the event of head-disk contact or lubricant pick-up by the head. To determine the PFPE-COC combinations which promote strong adhesion, an interaction energy study for various PFPEs on COCs is given.

2.1. Intramolecular PFPE Potential Energy

The intramolecular study [14] provided bond stretching, angle bending, and torsional potential energy parameters for PFPE molecules for use in atomistic MD simulations. In this paper, we start with a detailed ab initio treatment of the functional PFPE (Z-dol and Z-tetraol as benchmark examples) intramolecular force field as the foundation for larger scale models so that the fundamental physics of the microscopic scale is linked with the macroscopic. Here Z-dol and Z-tetraol are chosen for this analysis due to their prevalence in PFPE research; however, this method will be extended to other lubricants with various functionalities and molecular architecture. In order to investigate the oligomers with computationally intensive ab initio methods, model Z-dol (HOCH₂CF₂OCF₂CF₃) and Z-tetraol (HOCH₂HOCHCH₂OCH₂CF₂OCF₂CF₃) molecules are utilized which preserve the end-group structure.

The relaxed structures are calculated via geometry optimizations, followed by the intramolecular force field parameters using the method of Seminario [15]. We developed a code which uses LAPACK [16] routines for eigenvalue analysis of the ab initio second-derivative tensor (Hessian) and calculates stretching and bending parameters for harmonic potentials:

where and are the stretching and bending potentials, respectively, is the stretching force constant, is the bending force constant, is the equilibrium bond length, and is the equilibrium angle. The stretching interactions take place between two bonded atoms, while the bending accounts for the bending of the angle formed by two adjacent bonds. The torsional potential parameters, describing four bonded atoms, were calculated by generating the torsional energy profiles via a series of constrained geometry optimizations. The form of the torsional potential is

where is the torsional potential, are the torsional potential coefficients, is the truncation order, and is the torsional angle. These potential energy forms were selected due to their widespread use in popular MD codes.

Selected equilibrium bond lengths and angles and the corresponding Hessian-derived force constants are displayed in Table 1 for Z-dol. The study reveals the stiffness of short polar bonds such as H–O as well as the higher frequency of bond stretching modes in comparison to angle bending modes. In general, the bond and angle force constants we calculate are lower than those given in a force field such as the General AMBER Force Field [17]. These differences are due to the different methods used in the calculation. As shown in Table 2, torsional parameters for Z-dol were calculated for four torsional angles that would describe twisting along the chain: F–C–C–F, F–C–C–H, H–O–C–C, and H–C–C–O. The potential was parameterized to predict the energy barrier as well as the behavior of the profile near equilibrium. The F–C–C–F and F–C–C–H angles had the highest energetic barrier to rotation at ~20 kJ/mol indicating less twisting relative to the H–O–C–C and H–C–C–O (~15 and ~17 kJ/mol, resp.). Similar energy torsional behavior can be expected for Z-tetraol.

Table 1: Selected parameters for the intramolecular Zdol atomistic model.

Bonds

(kJ/mol Å²)

(Å)
Angles

(kJ/mol rad²)

(rad)

O–H
4850.4
0.9677
H–O–C
619.6
1.857

C–O
2656.1
1.404
O–C–H
720.5
1.883

C–H
3081.5
1.101
O–C–C
886.8
1.939

C–F
2689.2
1.354
H–C–H
454.7
1.895

C–C
2289.8
1.527
C–O–C
1329.3
2.077

Table 2: The torsional parameters calculated from the ab initio energy profiles for the five-term truncation of the potential energy function of (2). The coefficients have units of kJ/mol.

Atoms

F–C–C–F
9.986
−29.852
0.996
39.724
−0.791

F–C–C–H
10.798
−32.389
0.000
43.187
0.000

H–O–C–C
6.971
−20.913
0.000
27.884
0.000

H–C–C–O
9.303
−24.761
–5.309
33.390
4.338

The type of intramolecular analysis described here is necessary for accurately simulating an atomistic system with the numerous degrees of freedom of PFPEs. It provides the basis for atomistic molecular dynamics useful in validating higher scale models as will be described in a later section on integration.

2.2. Intermolecular PFPE Interactions

Although the functional PFPE end-groups are designed to have strong COC surface adhesion during disk rotation, the dynamic property of lubricant replenishment, or self-healing, capability must be addressed to ensure consistent surface coverage during intermittent head-disk contact and lubricant pick-up by the head. Superior lubricant replenishment requires rapid molecular mobility which is influenced by the strength of intermolecular interactions. In essence, lubricants with the strongest molecular interaction may have more sluggish self-healing performance. In determining the best lubricant to meet HDI performance criteria, we introduce the concept of PFPE nanoblends to obtain optimal material performance. This allows flexibility in designing lubricants which can adhere to the COC while simultaneously having superior self-healing properties. Thus, we evaluated interactions between lubricants of the same type (pure) and between lubricants of different forms (blend). The blended PFPEs with different forms could be a novel concept in designing the optimal lubricant. Our work on atomistic intermolecular interactions can be found in [18]. Several previous studies used molecular and mesoscale simulation techniques to model lubricant films at the HDI [19–23]. However, in this work we evaluated the intermolecular interaction as a function of end-group structure using ab initio theory to make inferences regarding the lubricant replenishment performance. End groups were truncated into model structures, as in the previous section, retaining the characteristic functionality. The summary of model structures used for the analysis is shown in Figure 2. Here, Z-dol, Z-tetraol, and dipropylamine (DDPA) are used as benchmark PFPE structures for the calculations allowing us to compare hydroxylated and nonhydroxylated functional forms. Dimers composed of two PFPE structures are relaxed, and a frequency analysis is performed to classify the stationary point as a minimum on the potential energy surface. The uncorrected interaction energy represented by the potential energy of the configuration is calculated as

where is the energy of the dimer and and are the energies of isolated monomers and , respectively. In simulating each case, a fixed position of the lubricant is used as a starting point, and a geometry relaxation was performed to obtain the energy at the equilibrium structure along with a frequency analysis to ensure that the result is neither a saddle point nor a transition state but rather a minimum on the potential energy surface.

Figure 2: Model PFPE and COC materials for atomistic scale study.

The finite basis sets used in the calculations cause each monomer to share basis functions with other monomers in the dimer calculation, which is known as basis set superstition error (BSSE). As a consequence, the interaction energy is overestimated by (3). Thus, to reduce the BSSE associated with our computational techniques, a counterpoise correction (CP) [24] is utilized for (3):

where is the energy of monomer at the dimer geometry and is the energy of monomer at the dimer geometry with the added basis functions of monomer , with an analogous treatment for monomer .

Figures 3 and 4 give the relaxed configurations for the PFPE dimers. The shortest intermolecular bond among the hydroxylated end groups is between the hydrogen and oxygen of the OH groups. This indicates that interaction between the monomers is influenced by hydrogen bonding where the positively charged hydrogen is attracted to the electronegative oxygen. For DDPA, the shortest intermolecular bond occurs between the propyl group’s hydrogen and the fluorinated group pointing toward an attraction between hydrogen and the electronegative fluorine.

(a)

(b)

(c)

Figure 3: Pure relaxed lubricant dimers for (a) Z-dol, (b) Z-tetraol, and (c) DDPA. Red: oxygen, purple: nitrogen, blue: fluorine, grey: carbon, and white: hydrogen.

(a)

(b)

(c)

Figure 4: Blended relaxed lubricant dimers for (a) Z-dol-Z-tetraol, (b) Z-dol-DDPA, and (c) Z-tetraol-DDPA. Red: oxygen, purple: nitrogen, blue: fluorine, grey: carbon, and white: hydrogen.

Table 3 shows the counterpoise corrected energies of interaction for the dimers. For the pure complexes corresponding to Figure 3, the Z-dol and Z-tetraol have the strongest interactions of −34.4 kJ/mol and −59.0 kJ/mol, respectively, where the negative value denotes a stable interaction. A much weaker interaction is observed for the DDPA dimer with a value of 3.97 kJ/mol. The shortest intermolecular bonds are among the hydroxylated end-groups for Z-dol and Z-tetraol occurring between the hydrogen and oxygen of the OH groups. For DDPA, the shortest intermolecular bond occurs between the hydrogens of the propyl groups. Compared to DDPA, the hydrogen bonding among the hydroxylated end groups gives them a greater capacity to form strongly bound dimer complexes. This indicates that the hydroxylated PFPEs would have sluggish mobility compared to the DDPA forms due to the strong intermolecular interactions. However, condensed phase molecular simulation must be performed to establish a more concrete linkage between the strength of end-group interaction and molecular diffusion, particularly in the presence of finite temperature.

Table 3: The interaction energies for pure and blended PFPE end groups. All units are in kJ/mol.

Pure

Blended

Z-dol
−34.4
Z-dol-Ztetraol
−35.6

Z-tetraol
−59.0
Z-dol-DDPA
5.23

DDPA
3.97
Z-tetraol-DDPA
−6.53

For the blended structures in Figure 4, the closest intermolecular bond for the Z-dol-Z-tetraol dimer is between the hydrogen and oxygen of the hydroxyl group similar to the pure Z-dol and Z-tetraol cases. This similarity in orientation between the pure and blended cases indicates that the mechanism of interaction is still driven by hydrogen bonding when the hydroxylated lubricants are blended. The Z-dol-DDPA complex in Figure 4(b) shows that the closest interaction takes place between fluorines on each monomer. The orientation of the dimer represents the case where the portion of the end group closest to the backbone interacts. Alternatively, as shown in Figure 4(c), the Z-tetraol-DDPA dimer has a close contact between the hydrogen of the first OH group and the DDPA fluorines. This type of orientation is expected to produce a more stable interaction than the Z-dol-DDPA orientation.

As given in Table 3, consistent with the configuration of Figure 4(a), the Z-dol-Z-tetraol dimer has the strongest interaction of −35.6 kJ/mol. This is also similar to the strength of the interaction for the Z-dol and Z-tetraol pure dimers providing further support for the point that hydrogen bonding is the key binding mechanism in both the pure and blended cases. Thus, there are no significant changes in the interaction when going from the pure case to the blended case for the hydroxylated end groups. As the orientation in Figure 4(b) indicates, the interaction for the Z-dol and DDPA is unstable with a value of 5.23 kJ/mol. In contrast, the dimer of Figure 4(c) has a stable interaction of −6.53 kJ/mol most likely due to the interaction between the fluorines and hydrogens. Thus, we have obtained both weakly repulsive and more strongly attractive dimer configurations for hydroxylated end groups with DDPA. In order to translate the strength of the end group functionality to lubricant self-healing performance, we are developing intermolecular force field parameters for molecular dynamics simulation where temperature effects can be evaluated.

We have examined intermolecular interactions among various model PFPE lubricants which augment our earlier studies of intramolecular interactions. For the pure lubricant complexes, we found superior intermolecular interaction between lubricant end groups with hydrogen bonding capabilities. However, the pure DDPA complex displays longer intermolecular bond lengths as well as a slightly repulsive interaction. The binary lubricant blends demonstrated strong interaction between blended hydroxylated lubricants and diminished strength of interaction when adding DDPA to hydroxylated lubricants. Therefore, we have identified the ability to tune lubricant interaction by introducing blends of various functional end-group types. This result indicates that self-healing properties may be tailored by selecting optimal blend ratios which can be investigated using optimization formulations. Our capabilities in determining force constants that govern both intermolecular and intramolecular degrees of freedom allow us to make a linkage between ab initio atomistic descriptions and molecular/mesoscale many-molecule models via coarse-graining procedures.

2.3. PFPE-COC Interactions

Although the interaction strength between PFPE functional groups is one component which governs lubricant behavior, the PFPE-COC interactions are also key, since lubricants need strong COC adhesion to prevent spin-off and resist head pick-up. The details of the PFPE-COC interaction analysis can be found in [25]. To simulate PFPE adhesion on conventional COC material, we construct a model DLC structure which includes both sp³- (diamond-like) and sp²- (graphitic) type carbon as in Figure 5. Figure 5 gives the relaxed geometries for Z-dol-, Z-tetraol-, and DDPA-DLC dimers. The configurations indicate a stronger attraction for the sp² bonded carbon compared to the sp³ form. The interaction energies are summarized in Table 4. All three PFPEs have large attractive interactions with the DLC where the highest value of −126.02 kJ/mol is observed for the DLC-Z-tetraol dimer. The closest contacts for Z-dol and Z-tetraol are between the PFPE and graphitic carbons, while for DDPA, the fluorines and DLC hydrogens. Based upon the interaction energies, it is predicted that Z-tetraol would adhere more strongly to DLC with potentially less spin-off. In addition to the chemical composition of the DLC, the dimensional characteristics could also contribute to enhanced adhesion since the structure the PFPE encounters is not smooth. One may expect a less smooth structure in the DLC case than in the case of pure graphene or diamond due to the processes used in manufacturing the sp²/sp³ hybrid.

Table 4: Counterpoise corrected interaction energies for relaxed PFPE-COC dimers. All units are in kJ/mol.

PFPE
COC

DLC
Diamond
Graphene

Z-dol
−104.2
−38.0
−80.2

Z-tetraol
−126.4
−63.8
−86.9

DDPA
−76.0
−39.9
−73.5

(a)

(b)

(c)

Figure 5: The relaxed PFPE-DLC dimers for (a) Z-dol, (b) Z-tetraol, and (c) DDPA. White = hydrogen, red = oxygen, grey = carbon, blue = fluorine, and violet = nitrogen.

Although the PFPE-DLC relaxations indicate that the PFPE is more strongly attracted to the graphitic form of carbon, we systematically examine the effect of carbon structure on PFPE adhesion by independently calculating interactions with diamond and graphitic carbon. To simulate PFPE interaction with pure diamond without periodic boundaries, we choose a simplified structure of sp³ bonded carbons with hydrogens terminating all dangling bonds. Figure 6 gives the relaxed PFPE geometries with pure diamond. Table 4 shows that the interactions for all three PFPE-diamond dimers are attractive. As with the PFPE-DLC dimers, Z-tetraol has the highest interaction energy of −63.8 kJ/mol with similar values of −38.0 kJ/mol for Z-dol and −39.9 kJ/mol for DDPA. For the diamond case, the Z-dol and Z-tetraol hydroxyl groups have the closest contact with the diamond hydrogens. For DDPA, the fluorines have the closest contact with the diamond hydrogens. All three PFPEs exhibit weaker interactions with the pure diamond than they did with the DLC structure. This supports the observation in the DLC relaxed structures that the PFPEs were interacting more strongly with the graphitic structures as evidenced by their positions.

(a)

(b)

(c)

Figure 6: The relaxed PFP-diamond dimers for (a) Z-dol, (b) Z-tetraol, and (c) DDPA. White = hydrogen, red = oxygen, grey = carbon, blue = fluorine, and violet = nitrogen.

We utilize a single graphene flake in order to represent the purely graphitic carbon. Figure 7 shows the relaxed geometries for the PFPEs with a graphene flake, and the interaction energies for the PFPEs with graphene are given in Table 4 where the Z-tetraol has a slightly stronger interaction. For all three cases, the interaction energies are larger than those for the PFPE-diamond cases. This again indicates that the PFPE end-group structures are more strongly attracted to the graphitic carbon. The hydroxyl group interacting with the graphene carbons gives the closest contact in the case of Z-dol, and the Z-tetraol hydrogens have the closest contact with graphene. DDPA’s propyl hydrogens have the closest graphene contact for the complex. For the graphitic carbon, the high electron density in the sp² bonds allows for transfer of charge to the PFPE, which most likely leads stronger bond than what is possible in the PFPE-diamond complex. The closest interatomic contacts are between the PFPE hydrogens and graphene carbons in Figure 7 supporting the existence of charge transfer between the hydrogens and sp²-bonded carbons. These results provide evidence that a higher sp² carbon fraction in COCs could produce stronger PFPE adhesion. Thus, this analysis demonstrates the power of molecular architecture to provide guidelines for achieving the desired material properties such as adhesion. The exact sp²/sp³ ratio which maximizes adhesion could be determined more precisely via optimization procedures.

(a)

(b)

(c)

Figure 7: The relaxed PFP-graphene dimers for (a) Z-dol, (b) Z-tetraol, and (c) DDPA. White = hydrogen, red = oxygen, grey = carbon, blue = fluorine, and violet = nitrogen.

Since graphene has the potential to be a paradigm-shifting overcoat due to the savings in thickness, we investigated in more detail the PFPE-graphene interaction. In the case of a graphene overcoat (GOC), a few layers may be necessary to adequately protect the media. The effect of the number of graphene layers on its properties is a topic that has been previously investigated in the literature [26, 27]. The interesting phenomenon in these studies is the dimensional crossover from a two-dimensional structure to a three-dimensional one. In our context, the dimensional crossover can be evaluated in terms of the effect of GOC thickness on PFPE interaction. Thus, we evaluated the effect of single- versus double-layer graphene on interactions with our model PFPEs: Z-dol, Z-tetraol, and DDPA. The interaction energies of each PFPE-GOC dimer in these test cases are examined.

Since the atomically thin graphene structure has never been utilized as a disk overcoat adjacent to PFPE lubricant, we studied, for the first time, the interaction between PFPE and a model graphene. The interaction energy results of Table 5 indicate that all of the end groups have stable interactions with the model graphene. However, DDPA’s interaction energy with the GOC is still much lower compared to Z-dol and Z-tetraol which is consistent with the intermolecular separations between the structures in Figure 8. To evaluate the effect of having a two-layer graphene surface on lubricant adhesion, an additional layer was added to the lubricant-GOC dimers as shown in Figure 9. The additional layer was added at the experimental graphene-layer separation of 3.35Å, and the interaction energy of the lubricant at this geometry was calculated. Table 5 shows that the lubricant adhesion improves with the introduction of the second layer, particularly for the Z-tetraol and DDPA forms. For the single-layer GOC, the Z-dol and Z-tetraol hydroxyl groups have the closest contacts, while the DDPA hydrogens are the closest for the DDPA-single GOC layer dimer. Similar close contacts are observed for the double-layer GOC cases. Thus, in addition to the extra protection to the magnetic media offered by an additional layer of graphene, our calculations indicate that more desirable lubricant adhesion can be achieved as well.

Table 5: The counterpoise corrected interaction energies for lubricant-GOC clusters for single-layer (denoted as ()) and double-layer (denoted as ()) graphene. All values are in kJ/mol.

PFPE
()
()

Z-dol
−7.75
−10.17

Z-tetraol
−6.11
−25.04

DDPA
−1.70
−16.37

(a)

(b)

(c)

Figure 8: (a) Z-dol-graphene equilibrium geometry, (b) Z-tetraol-graphene equilibrium geometry, and (c) DDPA-graphene equilibrium geometry. Red = oxygen, violet = nitrogen, blue = fluorine, grey = carbon, and white = hydrogen.

(a)

(b)

(c)

Figure 9: Two-layer graphene with (a) Z-dol, (b) Z-tetraol, and (c) DDPA. Red = oxygen, violet = nitrogen, blue = fluorine, grey = carbon, and white = hydrogen.

By our atomistic scale investigation, we have demonstrated a qualitative trend in end-group structure’s influence on adsorption on graphene. Our feasibility study indicates that graphene could replace existing COC when using conventional lubricants like Z-dol and Z-tetraol due to significant adhesive properties, and they may be superior to alternative lubricant materials like DDPA. However, Z-tetraol and DDPA exhibit the much stronger adhesion with the addition of a second layer. We further introduced the novel concept of a solid-state buffer layer between PFPE and graphene to add extra lubrication and protection of the GOC. Our study will eventually lead to rigorous investigation of the PFPE-PFPE and PFPE-COC interactions that can be tuned to develop atomistic design criteria in lubricant technology. Large size graphene sheets with significant defects including grain boundaries and its interaction with PFPEs are currently being investigated.

3. Integration of Atomistic and Molecular Scales

Although our atomistic scale study provides valuable information regarding the fundamental role of molecular architecture on material performance, the oligomeric structure of PFPE presents a challenge for efficient molecular/mesoscale computation of the extremely large number of atomic degrees of freedom in molecular simulations. Therefore, a method to reduce the number of degrees of freedom in the atomistic models of oligomeric PFPE is essential from the computational view point. Extensions beyond the atomistic scale are necessary to evaluate the lubricant/COC behavior in a manner more closely linked with realistic HDI materials. The power of the multiscale approach is the ability to base the higher scale descriptions upon the information gleaned from the lowest-scale descriptions.

Historically, reduced order modeling (ROM) [28] has been used to describe continuum/device scale models using fewer degrees of freedom as shown in Figure 10. In this paper, we describe how the ROM methodology techniques can be applied in bridging atomistic/molecular scale models of PFPE. In ROM, the full model (FM) is solved for various inputs to generate snapshots of the model output spanning a predefined domain. The dimension of the snapshot set is then reduced by projecting onto a fewer basis functions and expanding the output in terms of the reduced basis resulting in a reduced model (RM).

Figure 10: The full Z-dol model (above) is mapped onto the reduced Z-dol model (below). Red atoms = oxygen, blue atoms = fluorine, white atoms = hydrogen, and grey atoms = carbon.

The atomistic FM introduced here is defined via the following intramolecular Hamiltonian, :

Here, , , and are the coefficients of the bond, bending, and torsional potentials, respectively. The parameters , , , and are the instantaneous and equilibrium bond lengths and angles, respectively. represents a Legendre polynomial. These parameters for the FM simulation are generated via the method of Hessian eigenvalue analysis as previously described in Section 2.2. The snapshots are generated via FM molecular dynamics (MD) simulation which yields atom trajectories, , corresponding to the coordinates of an atom, , over time, . To reduce the dimension of the model, based on the snapshot data, we incorporate a mapping operator which maps clusters of atoms onto single RM sites. Our RM sites have a physical linkage to the FM. Noid et al. [29] showed that a center of mass mapping operator

is sufficient for consistency between an atomistic and coarse-grained (molecular) model in phase space. Here, we map onto , with the new dimension . Thus, with this reduction, we have the trajectories of only nm RM sites over time. From this analysis, for simplicity, we generate a new RM intramolecular stretching Hamiltonian, :

In our earlier analysis, MD of coarse oligomer models has been accomplished using a bead-spring model with only bond stretching constraints which can adequately capture chain conformational behavior at this larger molecular scale. The coefficients of are calculated by fitting the Boltzmann distribution

which describes statistics in the canonical ensemble, into the probability distribution of the data set obtained from mapping onto reduced space. Here, is the Boltzmann constant, is the absolute temperature, and is a normalization constant. We apply our ROM methodology to the lubricant PFPE Z-dol, where we reduce 42 atoms to 7 RM sites. Table 6 shows selected bond and angle parameters for the FM intramolecular Hamiltonian. After the ROM procedure, as illustrated in Figure 10, we obtain the 6 RM force constants corresponding to (7) where selected values are given in Table 6. Compared to the stretching parameters of the FM, the RM parameters have significantly smaller values indicating more flexible bonds between the RM sites. Since the angle bending parameters in the FM are consistently less than the stretching parameters, it is justified to describe the RM using the bead-spring approach since RM bending would be a relatively insignificant constraint. In summary, we have applied a traditionally continuum/device level ROM methodology to link atomistic/molecular scale models which can span the entire scale integration, as shown in Figure 10, yielding a full scale HDI design.

Table 6: Parameters for reduced model intramolecular Hamiltonian.

RM Bond
(kJ/mol Å²)
(Å)

End-backbone
244.09
2.52

Backbone-backbone
179.20
2.45

In the next section, we move our focus to coarse-grained descriptions of PFPEs established on RM Hamiltonians such as those described in this section. Based upon this lower resolution molecular scale model (bead-spring chains), we investigate the static and dynamic properties of PFPEs via both equilibrium and nonequilibrium MD. In this higher scale approach, the functionality is dependent upon bead type instead of the detailed atomistic architecture. In addition, the nanorheology and nanomechanics are explored to determine the effect of lubricant blend composition and chain structure on the molecular/mesoscale behavior that is critical in the HDI application. We also extend our atomistic calculation of various PFPE/COC combination properties, as was done in Section 2.3, to the molecular scale by utilizing the atomistic scale interaction results.

4. Molecular Simulation

4.1. Coarse-Grained Model

PFPE molecules were characterized via a coarse-grained MD (CGMD), which simplifies the detailed atomistic structures by categorizing groups of atoms with uniform bead size for backbone and end group based on the dimension of the oligomer. This approach, nevertheless, retains the essence of molecular structure and chain flexibility. To implement the atomistic level theory in MD models, we applied a coarse-graining procedure linking atomistic potentials and molecular/mesoscale simulations combining the parameters obtained from quantum mechanical calculation, as illustrated in Figure 11.

Figure 11: Illustration of coarse-graining procedure by using ROM for PFPE molecular systems.

A flat and smooth surface was assumed to model the COC. Coarse-grained bead-spring PFPE is composed of a finite number of beads with different physical or chemical properties. For simplicity, we assume that all the beads, including the end beads, have the same radius. Force field in this model is constructed by Lennard-Jones (LJ1) and van der Waals (LJ2) potentials for dispersive bead-bead and bead-wall interactions, respectively (Figure 12). For polar interactions, exponential potential functions (EXP1 and 2) were added to both bead-bead and bead-wall cases. For the bonding potential between adjacent beads in the chain, a finitely extensible nonlinear elastic (FENE) model was used. For example, PFPE Z-dol can be characterized differently from PFPE Z by assigning the end bead a polarity originating from the hydroxyl group in the chain end.

Figure 12: Potential energies for the bead-spring model. LJ1: Lennard-Jones potential, LJ2: van der Waals potential, EXP1 and EXP2: short-range polar potential, FENE: finitely extensible nonlinear elastic potential.

Neighboring beads are connected via a FENE spring as shown in the following equation:

where is the bond length, is the equilibrium bond length, is the maximum extensible range of the spring, and is the spring constant which quantifies the rigidity of the bond.

All beads, including the end beads in PFPE Z-dol, interact with each other by a pairwise, dispersive, truncated Lennard-Jones (LJ) potential as shown in the following equation:

where is the distance between beads, represents the well depth of interaction, and is the bead diameter.

For the dispersive bead- (db-)wall interaction, we used the attractive potential, :

Here, is the distance from the wall and is the potential depth.

Previous studies have indicated that the coupling between end groups and their binding to the solid surface is more like a hydrogen bonding interaction [30] from an ab initio calculation [31]. While other studies considered the possibilities of hydrogen transfer [32] and esterification as a result of annealing [33]. In the functional PFPE model, the polarity interaction is assumed to happen within a short range and is modeled as an exponential decay function. Equation (12) describes the potential function among end beads:

where represents the interaction strength between end beads (Table 7), is the critical distance between a pair of end beads, which is defined by the diameter of an endbead, and is the decay length in the exponential function, which must be small to meet the criteria of the short-range interaction.

Table 7: Endbead-endbead interaction () between PFPEs, that is, Z-dol, Z-tetraol, and DDPA (kJ/mol).

PFPE

Z-dol
Z-tetraol
DDPA

Z-dol
−34.46
−35.63
5.23

Z-tetraol
—
−59.03
–6.53

DDPA
—
—
–3.98

Further, end beads interact with the solid wall as

where is the critical distance of the interaction. Note that the larger the (Table 4), the stronger the interaction.

By varying molecular structures and inputting ab initio-based , atomistically architectured PFPEs are modeled in this simulation. Monolayer and multilayer films were investigated within the simulation box dimensions in the and directions of with NVT ensemble, where is the diameter of the bead. After the system equilibration, the simulations were performed for 40 000τ with the time step of 0.005τ, where is , and is the mass of the bead. Periodic boundary conditions were applied in the and directions.

4.2. Molecular Dynamics

The equation of motion can be expressed in the following Langevin equation form [34–36]:

Here, spans and spans , where and imply the number of PFPE molecules and monomers, respectively, is the mass of the beads, and is the potential energy of the system expressed as some of all potential energy given in (9)–(13).

is the frictional tensor. We assume it to be isotropic; that is, ( is unit tensor). is assumed to be Gaussian white noise, which is generated according to the fluctuation-dissipation theorem [35, 37]

Here, the angular brackets denote an ensemble average, which is the same as time average from the ergodic hypothesis. and are Kronecker delta and is Dirac’s delta function.

4.3. Static and Dynamic Properties

We examined the molecular conformations of the confined film via the anisotropic radius of gyration () in parallel () and perpendicular () components:

where , , and are the locations of the chain’s center of mass. The anisotropic radius of gyration as a function of distance from the wall (Figure 13) demonstrates that the chain conformation becomes more oblate as the wall distance decreases.

(a)

(b)

Figure 13: Radius of gyration as a function of the distance from the wall, , for : (a) parallel component and (b) perpendicular component for (■) , (●) , and (▼) . Error bars indicate one-half of the standard deviation and are shown only for .

The chains in the film extend more in a perpendicular direction than in a parallel direction, resulting in a layered structure. However, the endbead interaction does not significantly affect the conformation of both and , where oscillates around that of chains with .

The difference plotted in Figure 14 shows that maximizes near and and minimizes near and ( represents the mean perpendicular radius of gyration.). This conformational fluctuation is due to the coupling and orientation of the functional end beads.

Figure 14: for the chains () with functional end beads for (●) and (▼) . ( is the distance from the wall.) Error bars indicate one-half of the standard deviation and are shown only for .

Figure 15 shows the density profile as a function of the distance from the wall. The polymer chains () with functional end beads ( and ) exhibit pronounced peaks, while the nonfunctional end beads () are monotonic. The end beads are believed to couple near the peak. Because the end beads in the first layer attach to the wall due to the attractive functional interaction, the orientation of end beads is induced in the second layer, and likewise, the orientation propagates to subsequent layers. In the case of the chains located within the end beads’ density peak, the functional end beads also stay within the peak to couple with other end beads resulting in a flatter conformation. The chains located out of the peak unfold themselves so that their end beads can reach other end beads. The characteristic oscillation in the molecular conformation (Figure 14) and the density profile (Figure 15) is induced by functional end beads.

(a)

(b)

(c)

Figure 15: Density profiles of the end beads as a function of the distance from the wall, : (a) and (b) , and (c) .

To establish a qualitative relationship between orientation and the layer structure, we examined the number of layers at , and the endbead density is plotted in Figure 16 for . The bead density profile has three peaks: at the wall, between the second and third layers, and between the fourth and fifth layers. The adsorption of functional end beads results in an alternate ordering in the subsequent layers, that is, up in the second layer, down in the third, and so forth. Our result provides a direct interpretation of experimental surface energy data for PFPE films with functional end groups [30]. The nondispersive component of surface energy exhibited an oscillatory pattern with increasing film thickness and was shown to be approximately proportional to end-group density, as demonstrated in our previous study [30, 38]. Our simulation results also suggest that the density variation of the end groups is related to the characteristic behavior of the surface energy of PFPE films.

Figure 16: Density profiles of end beads and number of layers (, and ) as a function of the distance from the wall ().

Furthermore, functional end-group density variation could qualitatively explain the spreading profiles of PFPE films. The surface diffusion coefficient obtained from the experimental spreading profile depends strongly on the layer positions via end-group orientation [39]. That is, the diffusion coefficient between the second and third layers is smaller than that between the first and second layers because of the endbead coupling, which is consistent with the observation made by Ma et al. [40], who also provided a qualitative explanation for the anomalous peaks in thickness-dependent diffusion coefficients.

As temperature increases, the layer expands and the orientation of end beads smears. This is shown in Figure 17. The observed expansion of layer thickness is attributed mainly to the temperature dependence of the intermolecular interaction and excluded volume effect not due to bond stretching. Temperature dependence, which is pertinent to the annealing process of thin polymeric lubricant films, has been carefully examined by Hsia et al. [41].

Figure 17: Temperature effect on the density profiles of end beads () as a function of the distance from the wall (). (dot-dashed line), (dashed line), (dotted line), and (solid line).

The steady-state nanoscopic properties were examined, and the result agrees qualitatively with the simulation results obtained via the simple reactive sphere (SRS) model discussed. Our results provide a direct interpretation of experimental surface energy data for PFPE films with functional end groups [30, 38, 42, 43]. Stable films did not experience dewetting or film rupture. However, a rougher surface morphology was observed for smaller molecular weights (Figure 18(b)) and strong endbead functionality (Figure 18(c)). Visualizations of surface morphology could provide powerful tool for describing air bearing stability.

(a)

(b)

(c)

Figure 18: Morphology of PFPE films from MD simulations: (a) schematic of simulated surface roughness for (b) molecular weight dependence ( (upper), 15 (middle), and 20 (bottom)) and (c) endbead functionality dependence ( (upper) and (bottom)).

Another noteworthy observation is the segregation, or localization, of chain-ends at the film surface, which was found both in the nonfunctional and functional chain-ends. The endbead density obtained from the MD simulations is normalized with respect to the total number of beads and plotted in Figure 19 as a function of endbead functionality. The segregation of the nonfunctional chain-ends at the surface () assures that the driving force is entirely entropic. The chains can avoid the decrease in entropy by having the nonfunctional ends localized at the film surface where the number of configurations is reduced due to chain folding. The mildly functional end beads () still localize at the film surface, which suggests that the entropic contribution is more dominant, and, as a result, the free energy of the system is minimized by having the end beads localized. The strongest functional end beads, on the other hand, are totally depleted from the surface (). The contribution of surface enthalpy to the free energy is too large to be compensated by the entropic component. In this case, therefore, all of the functional end beads minimize the energy of their functional interaction by coupling with other end beads inside the film. In the same manner, the segregation of small molecular weight fraction at interfaces was observed from our binary MD simulations.

Figure 19: Normalized endbead densities as a function of distance from the surface () for (a) , (b) , (c) , and (d) .

For the dynamic property, we will examine the self-diffusion process of a tagged PFPE molecular center of mass from the simulation to achieve insight into the excitation of translational motion, that is, spreading and replenishment. The squared displacement of the center mass of a molecule or a bead is used as a measure of translational movement. The self-diffusion coefficient, , can be represented as a velocity autocorrelation function:

where denotes the velocity vector of the center of mass for a PFPE molecule. Figure 20 illustrates the results for the and (the components of the diffusion coefficients in the directions parallel and perpendicular to the solid surface) obtained from our preliminary MD simulations.

Figure 20: Self-diffusion coefficients () as a function of distance from the surface for nonfunctional () and functional chain-ends ().

We may also study the molecular relaxation process, which has been investigated by the autocorrelation function of normal modes for a linear polymer chain [44]. The relaxation spectrum can be analyzed by the Kohlrausch-Williams-Watts function [45]:

Here, represents the elementary relaxation time and the stretching parameter is descriptive for the nonexponential decay process (or distribution of relaxation times), or it physically quantifies the confinement effect; that is, the larger the is, the stronger the confinement becomes. Figure 21 demonstrates the contribution of the solid surface to the relaxation process for end-functional PFPEs () by comparing the bulk () with nanofilm (). A longer is observed for the nanofilm than for the PFPE bulk, but the MW dependence, that is, , is weakened. For instance, for the bulk as the Rouse model predicts [46], but it decreases to for the nanofilm. In addition, for the nanofilm increased slightly (about 0.15) with MW as seen from the solid line in Figure 22. Note that the end-group ratio, , decreases for higher MW and that the functional end groups interact strongly with the solid surface. Consequently, the higher the end-group ratio in a PFPE chain is, the stronger the PFPE-solid surface interaction results. Compared with lower-MW PFPEs, the effective interaction with the surface reduces for higher MW PFPEs.

Figure 21: as a function of (: the number of beads per molecule) for the PFPE () in the bulk state () and in the nanofilm () with .

Figure 22: The stretching parameter and end-group ratio as a function of for the PFPE () in the bulk state () and in the nanofilm () with .

4.4. Various PFPEs and COCs

For the future HDI design, which increases the stability while reducing the configuration scale to obtain the higher sensitivity of the head and magnetic layer interactions, we introduced various COC materials including graphene. Due to the advantages of graphene including the single atomic thickness of graphene and superior mechanical and thermal characteristics, HDI with graphene will provide significant reduction in HMS along with high reliability if we optimize the lubricant structures (Figure 23). Various PFPEs properties on the selected COCs were also investigated and compared. As shown in Figure 24, all PFPEs (Z-dol, Z-tetraol, and DDPA) show a flat conformation () on three different carbon surfaces. By comparing the end-group functionalities listed in Tables 4 and 7, Figure 25 illustrates that the strength of interactions among end groups and the surface barely affects the perpendicular molecular conformation as long as the total number of beads (i.e., molecular weight) remains constant. Particularly, similar values for all cases, which mainly determine the film thickness in a monolayer lubricant system, indicate that intramolecular factors such as chain flexibility, molecular weight, or molecular structural complexity may affect film thickness and conformation more than intermolecular factors such as end-group functionality or COC material.

Figure 23: HDI configuration and graphene for the COC layer.

Figure 24: The ratio of to for each component of the material matrix.

Figure 25: The perpendicular component of radii of gyration () for each component of the material matrix.

Although the molecular conformations were identical, we observed that the end-group clusters on the surface vary for different PFPEs. Figure 26 shows the bottom view of Z-dol, Z-tetraol, and DDPA, where the end groups are located on the surface. Both Z-dol and Z-tetraol (Figures 26(a) and 26(b)) exhibit clusters and network structures, and the end groups of Z-dol form smaller clusters, while DDPA lubricant film, which has the lowest end-group interaction energy (–3.98 kJ/mol) between end groups (Table 7), does not display any clustering. In addition to the PFPEs on the diamond surface in Figure 26, we also found that PFPEs on the other surfaces including novel graphene surface (i.e., DLC and graphene) form similar end-group distributions as those on the diamond surface.

(a)

(b)

(c)

Figure 26: Bottom view snapshots of (a) Z-dol (pink: end-group), (b) Z-tetraol (red: end-group), and (c) DDPA (blue: end-group) on diamond overcoat surface (light blue: backbone).

From this result, we observed that the clusters formed in PFPE films depend on the interactions between end groups (Z-tetraol > Z-dol ≫ DDPA) while the interaction between carbon surfaces and the end groups merely affects the size of clusters.

Figures 27 and 28 compare self-diffusion coefficients of each component for Z-dol, Z-tetraol, and DDPA on the diamond, DLC, and graphene. The comparison between and values indicates that is the dominant diffusion process since perpendicular diffusion is confined by the monolayer thickness of PFPE films. The of PFPEs is ranked in reverse order compared to the interaction strength between end groups (Z-tetraol < Z-dol < DDPA). Also, the self-diffusion is affected by the carbon surface structure. Z-dol exhibits similar mobility on all carbon surfaces, while Z-tetraol shows minimal mobility which is sensitive to the carbon surface. DDPA, consisting of a different type of end group from Z-dol and Z-tetraol (hydroxyl end groups), shows significantly higher self-diffusion coefficients on the diamond and graphene structure COC. When the molecular conformation of PFPEs is considered, this result indicates that Z-tetraol may not be the best choice for the COC beyond the conventional materials, since a significant disadvantage in the self-healing ability exists. However, we see that the mobility of DDPA, which is similar to Z-dol on the conventional DLC overcoat, significantly increases on the graphene as well as diamond.

Figure 27: Parallel component of self-diffusion coefficients () for each component of the material matrix.

Figure 28: Perpendicular component of self-diffusion coefficients () for each component of the material matrix.

Here, we found that the changes in the monolayer film conformation affected by the perpendicular profile of PFPE molecular structures on the surface are insensitive to the kind of end groups or carbon surfaces. However, we observed that the end-group distributions on the surface depend on the end-group-end-group interactions. Also, different chemical types of end groups on PFPEs vary the mobility of PFPEs with respect to the carbon surface characteristics. In addition to the various linear functional PFPEs examined thus far, we investigated different types of PFPEs, which are atomistically/molecularly architectured functional lubricant molecules on the existing and novel COC layers (Figure 29). Since the multiscale methodology provides the advantage of scale integration, complicated combinatorial problems of PFPFs and COCs and numerous parameters can be estimated, and the optimal lubricant design criteria for the future HDI will be selected by considering the performance parameters from the bottom scale.

Figure 29: Different types of PFPEs with atomistically/molecularly architectured functional molecular structures on the various COC layers and output parameters.

Star-like molecules (TA-30 and QA-40), which prevent lubricant decomposition at high temperature by eliminating the acid labile unit (–OCF₂O–), were also investigated via coarse-grained model (Figure 30). It has been reported that an additional functional group in ZTMD affects the film conformation by anchoring the backbone [47, 48], which is adsorbed on the surface so that the structure results in a smaller film thickness compared to the linear PFPEs. Tani et al. compared branched PFPEs (e.g., TA-30 and D-4OH) with Z-tetraol by using the atomic force microscopy (AFM) and the molecular dynamics, and they found that TA-30 shows the highest surface coverage and lowest molecular height among the examined PFPEs [49]. Branches with functional end groups in star-like PFPEs also form backbones, anchored on the surface; however, TA-30 and QA-40 may exhibit different film conformations from the ZTMD due to the flexible branch backbone. To understand this, we specifically examined structural effects of anchoring arms on TA-30 and QA-40 on the film thickness. As shown in Figure 31, star-like PFPEs form similar end-group anchored structures while QA-40 exhibits more severe network structure due to the additional functional arm, which may provide the film stability as well as diminishing mobility (bottom view). Since QA-40 forms larger clusters inside and at the bottom of the film, less functional end groups on the surface of the film represent the reduction of film surface energy, improving the stability in the near-contact HMS system. The flexibility of arms and end-group interactions results in less flat molecular conformation as well as thicker film as shown in Figure 32.

Figure 30: Coarse-grained bead-spring models of linear PFPEs (Z-dol, Z-tetraol, and ZTMD) and star-like PFPEs with three and four arms (TA-30 and QA-40), where magenta and red beads represent hydrogenated functional groups.

Figure 31: Snapshots of multilayer and monolayer TA-30 and QA-30: side and bottom views.

Figure 32: Schematic description of molecular conformations of ZTMD, TA-30, and QA-40 on COC.

4.5. Nanorheology and Mechanics

In addition to the inherent PFPE/COC properties, we further explored the nonequilibrium dynamics of HDI operation and PFPE/COC responses. By adopting the coarse-grained model described in previous sections to the nonequilibrium MD, we have examined the nanorheological and nanomechanical responses for the nonfunctional PFPEs as a function of wall separation/film thickness and the molecule-solid surface interaction. We also examined the effect of end-group functionality and further compared the relaxation time of nanometer thick confined PFPE films with PFPE bulk.

4.5.1. Nanorheology

There has been great interest in resolving the rheological and tribological issues of PFPE lubricants during the intermittent contact between the head and the disk. We first examined the confinement effect. Figure 33 plots the shear-rate-dependent viscosity of nonfunctional PFPE () nanofilms with three different wall separation lengths (), where and are endbead functionalities. As , the shear viscosity is slightly larger than the bulk value. However, as the wall separation gets smaller, for example, or even , which contains only 1-2 layers of PFPEs, the film shear viscosity increases significantly, indicating a much stronger external confinement to the examined films. Thus, the confinement depends nonlinearly on the wall separation.

Figure 33: The shear viscosity of nonfunctional PFPE () nanofilms with three different wall separations (, , and ) as a function of shear rate.

We found that the surface confinement plays a vital role in determining the shear viscosity of PFPE nanofilms. Nonfunctional PFPE nanofilms confined between two solid surfaces with a separation of under different surface interaction strengths were simulated. The calculated shear viscosity at the shear rate of as a function of end-group functionalities is shown in Figure 34, where a second degree of polynomial dependence is found. Note that this dependence may vary with the wall separation and shear rate.

Figure 34: The shear viscosity of nonfunctional PFPE () films with the thickness of as a function of surface interactions () at the shear rate of . Note that 230 functional PFPEs with are confined in the unit cell of .

Using the conventional rheometer, Kono et al. [50] experimentally examined the rheological properties of PFPE bulk systems for the shear rate range of 10⁻⁶~10³ sec⁻¹. The bulk rheology examines the dynamic response of PFPE system under deformation (e.g., shear or elongation) while the nanorheology for PFPE, for the first time, provides additional information about PFPE relaxation (or shear viscosity) due to the nanoscale confinement as well as interaction between surface and confined molecules. Due to the ease in assigning external conditions, that is, shear/normal stress into the simulations, molecular dynamics (MD) turns out to be an alternative approach to understand the above-mentioned issues. Here, by modifying our previous simulation model [22] and incorporating two confined walls, we studied the rheological responses of PFPE nanofilm as a complementary tool to estimate the nanoscale lubricant tribological performance.

We have examined the viscosity dependence on the endgroup functionality. Figure 35 plots the shear-rate-dependent viscosity of both nonfunctional () and functional () PFPE films confined between two solid surfaces with a separation of . Here, the “shear thinning” behavior is observed as the shear rate increases. In addition, when compared to the shear viscosity for the bulk PFPE system (O(~10¹)) [51], the solid surfaces have significantly confined dynamic behavior of nonfunctional PFPE films (square symbols), which leads to a huge increase (O(~10⁴)) in the shear viscosity at the low shear rate region. On the other hand, the confined functional PFPE shows less increase in the steady shear viscosity than the bulk state (O(~10³)) [51]. A similar increase (O(~10²)) in the steady shear viscosity for the functional PFPE films is also observed (circled symbols). Therefore, we found that although functional end-group coupling can slow down the PFPE dynamics, the surface confinement plays a more dominant role in the rheological properties of PFPE nanofilms.

Figure 35: The shear viscosity of both nonfunctional () and functional () PFPE films confined between two solid surfaces with a separation of and with under .

4.5.2. Nano-Mechanics

During the HDD operation, the intermittent contact between the slider and lubricant surface inevitably results in lubricant pick-up, where the lubricant molecules could be partially transferred to the air bearing surface due to the horizontal air shear stress as well as the vertical head-disk interaction and eventually affect the tribological performance of the HDI integration, that is, lubricant depletion and head vertical displacement. Therefore, nanomechanics of the PFPE films, such as “compression” and “tension” along the normal coordinate, are of technical significance in the HDI integration. In addition, nanomechanics of PPFE films allows us to explore several interesting nano-lubrication fundamentals, that is, film surface tension, disjoining pressure, and work of adhesion and separation, which have been attempted via both experimental [52] and theoretical [53] tools. Here, using molecular dynamics (MD) simulation [54], we examined the nanomechanics of PFPE films, including “compression” and “tension,” where the lubricant end-group functionality is found to be a key determining factor in the event of HDI contact. Based upon the calculated normal stress profile, we further investigated the rheological property of PFPE nanofilms via the -modes Maxwell model.

In our “thought nanomechanics experiment,” 150 functional PFPE molecules () with the chain length of 20 beads per molecule were randomly “coated” onto two solid surfaces facing each other with the horizontal dimension of and a vertical separation of 25σ. The simulated nanomechanics snapshots are captured in Figure 36. The MD simulation was first executed to relax PFPE films to the equilibrium state (Figure 36(a)), where functional end beads in red are either localized near the solid surface or clustered. The “compression” mode was then carried out by moving down the top surface at a constant speed of 0.002σ/τ. As the separation became less than the LJ interaction range, interfacial molecules began to interact with each other (Figure 36(b)) and eventually merged into single nanoconfined film (Figure 36(c)). During the “tension” process, while the top surface was retracted with the same constant speed, a fluid bridge was found in between (Figure 36(d)), became thinner and thinner as the wall separation increased (Figure 36(e)), and eventually separated into two individual films. A similar nanomechanics simulation was also performed for nonfunctional PFPE nanofilm as shown in Figure 37, which showed significant difference compared to the functional PFPE nanofilm during the “tension” process. With the same time interval for data output, no apparent fluid bridge was observed for nonfunctional film. Also, the wall separation at the fluid bridge “break” point for nonfunctional PFPE film (~24σ) is much smaller than that of functional PFPE film (~33σ).

Figure 36: The simulated nanomechanics for functional PFPE () films: (a) two well-equilibrated films approach each other, so-called “compression”; (b) two films are in contact; (c) two films completely contact and merge into one confined film; (d) the top surface is retracted, so-called “tension”; and (e) the confined film is being separated into two parts in the nearly end of the “tension.” Note that light blue indicates the backbone beads and red represents the functional end beads.

Figure 37: The simulated nanomechanics for nonfunctional PFPE () films: (a) two well-equilibrated films approach each other, so-called “compression”; (b) two films are in contact; (c) two films completely contact and merge into one confined film; (d) the top surface is retracted, so-called “tension”; and (e) the confined film is being separated into two parts in the nearly end of the “tension.”

The corresponding normal stresses, ( is the coordinate perpendicular to disk surface), during “compression” and “tension” processes were further calculated for both nonfunctional and functional PFPE films. It is shown in Figure 38 that the film normal stress remains constant as the top surface moves downwards, but sharply drops when interfacial molecular attraction occurs and eventually reaches the minimum. The normal stress may be significantly increased afterwards with even a minor compressive strain due to the strong intermolecular repulsion force, that is, the potential. In fact, our “compression” mechanics study suggests ~1σ clearance gain for functional film than nonfunctional film before the HDI interaction takes place, indicating less lubricant pick-up event for functional PFPEs. During the “tension” mechanics, an irreversible normal stress profile (hysteresis) was observed, especially for the functional PFPE film, which could be relevant to the viscoelasticity and stiction of PFPE nanofilm and also depends on the retraction rate of the top surface.

Figure 38: The normal stress of both nonfunctional (dashed lines) and functional (solid) PFPE films with during the simulated nanomechanics, that is, “compression” and “tension.”

As one of the most classical viscoelastic models for a polymer system, the Maxwell model with multiple relaxation modes is employed here to understand the viscoelastic properties of PFPE nanofilms. The Maxwell model can be represented by a purely viscous dash pot with a viscosity of and a purely elastic spring with the spring constant of connected in series, as shown in Figure 39(a).

(a)

(b)

Figure 39: The schematics of the (a) single mode Maxwell and (b) -modes () Maxwell model. Note that is the spring constant and is the viscosity of the dash pot.

The relaxation process of polymeric/oligomeric systems usually cannot be described by single relaxation time; therefore, -modes Maxwell model is typically preferred (Figure 39(b): ), where the relaxation moduli are given by

Here, and () are the discrete spring constant and the associated relaxation times of the material, respectively. This idea can be translated from the time domain to the space domain in that the retraction rate () of the top wall remains constant. Therefore, (19) is employed to fit the simulated film normal stress with . This fit is only limited to the tension process as indicated by the blue symbols in Figures 40 and 41:

Here, “0” stands for the reference point, where two films are completely merged into one confined film before the repulsion occurs inside the film.

Figure 40: The nonlinear square fit of normal stress data using (19) for nonfunctional PFPE () films.

Figure 41: The nonlinear square fit of normal stress data using (19) for functional PFPE () films.

From the “tension” normal stress profiles, we observed that the PFPE molecules might experience two steps of relaxation, especially for the functional PFPEs, showing a sharp slope in the beginning of the “tension” process while the slope gets compromised once the wall separation is larger than 20σ in Figure 41. Therefore, is set to be 2 in (19). The fitted elastic moduli and via the -modes Maxwell model for both nonfunctional and functional PFPE films are listed in Table 8, where only one physical mode is observed for the nonfunctional PFPE film, related to the relaxation of the backbone beads, while a second mode with weaker elastic modulus but longer relaxation time is captured for the functional PFPE film, related to the relaxation of the functional end beads. Therefore, we believe that in the initial state of the “tension” process, the relaxation of backbone beads is dominant. However, the second stage of “tension,” which shows the apparent long-tail behavior, is mainly ascribed to functional end beads.

Table 8: The elastic moduli and relaxation times for the -modes Maxwell model for both nonfunctional and functional PFPE films.

PFPE
Nonfunctional
Functional

mode
1
2
1
2

(ε/σ ⁴)
0.29
0.09
0.44
0.12

(σ)
1.13
1.1 × 10¹⁹
1.68
18.36

(τ)
565
Unphysical
840
9180

5. Summary and Discussion

In this paper, we have described a multiscale framework for modeling HDI materials. Beginning with the atomistic scale, we investigated the effect of detailed atomistic architecture on intramolecular and intermolecular PFPE degrees of freedom. The bottom-scale models allow us to investigate the strength of atomistic interactions between various PFPE/COC material pairs motivated by the need to engineer mobile lubricants with strong COC adhesion. Via a bottom to top approach, the atomistic information is passed to the molecular scale via coarse-graining procedures. This coarse-graining involves developing Hamiltonians with reduced degrees of freedom based on atomistic conformations and interaction strengths. At the molecular scale, molecular dynamics allows for the investigation of bulk and surface PFPE properties for a range of functional group compositions. These properties include diffusion, radii of gyration, film morphology, and nanorheological performance. Thus, we have demonstrated a systematic methodology for probing the HDI; however, this approach can be applied to other physical systems with phenomena occurring over a range of time and length scales. Following this same methodology, extensions will be made to other critical issues in the performance of the HDI system. As discussed earlier, advanced COC materials such as graphene have the potential to provide significant savings in the COC layer thickness. However, factors such as the thermo-mechanical properties in the harsh HDI environment must be well understood. Thus, we apply our multiscale approach to investigating the mechanical performance of graphene under stress with various defects and grain boundaries via ab initio and molecular dynamics calculations [55]. Other COC configurations are also being considered with graphene reinforced by fullerenes in a novel solid-state buffer layer as shown in Figure 42.

Figure 42: A schematic of the HDI with fullerene buffer layer between the PFPE lubricant and graphene. The motion of the fullerene on graphene can provide additional lubrication.

Due to significant interest in heat-assisted magnetic recording (HAMR), we also study heat transfer in graphene and other nanoscale films via MD and the mesoscale lattice Boltzmann method (LBM) [56–58]. Consistent with our multiscale approach, the issue of HAMR’s effect on HDI materials is also investigated at the atomistic scale using density functional theory-based ab initio molecular dynamics (AIMD). The AIMD simulation is in the context of HAMR conditions and is a natural complement to the higher scale classical MD/LBM studies where the best performing materials at the bottom, high resolution scale are further evaluated at higher scales.

In our integration model, we obtained the mesoscale properties from the coarse-grained model, which enables the bottom-up approach. This implies that the mesoscale model can be systematically derived from ab initio simulations via sequential coarse graining. LBM is the mesoscale/continuum theory, which can be incorporated to the highest coarse-grained model (e.g., SRS model [59]) since the LBM is on particle-based assumption [60–65]. Therefore, LBM can be used as a fundamental building block for constructing mesoscale/continuum models descriptive for lubricant as well as air bearing system (ABS). Since the LBM has the nature of transient physics and simple algorithm for the parallel calculation, thermal/mechanical models will be incorporated with the integrated model also covering adsorption/desorption kinetics, evaporation, and thermal oxidative degradation of PFPE lubricant at elevated temperature.

Acknowledgments

This work is supported by Seagate Technology, and the publication fee was paid by HGST, a Western Digital Company.

References

D. Kim, P. S. Chung, P. Jain, S. H. Vemuri, and M. S. Jhon, “Multiscale modeling of head disk interface,” IEEE Transactions on Magnetics, vol. 46, no. 6, pp. 2401–2404, 2010.

L. J. Chen, H. J. Qian, Z. Y. Lu, Z. S. Li, and C. C. Sun, “An automatic coarse-graining and fine-graining simulation method: application on polyethylene,” The Journal of Physical Chemistry B, vol. 110, no. 47, pp. 24093–24100, 2006.

J. Fish, “Bridging the scales in nano engineering and science,” Journal of Nanoparticle Research, vol. 8, no. 5, pp. 577–594, 2006.

O. Al-Khayat and H. P. Langtangen, “Computational aspects of multiscale simulation with the lumped particle framework,” Communications in Computational Physics, vol. 12, pp. 1257–1274, 2012.

A. Lyubartsev, Y. Tu, and A. Laaksonen, “Hierarchical multiscale modelling scheme from first principles to mesoscale,” Journal of Computational and Theoretical Nanoscience, vol. 6, no. 5, pp. 951–959, 2009.

F. Muller-Plathe, “Coarse-graining in polymer simulation: from the atomistic to the mesoscopic scale and back,” ChemPhysChem, vol. 3, no. 9, pp. 754–769, 2002.

A. B. Mhadeshwar and D. G. Vlachos, “Hierarchical multiscale mechanism development for methane partial oxidation and reforming and for thermal decomposition of oxygenates on Rh,” Journal of Physical Chemistry B, vol. 109, no. 35, pp. 16819–16835, 2005.

T. Murtola, A. Bunker, I. Vattulainen, M. Deserno, and M. Karthunen, “Multiscale modeling of emergent materials: biological and soft matter,” Physical Chemistry Chemical Physics, vol. 11, pp. 1689–1892, 2009.

N. Sheng, M. C. Boyce, D. M. Parks, G. C. Rutledge, J. I. Abes, and R. E. Cohen, “Multiscale micromechanical modeling of polymer/clay nanocomposites and the effective clay particle,” Polymer, vol. 45, no. 2, pp. 487–506, 2004.

S. P. Xiao and T. Belytschko, “A bridging domain method for coupling continua with molecular dynamics,” Computer Methods in Applied Mechanics and Engineering, vol. 193, no. 17–20, pp. 1645–1669, 2004.

Q. Shi, S. Izvekov, and G. A. Voth, “Mixed atomistic and coarse-grained molecular dynamics: simulation of a membrane-bound ion channel,” Journal of Physical Chemistry B, vol. 110, no. 31, pp. 15045–15048, 2006.

K. Kremer and F. Müller-Plathe, “Multiscale simulation in polymer science,” Molecular Simulation, vol. 28, no. 8-9, pp. 729–750, 2002.

H. Wang, C. Junghans, and K. Kremer, “Comparative atomistic and coarse-grained study of water: what do we lose by coarse-graining?” The European Physical Journal E, vol. 28, no. 2, pp. 221–229, 2009.

R. Smith, P. Seung Chung, J. A. Steckel, M. S. Jhon, and L. T. Biegler, “Force field parameter estimation of functional perfluoropolyether lubricants,” Journal of Applied Physics, vol. 109, no. 7, Article ID 07B728, 2011.

J. M. Seminario, “Calculation of intramolecular force fields from second-derivative tensors,” International Journal of Quantum Chemistry, vol. 60, no. 7, pp. 1271–1277, 1996.

E. Anderson, Z. Bai, C. Bischof et al., LAPACK Users’ Guide, SIAM, Philadelphia, Pa, USA, 1999.

J. Wang, R. M. Wolf, J. W. Caldwell, P. A. Kollman, and D. A. Case, “Development and testing of a general Amber force field,” Journal of Computational Chemistry, vol. 25, no. 9, pp. 1157–1174, 2004.

R. L. Smith, P. S. Chung, S. H. Vemuri, G. Y. Yeom, L. T. Biegler, and M. S. Jhon, “Atomistic simulation method in head-disk interface of magnetic data storage systems,” Journal of Applied Physics, vol. 111, no. 7, Article ID 07B717, 3 pages, 2012.

P. S. Chung, H. Park, and M. S. Jhon, “The static and dynamic responses of binary mixture perfluoropolyether lubricant films- molecular structural effects,” IEEE Transactions on Magnetics, vol. 45, no. 10, pp. 3644–3647, 2009.

Q. Guo, P. S. Chung, M. S. Jhon, and H. J. Choi, “Nano-rheology of single unentangled polymeric lubricant films,” Macromolecular Theory and Simulations, vol. 17, no. 9, pp. 454–459, 2008.

H. Chen and M. S. Jhon, “Relationship between surface coverage and end group functionality of molecularly thin perfluoropolyether films,” Journal of Applied Physics, vol. 103, no. 7, Article ID 07F536, 3 pages, 2008.

P. S. Chung, H. Chen, and M. S. Jhon, “Molecular dynamics simulation of binary mixture lubricant films,” Journal of Applied Physics, vol. 103, no. 7, Article ID 07F526, 3 pages, 2008.

Q. Guo, P. S. Chung, H. Chen, and M. S. Jhon, “Molecular rheology of perfluoropolyether lubricant via nonequilibrium molecular dynamics simulation,” Journal of Applied Physics, vol. 99, no. 8, Article ID 08N105, 3 pages, 2006.

R. M. Balabin, “Communications: intramolecular basis set superposition error as a measure of basis set incompleteness: can one reach the basis set limit without extrapolation?” Journal of Chemical Physics, vol. 132, no. 21, Article ID 211103, 4 pages, 2010.

R. Smith, P. S. Chung, S. H. Vemuri, L. T. Biegler, and M. S. Jhon, “Atomistically tuning lubricant adhesion on carbon overcoat surface,” IEEE Transactions on Magnetics, vol. 48, no. 11, pp. 4273–4276, 2012.

S. Ghosh, W. Bao, D. L. Nika et al., “Dimensional crossover of thermal transport in few-layer graphene,” Nature Materials, vol. 9, no. 7, pp. 555–558, 2010.

W. R. Zhong, M. P. Zhang, B. Q. Ai, and D. Q. Zheng, “Chirality and thickness-dependent thermal conductivity of few-layer graphene: a molecular dynamics study,” Applied Physics Letters, vol. 98, no. 11, Article ID 113107, 3 pages, 2011.

A. Agarwal, L. T. Biegler, and S. E. Zitney, “Simulation and optimization of pressure swing adsorption systems using reduced-order modeling,” Industrial and Engineering Chemistry Research, vol. 48, no. 5, pp. 2327–2343, 2009.

W. G. Noid, J. W. Chu, G. S. Ayton et al., “The multiscale coarse-graining method. I. A rigorous bridge between atomistic and coarse-grained models,” Journal of Chemical Physics, vol. 128, no. 24, Article ID 244114, 11 pages, 2008.

S. Izumisawa and M. S. Jhon, “Stability analysis of ultra-thin lubricant films with chain-end functional groups,” Tribology Letters, vol. 12, no. 1, pp. 75–81, 2002.

R. J. Waltman, G. W. Tyndall, and J. Pacansky, “Computer-modeling study of the interactions of Zdol with amorphous carbon surfaces,” Langmuir, vol. 15, no. 19, pp. 6470–6483, 1999.

P. H. Kasai and A. M. Spool, “Z-DOL and carbon overcoat: bonding mechanism,” IEEE Transactions on Magnetics, vol. 37, no. 2, pp. 929–933, 2001.

R. J. Waltman, D. J. Pocker, and G. W. Tyndall, “Studies on the interactions between ZDOL perfluoropolyether lubricant and the carbon overcoat of rigid magnetic media,” Tribology Letters, vol. 4, no. 3-4, pp. 267–275, 1998.

T. Aoyagi, J. Takimoto, and M. Doi, “Molecular dynamics study of polymer melt confined between walls,” Journal of Chemical Physics, vol. 115, no. 1, pp. 552–559, 2001.

M. S. Jhon, G. Sekhon, and R. Armstrong, “The response of polymer molecules in a flow,” in Advances in Chemical Physics, I. Prigogine and S. A. Rice, Eds., vol. 66, pp. 153–211, John Wiley, New York, NY, USA, 1987.

Q. Guo, S. Izumisawa, D. M. Phillips, and M. S. Jhon, “Surface morphology and molecular conformation for ultrathin lubricant films with functional end groups,” Journal of Applied Physics, vol. 93, no. 10, pp. 8707–8709, 2003.

K. Kremer and G. S. Grest, “Dynamics of entangled linear polymer melts:  a molecular-dynamics simulation,” Journal of Chemical Physics, vol. 92, no. 8, pp. 5057–5086, 1990.

S. Izumisawa and M. S. Jhon, “Stability analysis and molecular simulation of nanoscale lubricant films with chain-end functional groups,” Journal of Applied Physics, vol. 91, no. 10, p. 7583, 2002.

S. J. Vinay, D. M. Phillips, Y. S. Lee et al., “Simulation of ultrathin lubricant films spreading over various carbon surfaces,” Journal of Applied Physics, vol. 87, no. 9, pp. 6164–6166, 2000.

X. Ma, J. Gui, L. Smoliar et al., “Spreading of perfluoropolyalkylether films on amorphous carbon surfaces,” Journal of Chemical Physics, vol. 110, no. 6, pp. 3129–3137, 1999.

Y. T. Hsia, Q. Guo, S. Izumisawa, and M. S. Jhon, “The dynamic behavior of ultrathin lubricant films,” Microsystem Technologies, vol. 11, no. 8-10, pp. 881–886, 2005.

G. W. Tyndall, R. J. Waltman, and D. J. Pocker, “Concerning the interactions between Zdol perfluoropolyether lubricant and an amorphous-nitrogenated carbon surface,” Langmuir, vol. 14, no. 26, pp. 7527–7536, 1998.

G. W. Tyndall, T. E. Karis, and M. S. Jhon, “Spreading profiles of molecularly thin perfluoropolyether films,” Tribology Transactions, vol. 42, no. 3, pp. 463–470, 1999.

K. Binder, Monte Carlo and Molecular Dynamics Simulations in Polymer Science, Oxford University Press, Oxford, UK, 1995.

Q. Guo, S. Izumisawa, M. S. Jhon, and Y. T. Hsia, “Transport properties of nanoscale lubricant films,” IEEE Transactions on Magnetics, vol. 40, no. 4, pp. 3177–3179, 2004.

P. G. de Gennes, Scaling Concepts in Polymer Physics, Cornell University Press, Ithaca, NY, USA, 1979.

X.-C. Guo, B. Knigge, B. Marchon, R. J. Waltman, M. Carter, and J. Burns, “Multidentate functionalized lubricant for ultralow head/disk spacing in a disk drive,” Journal of Applied Physics, vol. 100, no. 4, Article ID 044306, 2006.

X.-C. Guo, T. Karis, H. Deng, Q. Dai, J. Burns, and R. J. Waltman, “Fomblin multidentate lubricants for ultra-low magnetic spacing,” IEEE Transactions on Magnetics, vol. 42, no. 10, pp. 2504–2506, 2006.

H. Tani, K. Sakamoto, and N. Tagawa, “Conformation of ultrathin pfpe lubricants with different structure on magnetic disksdirect observation and MD simulation,” IEEE Transactions on Magnetics, vol. 45, no. 11, pp. 5050–5054, 2009.

R. N. Kono, S. Izumisawa, M. S. Jhon, C. A. Kim, and H. J. Choi, “Rheology of perfluoropolyether lubricants,” IEEE Transactions on Magnetics, vol. 37, no. 4 I, pp. 1827–1829, 2001.

H. J. Choi, Q. Guo, P. S. Chung, and M. S. Jhon, “Molecular rheology of perfluoropolyether lubricant via nonequilibrium molecular dynamics simulation,” IEEE Transactions on Magnetics, vol. 43, no. 2, pp. 903–905, 2007.

C. M. Mate and V. J. Novotny, “Molecular conformation and disjoining pressure of polymeric liquid films,” The Journal of Chemical Physics, vol. 94, no. 12, pp. 8420–8427, 1991.

S. Izumisawa and M. S. Jhon, “Calculation of disjoining pressure for lubricant films via molecular simulation,” IEEE Transactions on Magnetics, vol. 42, no. 10, pp. 2543–2545, 2006.

Q. Guo, The static and dynamic properties of nano-structured thin oligomeric films [Ph.D. thesis], Carnegie Mellon University, Pittsburgh, Pa, USA, 2006.

Y. I. Jhon, S. E. Zhu, J. H. Ahn, and M. S. Jhon, “The mechanical responses of tilted and non-tilted grain boundaries in graphene,” Carbon, vol. 50, no. 10, pp. 3708–3716, 2012.

W. T. Kim, D. Kim, S. H. Vemuri, S. C. Khang, P. S. Chung, and M. S. Jhon, “Multicomponent gas mixture air bearing modeling via lattice Boltzmann method,” Journal of Applied Physics, vol. 109, no. 7, Article ID 07B759, 3 pages, 2011.

H. M. Kim, J. H. Kang, and M. S. Jhon, “Hydro-kinetic approach in non-Newtonian lattice Boltzmann flow simulation,” Journal of the Korean Physical Society, vol. 58, article 444, 2011.

D. Kim, H. M. Kim, M. S. Jhon, S. J. Vinay, and J. Buchanan, “A characteristic non-reflecting boundary treatment in lattice Boltzmann method,” Chinese Physics Letters, vol. 25, no. 6, article 1964, 2008.

X. Ma, C. L. Bauer, M. S. Jhon, J. Gui, and B. Marchon, “Monte Carlo simulations of liquid spreading on a solid surface: effect of end-group functionality,” Physical Review E, vol. 60, no. 5, pp. 5795–5801, 1999.

D. M. Phillips, A. S. Khair, and M. S. Jhon, “Mathematical simulation of ultra-thin polymeric film spreading dynamics,” IEEE Transactions on Magnetics, vol. 37, no. 4 I, pp. 1866–1868, 2001.

D. M. Phillips and M. S. Jhon, “Dynamic simulation of nanoscale lubricant films,” Journal of Applied Physics, vol. 91, no. 10, p. 7577, 2002.

W. T. Kim, M. S. Jhon, Y. Zhou, I. Staroselsky, and H. Chen, “Nanoscale air bearing modeling via lattice Boltzmann method,” Journal of Applied Physics, vol. 97, no. 10, Article ID 10P304, 3 pages, 2005.

Q. Guo, L. Li, Y.-T. Hsia, and M. S. Jhon, “A spreading study of lubricant films via optical surface analyzer and molecular dynamics,” IEEE Transactions on Magnetics, vol. 42, no. 10, pp. 2528–2530, 2006.

Q. Guo, L. Li, Y.-T. Hsia, and M. S. Jhon, “Stability analysis of ultrathin lubricant films via surface energy measurements and molecular dynamics simulations,” Journal of Applied Physics, vol. 97, no. 10, Article ID 10P302, 3 pages, 2005.

H. M. Kim, D. Kim, W. T. Kim, P. S. Chung, and M. S. Jhon, “Langmuir slip model for air bearing simulation using the Lattice Boltzmann method,” IEEE Transactions on Magnetics, vol. 43, no. 6, pp. 2244–2246, 2007.