Tribology: Friction and Lubrication of Hematite

Tribology Friction Boundary lubrication Hematite ReaxFF Molecular dynamics NEMD

This tutorial will teach you how to:

• build a lubricated contact from a hydroxylated hematite slab,

• prepare a hexadecane/stearic-acid mixture between two surfaces,

• define fixed, thermostatted, and driven regions,

• equilibrate the contact under a normal load,

• perform a non-equilibrium molecular dynamics (NEMD) sliding simulation, and

• estimate an apparent friction coefficient from the tangential and normal reaction forces.

Note: The hands-on calculation is deliberately small and short so that the complete workflow can be explored on a laptop. It is intended to teach the setup. The production results shown below use larger systems, longer trajectories, several loads, and independent replicas.

Production results

The workflow described in this tutorial was extended to compare three contacts: dry hydroxylated hematite, hematite lubricated by hexadecane, and hematite lubricated by a hexadecane mixture containing 25 mol% stearic acid.

For each trajectory, an apparent friction coefficient was calculated from the time-averaged reaction forces on the driven upper slab,

\[\mu_{\mathrm{app}} = \frac{|\langle F_t \rangle|}{|\langle F_z \rangle|}.\]

At the higher loads, 10, 12.5, and 15 nN, the normal reaction was reasonably well controlled. The three load-specific values were therefore averaged to give one simple descriptor for each formulation.

Table 13 Mean apparent friction coefficient from the high-load production calculations.

Formulation	μ at 10 / 12.5 / 15 nN	Mean apparent μ
Dry hematite	4.391 / 3.632 / 2.976	3.666
Hexadecane	0.342 / 0.402 / 0.297	0.347
Hexadecane + 25 mol% stearic acid	0.317 / 0.294 / 0.340	0.317

Note

The main result is the large friction reduction produced by the lubricant film. The difference between hexadecane and the stearic-acid mixture is much smaller: adding 25 mol% stearic acid lowers the mean apparent friction by about 9% in this dataset. One individual load reverses this ordering, so this should be interpreted as a modest average trend rather than a guaranteed reduction at every load.

The values above are apparent friction coefficients for the simulated conditions, not universal material constants. In particular, nanoscale adhesion, finite contact size, confinement, sliding velocity, and surface structure can all influence the ratio.

Fig. 41 Mean apparent friction coefficient for dry hematite, hexadecane, and hexadecane + 25 mol% stearic acid. The bars are averages over 10, 12.5, and 15 nN.

Introduction

Boundary lubrication occurs when two solid surfaces are separated by a molecular film whose structure and interaction with the surfaces control the shear response. Long-chain hydrocarbons can reduce direct solid-solid contact, while polar organic friction modifiers such as stearic acid can interact more strongly with an oxide surface.

In a tribology NEMD calculation, a normal load presses the two surfaces together and one surface is driven laterally. The resulting tangential and normal reaction forces can be used to estimate friction. Here the upper hematite slab slides along the y direction while it remains free to respond to the applied load along z.

The complete production workflow used the following sequence:

• pre-optimize the hydroxylated hematite slab,

• build two opposing surfaces and insert the lubricant,

• pre-optimize the complete interface,

• equilibrate under the normal load without sliding, and

• start the NEMD production run from the equilibrated coordinates and velocities.

The GUI calculation below follows the same sequence but uses only one mixture, one load, and one trajectory. The production mixture contains 25 mol% stearic acid, whereas the compact GUI example uses eight hexadecane molecules and two stearic-acid molecules, corresponding to 20 mol% stearic acid.

Part A: Prepare the Hydroxylated Hematite Slab

Import the prepared slab

The starting structure is Fe2O3_0001_H2O_rectangular.xyz. It contains 72 atoms and a rectangular two-dimensional cell of approximately 5.05 × 8.75 Å.

The slab was prepared beforehand from an α-Fe₂O₃ (0001) surface. Conceptually, the procedure consists of cutting the (0001) surface from bulk hematite, replicating the in-plane unit cell to obtain a rectangular supercell, hydroxylating the exposed surfaces, and relaxing the resulting slab. Surface construction and hydroxylation are not the focus of this tutorial, so we start from the prepared XYZ file. See Crystals and Surfaces for instructions on building crystals and slabs.

1. Open AMSinput: SCM → New Input

2. Choose File → Import Coordinates (System) and select Fe2O3_0001_H2O_rectangular.xyz

3. Confirm that the system has Slab periodicity

4. Open Model → Lattice and verify the rectangular in-plane cell

5. Rotate to a side view with Cmd/Ctrl+2 and inspect the two exposed hydroxylated surfaces

Pre-optimize the slab

Before building the contact, perform a short ReaxFF optimization of the small slab. In the production workflow this step is used to relax the prepared surface and its in-plane lattice before replication.

1. Switch to the ReaxFF engine: →

2. Set Force field to FeOCHCl-ox.ff

3. Set Task to Geometry Optimization

4. Go to Details → Geometry Optimization

5. Set Optimize lattice to Yes

6. Go to Model → Geometry Constraints and PES Scan

7. Constrain the cell to remain rectangular by checking the XY component of the Freeze Strain matrix

8. Save the job as hematite_slab_relax and run it

9. When the job finishes, update AMSinput to the optimized geometry

Fig. 42 Hydroxylated hematite slab after geometry optimization.

Part B: Build the Lubricated Contact

Create a small two-surface model

The production calculations use a 4 × 4 supercell. Here we use a 2 × 2 replication of the small slab. The final system contains 1,088 atoms and leaves enough space for the lubricant mixture.

1. Choose Edit → Crystal → Generate Super Cell

2. Repeat the slab 2 times along the first lattice vector and 2 times along the second

3. Select all atoms of the replicated slab, then copy and paste them to create the upper slab

4. Translate the upper slab along +z until the two inner surfaces are separated by about 40 Å. To do this, choose Edit → Tune Geometry. In the Translate atoms, direction section, set Step to 56 and z to 1, then apply the translation by clicking the small triangle pointing to the right

5. Wrap the atoms in the supercell with Edit → Crystals → Map Atoms To (0..1)

Inspect the contact and check that the lowest H atom in the upper slab is approximately 40 Šfrom the highest H atom in the lower slab. The resulting tutorial surface area is approximately 175 Ų.

Fig. 43 Two opposing hematite slabs separated by an approximately 40 Å gap.

Define the solid regions

The solid is assigned to three regions. The lowest 25% of the lower slab is fixed, while the remaining 75% of the lower slab and the entire upper slab are thermostatted. The entire upper slab also receives the normal load and sliding velocity.

Create the following three regions with Model → Regions:

Table 14 Solid regions used in the tutorial.

Region	Approximate part of slab	Role
bottom_fixed	lowest 25% of bottom slab	fixed mechanical support
thermo	remaining 75% of bottom slab and entire top slab	thermostatted solid atoms
top_drive	entire top slab	receives the total load and net sliding velocity

1. In the side view, select the lowest 25% of the lower slab and create bottom_fixed

2. Select the remaining 75% of the lower slab together with the entire upper slab and create thermo

3. Select the entire upper slab and create top_drive

4. Select each region in the Regions window and visually verify its membership

Note

The overlap between thermo and top_drive is intentional: every atom in the upper slab belongs to both regions. AMS automatically excludes the net velocity imposed by Apply velocity from thermostatting, so the thermostat acts on the thermal motion of the driven slab rather than opposing its prescribed sliding motion.

This simplified scheme thermostats all mobile solid atoms, including those at the contact. It should therefore be applied consistently to every formulation being compared, because changing the thermostatted region can affect the measured friction response.

Warning

bottom_fixed must not overlap thermo or top_drive. Only the overlap between thermo and top_drive is intended.

Fig. 44 Three-region assignment for the two hematite slabs. The upper slab belongs to both thermo and top_drive.

You can save the job as stearic_preopt.

Add the hexadecane/stearic-acid lubricant

Pack an 80:20 molecular mixture between the two surfaces. For the small tutorial cell, it is simpler to specify the molecule counts directly: eight hexadecane molecules and two stearic-acid molecules.

First, switch to three-dimensional periodicity to enable packing in a 3D box.

1. In the Main panel, switch Periodicity to Bulk

2. Set the third component of the c lattice vector in Model → Lattice to 76 Å

This defines the packing volume between the two inner surfaces. Then pack the molecules.

1. Open Builders → Packed Molecules

2. Choose the option to specify the number of molecules

3. Add hexadecane with SMILES CCCCCCCCCCCCCCCC and set the number of molecules to 8

4. Add a second component with the + button

5. Add stearic acid with SMILES CCCCCCCCCCCCCCCCCC(=O)O and set the number of molecules to 2

6. Assign both molecular components to a region named fluid

7. Generate the molecules and inspect the structure from the side

8. If necessary, regenerate the molecules until you obtain an acceptable configuration

Note

Ten molecules are a compromise for a lightweight tutorial. Eight hexadecane molecules and two stearic-acid molecules give exactly 20 mol% stearic acid, but the initial density is only approximate. The production workflow uses a larger surface and density-based Packmol packing.

Fig. 45 Eight hexadecane molecules and two stearic-acid molecules between the hematite surfaces.

Close the builder and switch Periodicity back to Slab. Save the job and continue to the optimization.

Pre-optimize the complete interface

Packmol creates a reasonable starting arrangement, but local contacts can still be strained. Perform a short geometry optimization before starting molecular dynamics.

1. Keep the ReaxFF force field FeOCHCl-ox.ff

2. Set Task to Geometry Optimization

3. In Details → Geometry Optimization, uncheck Optimize lattice

4. Set Convergence to Custom

5. Set the Gradient and Energy convergence criteria to 0.05 eV/Å and 0.0001 eV, respectively

6. Select the region bottom_fixed with Select → Select Region → bottom_fixed

7. Open Model → Geometry Constraints and PES Scan and add a fixed-position constraint for the selected region

8. Save the job and run it

9. Inspect the optimization in SCM → Movie

10. After confirming that the final geometry is reasonable, update AMSinput to that geometry

Tip

If the optimization immediately develops very large forces because two molecules overlap, return to the packing step and generate a better initial configuration. If the optimization reaches the maximum number of iterations, inspect the final geometry carefully. It may still be suitable for load equilibration if it no longer contains unphysical contacts.

Part C: Equilibrate Under the Normal Load

The lubricant should first equilibrate under the normal load without sliding. This allows the film thickness and normal reaction to settle before shear is introduced.

The production series used 10, 12.5, and 15 nN on a much larger surface. To keep approximately the same normal pressure in the 2 × 2 tutorial cell, use a total load of about 3 nN. This is close to the pressure of the 12.5 nN production calculation after scaling by contact area.

Set up the load-equilibration MD

1. Update AMSinput to the pre-optimized interface

2. Set Task to Molecular Dynamics

3. In Model → MD, set Time step to 0.25 fs

4. Set Number of steps to 40000 (10 ps)

5. Set the initial velocities to Random at 300 K

6. Keep bottom_fixed fixed

7. In Thermostat, add one NHC thermostat at 300 K with Tau = 100 fs and set Region to thermo

8. Open Model → MD → NEMD Velocities and Forces and add an Apply force entry

9. Set Region to top_drive and Force to 0.0 0.0 -3e-9 N

10. Keep Per atom unchecked so that 3 nN is the total load on top_drive

11. In Model → MD → MD trajectory options, confirm that Write velocities is enabled

12. Save the job as stearic_equil and run it

Note

There is no lateral driving velocity in this stage. The upper slab must remain free to move in z under the normal load.

Inspect the equilibration

1. Open the running or completed job in SCM → Movie

2. View the contact from the side (Cmd/Ctrl+1) and follow the inter-slab separation

3. Model → Add Graph and plot the temperature with MD Properties → Temperature

4. Select one atom from each slab with approximately the same x and y coordinates

5. Model → Add Graph and plot the separation with Properties → Distance, Angle, Dihedral

6. Confirm that the lubricant remains confined between the surfaces

7. Check that the gap is not undergoing an obvious monotonic collapse or expansion at the end of the run

For a research calculation, continue the equilibration until the film thickness and normal reaction are statistically stationary. The 10 ps run here is only for learning the setup.

Fig. 46 Energy, temperature, and inter-slab distance during the equilibration run.

Part D: Run the Sliding NEMD Calculation

The NEMD stage starts from the equilibrated coordinates and velocities, keeps the same normal load, and imposes a lateral velocity on top_drive.

Set the sliding velocity

1. When stearic_equil finishes, update AMSinput to its final geometry

2. Keep Task set to Molecular Dynamics

3. Set Initial velocities to From File and select the ams.rkf file from stearic_equil

4. Keep the bottom_fixed constraint, the 300 K NHC thermostat on thermo, and the force in the negative z direction on top_drive

5. Confirm that Time step is 0.25 fs, then set Number of steps to 100000 (25 ps)

6. In Model → MD → NEMD Velocities and Forces, add an Apply velocity entry

7. Set Region to top_drive

8. Set Velocity to 0.0 0.001 0.0 Angstrom/fs

9. Set Components to XY

Only the net x and y velocity components of top_drive are controlled. The upper slab can therefore still move along z in response to the normal force. Although the upper slab also belongs to thermo, AMS excludes the imposed net velocity from thermostatting.

Set the trajectory output and run

1. Open Model → MD → MD trajectory options

2. Enable Write molecules and Write velocities

3. Optionally enable Write engine gradients if you want to analyze per-atom forces

4. Save the job as stearic_nemd

5. Start the calculation with File → Run

Note

0.001 Angstrom/fs corresponds to 100 m/s, much faster than a macroscopic sliding experiment. High velocities are commonly used in atomistic NEMD because accessible simulation times are short. The resulting friction should be interpreted for the simulated conditions.

Fig. 47 Normal-force and sliding-velocity settings for the NEMD calculation.

Part E: Inspect the Forces and Estimate Friction

Inspect the sliding trajectory

1. Open stearic_nemd in SCM → Movie

2. Play the trajectory and confirm that the upper slab moves along +y

3. Confirm that the bottom support remains fixed while the upper slab can move along z

4. Check that the lubricant stays between the two surfaces

5. Inspect the end of the trajectory for obvious drift in the gap or contact structure

Inspect the tangential and normal reaction forces

When Apply velocity is active, AMS automatically writes the instantaneous, mean, and standard deviation of the net engine force on the driven region to the trajectory. These quantities exclude the externally applied force. The tangential component opposes sliding, while the z component provides the normal reaction used in the apparent friction ratio.

1. In AMSmovie, choose Graph → Add Graph

2. Inspect the MD mean engine-force components for top_drive

3. Plot the y component, which is the tangential reaction for this calculation

4. Plot the z component, which gives the normal reaction; after equilibration, its magnitude should be close to the applied 3 nN load

5. Average a late-time portion of both signals; for this short tutorial the second half is a reasonable first inspection

6. Calculate the apparent friction coefficient from the averaged components

For this run, averaging over the second half of the trajectory gives \(\langle F_y \rangle = -0.624454\ \mathrm{nN}\) and \(\langle F_z \rangle = 3.10325\ \mathrm{nN}\). Therefore,

\[\mu_{\mathrm{app}} = \frac{|\langle F_y \rangle|}{|\langle F_z \rangle|} = \frac{0.624454}{3.10325} = 0.201.\]

Important

Do not expect the 25 ps, 1,088-atom, single-trajectory tutorial result to reproduce the production result for the stearic-acid formulation. The force signal is noisy, and friction is a statistical observable. The purpose of this calculation is to verify the workflow and identify the quantities that must be converged in a production study.

Fig. 48 Tangential and normal reaction-force traces for the driven upper slab.

Part F: Extend the Calculation to Production Quality

Once the small calculation behaves sensibly, increase the spatial and statistical sampling. The production calculations summarized at the beginning used the same physical workflow with larger systems and longer runs.

A useful progression is:

• increase the surface to approximately 20 × 35 Å,

• equilibrate each load for approximately 100 ps,

• run approximately 1 ns of sliding NEMD,

• repeat at 10, 12.5, and 15 nN,

• use at least two independent replicas per load,

• compare dry hematite, hexadecane, and hexadecane + 25 mol% stearic acid, and

• calculate friction from a statistically settled late-time interval before combining the high-load values.

The automated Python workflow follows these steps and is useful once the GUI calculation has been understood.

Download the production workflow: tribology_case.py

Summary

In this tutorial, you built a small boundary-lubricated hematite contact, added a hexadecane/stearic-acid film, equilibrated it under a normal load, and imposed sliding with ReaxFF NEMD. The same workflow can be extended to larger cells and longer simulations to compare lubricant formulations.

For the production calculations used here, hexadecane reduces the mean apparent friction by roughly an order of magnitude relative to dry hematite, while adding 25 mol% stearic acid gives a smaller additional reduction on average.