Atomic insights into thickness-dependent deformation mechanism and mechanical properties of Ag/PMMA ultra-thin nanofilms

Abstract

In this work, the nanoindentations on bilayer composite nanofilms composed of metal Ag and polymer PMMA were simulated using molecular dynamics. The effects of the thickness of Ag and PMMA on the elastic moduli of the composite films were analyzed from Hertz contact theory, dislocation evolution and atomic migration. The results show that the maximum penetration depth that the Hertz model could well describe is about 6 Å, and this limiting value is almost independent on the film thickness. The deformation mode of the Ag films gradually changes from bending mode to indentation mode with an increase in Ag thickness, which improves the elastic modulus of the composite films. The rule of mixtures could give a theoretical prediction about the elastic modulus of the composite film close to the nanoindentation, and Hertz theory could also be used as long as the thickness of Ag films exceeded a certain value. The introduction of a PMMA layer impedes the development of dislocation in the Ag layer and improves the elastic limit of the composite films. This work provides an important basis for experimentally measuring the overall elastic modulus of metal/polymer composite film based on nanoindentation or extracting the elastic modulus of metal film from the overall indentation response of the composite film.

---

1. Introduction

During the manufacturing of some functional devices, such as solar cells nanosensors,¹ micro-electromechanical systems (MEMS)² and biological medicine devices,^3,4 film deposition processes are commonly involved, generating complex stacks of thin films.⁵ These stacks of thin films may combine different families of materials with dissimilar properties, including metals, polymers, and inorganic nonmetallic materials,^6,7 and utilize the excellent material properties or functionalities of each component. For example, Liu et al.¹ developed multi-functional nano-wrinkled composite films composed of Al/Ag/Si/ITO/PDMS for enhancing energy conversion from sunlight and raindrops. Li et al.⁸ designed a flexible LED array device made of Si-microwires and ZnO nanofilm, which are mechanically flexible, stable and lightweight, as well as having excellent optoelectronic properties. Zeng et al.⁹ irreversibly synthesized an in-plane-bonded two-dimensional polymer material (2D polyaramid) under solution phase conditions, and it shows great potential in lightweight nanofiltration composites due to its high modulus, stiffness and yield strength, and even ultra-low gas permeability. Shakil et al.¹⁰ found that compared with room temperature, the hardness, reduced modulus and yield strength of GO-SF-CNC (graphene oxide, silk fibroin and cellulose nanocrystal) nanocomposites dramatically increase at 80 °C without any significant degradation in mechanical or creep behavior, so they can be used in flexible electronics, thermal interface materials, and energy storage devices with long-term exposure to different temperatures. Gu et al.¹¹ presented stretchable conductors made from multiple layers of gradient assembled polyurethane (GAP) comprising gold nanoparticles capable of self-assembly under strain, which delivers stable power with a high-rate capability under high strain. Bae et al.¹² developed an unconventional 3D hybrid nanocomposite by embedding a uniform 3D alumina nanoshell into an epoxy-siloxane molecular hybrid (ESMH) matrix, obtaining metal-like hardness, high flexibility and transparency, and fine wear impact resistance, which might be useful for protective films for flexible optoelectronic devices. Therefore, by rationally designing nanolayered heterogeneous structures, composite films could provide multifunctional properties that a single material cannot provide.

Although the methods of preparing and characterizing thin films have made great progress in recent years, current thin film technology is still facing a series of multidisciplinary challenges in order to meet diverse performance requirements.^13–16 One of these challenges is measurement of the material properties of the thin film stack due to its composite response.^17,18 The yield or fracture of micro–nano thin film devices when subjected to external loads are strongly related to a combination of mechanical properties and geometries, especially the thickness of each film layer.^19–21

For nanoscale films or coating materials, the relevant literature^22,23 has shown that it is extremely challenging to experimentally measure the mechanical properties of nanoscale films by traditional tensile and bending testing methods. Nanoindentation has been widely adopted for the characterization of mechanical properties of thin film materials, where one typically starts by measuring the indentation load–displacement relationship of a selected thin film, and then fits the result to an analytical model to extract the elastic modulus and strength of the material. In their important applications, thin film materials are almost always attached to a substrate via vdW interaction. The overall indentation load–displacement relationship of a thin film on a substrate could be measured via instrumented nanoindentation. Then the reduced modulus could be obtained by either nonlinearly fitting the initial elastic loading segment to the Hertz elastic contact theory,²⁴ or by linearly fitting the slope of the initial unloading segment using the Oliver–Pharr method.²⁵

With regard to composite films, it is also challenging to conduct tensile or bending tests experimentally because nanofilms are commonly deposited on rigid substrates. And the theoretical estimation approach may be limited by the unknown mechanical properties or unknown thicknesses of each component. By contrast, nanoindentation has been proven to be a more feasible experimental technique to determine the elastic modulus of multilayer heterogeneous structures. For example, Tang et al.²⁶ determined the elastic modulus of seven layers of alternating Al and SiC deposited on an Si substrate using the nanoindentation method. The variation of the elastic modulus, measured from the indentation unloading, as a function of layer thickness and indentation depth has been investigated. Zhang et al.²⁷ systematically investigated the length scale dependent hardness and elastic modulus of nanostructured Cu/Cr multilayer films from 5 nm to 250 nm by using a nanoindentation technique with a Berkovich indenter. Zeng et al.²⁸ used the nanoindentation method based on atomic force microscopy to determine the Young's modulus and yield strength of two-dimensional metal–organic framework (MOF) nanosheet materials with atom-sized pores and thickness ranging from 10 nm (dozens of monolayers) to 400 nm (stacked multilayer). Tang et al.²⁹ studied the hardness and elastic modulus of Au/NiCr/Ta soft metal multilayers by using a nanoindentation continuous stiffness measurement technique. Ghidelli et al.³⁰ measured the elastic modulus and residual stress of Au-TiW bilayer films by a nanoindentation method.

The Hertz model was successfully used to describe frictionless, non-adhesive contacts between an indenter and the surface of a half-space material.³⁵ However, for many cases of interest (microelectronics, coatings, electromechanical devices, composites, etc.), the Hertz elastic contact model has been extended to study the mechanical response of multilayer composite films or heterogenous materials.^26,31–34 For example, Zhang et al.³² studied the elastic modulus of porous metallic glass using nanoindentation molecular dynamic simulations combined with the classical Hertz elastic contact model. Feng et al.³¹ studied the nanoindentation process of WC-Co composite materials by molecular dynamics simulation, evaluated the hardness through the Oliver–Pharr method and estimated the elastic modulus according to Hertz theory. Humood et al.³³ fabricated bilayer and quadri-layer film stacks by layer-by-layer assembly and performed nanoindentation experiments to obtain the elastic modulus, hardness, and elastic recovery of the PEI/MMT composite film using the Oliver–Pharr method.

Correlating the indentation response to the mechanical properties of the composite nanofilms is of great importance. However, when thin film stacks are nanometer scale in thickness, the microstructure of materials cannot be ignored because atoms are essentially arranged in a discrete manner. Therefore, previous theoretical or numerical studies of the nanoindentation on thin film stacks based on the continuum hypothesis may not be able to fundamentally explain the nanofilm deformation mechanism. In addition, for nanoindentation experiments, it is currently very difficult and expensive if not impossible to prepare nanoscale multi-layer films with a smooth surface and uniform thickness. At the same time, it is also difficult to ensure the accuracy and repeatability of the indentation results due to the existence of complex factors, such as adhesion, friction and electrostatic interaction, in the experiment.³⁴ Therefore, molecular dynamics (MD) simulation, which can not only modify the calculation model and simulation conditions at will, but also track the microstructural evolution of thin films, provides an effective approach to study the micro–nano deformation mechanism of thin films at the nanoscale and interpret the data obtained from nanoindentation experiments. The goal of this work is to use MD simulations to investigate the deformation mechanism of a metal/polymer bilayer under an indentation load and to help interpret the indentation data.

In this paper, Ag and poly(methyl methacrylate) (PMMA) polymer, both of which are commonly used in a variety of science and engineering fields, were selected as the components to form the composite nanofilms. MD simulations of the nanoindentations on the bilayer composite films (Ag as the upper layer and PMMA as the lower layer) supported by the Si substrate were performed in order to investigate the effect of the metal/polymer interface and film thickness combinations on the Young's modulus of the composite film and dislocation in the Ag layer. In Section 2, we introduce the MD simulation method, including the potential function, model setup and numerical model parameters used for the nanoindentation simulation. In Section 3, after verifying the model parameters, the effects of Ag film thickness and PMMA film thickness on the mechanical properties of the composite films are studied by Hertz contact theory, dislocation analysis and atomic migration distribution mapping. The rule of mixtures is used to provide a theoretical estimation of the overall elastic modulus. Finally, in Section 4, we summarize the results and present the conclusions and prospects.

2. Method

2.1 Potential function and model setup

In this study, a series of molecular dynamics (MD) simulations were performed to simulate the nanoindentations on metal Ag, polymer PMMA, and their composite films through Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) code.³⁵ All these thin films are rigidly supported by the silicon substrate. Four different atomic potentials, the embedded-atom method (EAM),³⁶ the consistent valence force field (CVFF),³⁷ Tersoff/mod,³⁸ and Lennard–Jones (LJ 12-6) potentials are adopted in this study. The EAM atomic potentials, which are obtained by fitting the experimental and first-principles data and can well reproduce the thermal and mechanical properties of Ag, were used to describe the interaction between Ag atoms. The CVFF potential, which is suitable for many organic and polymer systems, was used to precisely describe the molecular interactions within PMMA in this study. The total potential energy includes the bonding interaction terms (bond stretching energy, bond angle bending energy, dihedral angle twist energy, energy of cross terms) and intermolecular non-bonding interactions (van der Waals electrostatic interactions). The Tersoff/mod potentials were optimized by Kumagai et al.³⁸ on the basis of the Tersoff potential function for the elastic constant and melting point of diamond-structured silicon, which was applied for simulation of the silicon substrate in this study. The van der Waals interaction between adjacent layers of multilayer films was described by the well-established LJ 12-6 potential, whose expression is:where σ is the collision diameter at which E_LJ = 0 and ε is the depth of the potential well; the cutoff radius is uniformly set to 12.5 Å. The σ and ε parameters between each atom are listed in Table 1, while the σ and ε parameters between different atoms are fitted by the following rules:³⁹. These potentials have been widely adopted to study the mechanical properties of a variety of materials under different loading conditions, such as uniaxial tension, nanoindentation, and shock compression,^40–42 which provides a necessary guarantee of the reliability of the simulation results.

Table 1 Lennard–Jones parameters of the non-bonded interactions

Atom type	Ag	PMMA					Si
		H	C^1–3	O	C	C′
σ (kcal mol^−1)	4.56	0.038	0.039	0.288	0.16	0.148	0.401
ε (Å)	2.955	2.45	3.875	2.86	3.475	3.617	3.826
Reference	43	37					44

---

Considering the complexity of the multilayer structure and the accuracy of the calculation results, two strategies were selected for modeling the regular crystal structures (FCC Ag and diamond Si) and amorphous polymer PMMA, respectively. For the case of Ag thin films, the face-centered cubic unit cell with a lattice constant of 4.09 Å was replicated to generate Ag films with a surface area of 24 × 24 nm and thicknesses of 2, 4, 6 and 8 nm. The X, Y, Z axes correspond to the crystal orientations of [100], [010] and [001], respectively, and the upper surface along the Z axis is subjected to nanoindentation in the simulations. For the case of the Si substrate, the diamond-structured Si unit cell with a lattice constant of 5.43 Å was replicated to generate the Si substrate with a volume of 24 × 24 × 3 nm. The Si substrate was divided into two layers: the lowest bottom layer of the Si substrate called the fixed layer is a group of atoms with a thickness of 1 nm and its role is to prevent the model from moving during relaxation and indentation processes; while the remaining atoms in the upper part of the Si substrate is called the thermostatic layer, whose temperature is kept at a constant 300 K using a temp/rescale method to dissipate the Joule heat and balance the temperature difference between the Newton layer and the fixed layer.

The initial configuration of the polymer PMMA in the composite film was created with the aid of the Amorphous Cell module in Materials Studio software, where the confine layer function was used to construct two-dimensional films (without periodicity in the Z-direction). Firstly, a PMMA single chain composed of 50 monomers was generated; then 150, 300, 450 and 600 PMMA chains were packed into the simulation cells with thicknesses of 2, 4, 6, 8 nm, respectively. The initial density was set to 1.10 g cm⁻³. However, due to the overlapping or hole phenomenon in the atomic distribution of the initial configuration of the polymer film, such a model is usually in a high-energy state and cannot be used for subsequent nanoindentation simulations. Therefore, following previous work,^45,46 in LAMMPS an annealing process was used to eliminate the inhomogeneity of atomic distribution within the PMMA model. The annealing process was as follows: (1) virtual and rigid walls were inserted parallel to the X–Y plane on the upper and lower surfaces of the PMMA, where the distance between the upper and lower bounds in the Z-direction was the same as the initial thickness of the film and remained unchanged. This operation can keep the roughness of the PMMA surface relatively small during the relaxation process, thereby minimizing the impact on the contact zero point of the indenter and the contact area in the subsequent nanoindentation simulation. (2) The PMMA films were relaxed for 50 ps with a temperature of 450 K, which is higher than the glass transition temperature of PMMA (T_g ≈ 380 K⁴⁷). Therefore, after sufficient time for the relaxation process, molecular chains with sufficient mobility achieve a stable and uniform distribution eliminating the defects. (3) In the next step, the PMMA films were cooled from 450 K to 300 K at a constant rate of 3 K ps. (4) Then the PMMA films were relaxed for 50 ps at 300 K. All the above simulations were performed under the NVT ensemble. (5) Finally, the ensemble was switched from NVT to NPT and maintained at 300 K and 0 Pa for 50 ps. After all the initial configurations of PMMA had undergone a complete annealing process, their densities remained between 1.06 and 1.09 g cm⁻³ and the pressure was maintained between 100 and 200 atm during the last 30 ps of relaxation time. It can be deduced that the PMMA thin films were in an equilibrium state.

After each layer of film was prepared, films of different thickness were combined into composite films in the sequence Ag, PMMA and Si substrate. As an example, Fig. 1 illustrates the numerical model of the composite film system with 4 nm of Ag and 4 nm of PMMA. Periodic boundary conditions were imposed in the X- and Y-directions, while the vacuum region at the upper boundary was used to simulate the free surface, and the lower boundary was supported by the fixed Si substrate. The relaxation process of the composite films was carried out at 300 K in the NVT ensemble for 1 ns (time step of 1 fs) to eliminate the free-floating of the film in the Z-direction caused by the introduction of the vacuum layer. After sufficient relaxation of all models, nanoindentation simulation could be performed for the multilayer film system under the same NVT ensemble at 300 K.

Fig. 1 Schematic illustration of the composite film model with 4 nm of Ag film and 4 nm of PMMA film.

2.2 Indentation setup

Many previous studies^48–50 have reported that indentation results are greatly affected by the indenter radius, indentation speed and penetration depth. Therefore, the simulation conditions must be chosen carefully to exclude the influence of irrelevant variables. A virtual and rigid spherical indenter was used in all nanoindentation simulations, the advantage of which is that the virtual indenter has a completely repulsive and smooth surface, so the adsorption and friction phenomena that often occur in experiments can be ignored, simplifying the complexity of the problem. The force applied by the indenter during nanoindentation is defined as follows:where K is the specified force constant and is chosen to be 10 eV Å⁻³, r is the distance of an atom to the indenter center, and R is the indenter radius. In this study, the thickness of the composite films was between 4 nm to 16 nm, so the indenter radius was set as 8 nm, which ensures enough atomic contacts during the indentation process so the indentation results can more intrinsically reflect the mechanical properties of the sample itself. The penetration depth was limited to 15% of the total thickness of the composite films.⁵⁰

The complete nanoindentation simulation process involves three stages: loading at a constant velocity, holding at the deepest penetration depth for a short time, and reverse unloading at the same velocity. The loading velocity was kept at 5 m s⁻¹, which lies well below the speed of sound and the dislocation glide velocity,⁵¹ so the entire nanoindentation simulation is equivalent to a quasi-static process, and the time-dependent creep phenomenon is not involved in this study. From the subsequent simulation results, it can also be seen that the mean value of the indentation force during the holding process is approximately the same as that at the end stage of loading and the initial stage of unloading, which is strong proof of a quasi-static process.

2.3 Elastic modulus

It is well known that Hertz elastic contact theory has been widely used to study the mechanical contact behavior of metals, polymers and other materials in the initial elastic stage, which can obtain the most basic mechanical parameters, such as elastic modulus and hardness. Meanwhile, compared with the traditional Oliver−Pharr method based on the unloading part of the indentation load–depth curve, some studies^48,52 also proved that the Hertz model could more accurately reflect the mechanical properties of the measured materials.

In the case of a spherical indenter, the relation between indentation force P and penetration depth h can be expressed as:

where E_r is the reduced modulus and R is the radius of the indenter. Then the reduced modulus obtained through eqn (3) has the following relationship with the mechanical properties of indenter and films:where E_s and E_i are the Young's modulus of the film and the indenter, respectively; and ν_s and ν_iare the Poisson's ratio of the film and the indenter, respectively. As mentioned above, the indenter is rigid and non-deformable, which means that E_i tends to infinity. So, the second term in eqn (4) also tends to zero and can be ignored. In addition, according to previous research, ν_Ag and ν_PMMA are 0.34 and 0.337,^53,54 respectively. Therefore, in this work the Poisson's ratio of the composite films is set to 0.34 for simplicity, and eqn (4) is finally simplified to:

E_s = 0.8844E_r (5)

In this study, firstly, the nanoindentation load–displacement curves within the initial elastic contact segment are nonlinearly fitted according to eqn (3), and the reduced modulus E_r is used to obtained the elastic modulus E_N of the composite film by combining with eqn (5). Since the bilayer nanofilms are anisotropic, the modulus measured by the nanoindentation technique is nominally the Young's modulus perpendicular to the plane of the film. However, given the geometry of the indenter, the modulus when actually measured is a weighted average over all directions. Although the indentation technique cannot provide the components of the modulus along two orthogonal directions of layered film, for many practical applications such as flexible electronics, knowledge of the “averaged” Young's modulus or the relative elasticity to some known material is also desirable.

3. Results and discussion

3.1 Numerical model validation

Prior to the nanoindentation simulation for the Ag/PMMA composite films, in order to verify the model parameters of the two film materials, nanoindentation simulations of pure Ag and pure PMMA films of different thicknesses were performed and compared with previous simulations and experimental results.

Taking an Ag film with a thickness of 8 nm as an example, Fig. 2(a) shows the curve of the indentation force as a function of the penetration depth during the loading and unloading processes. It can be seen from this figure that the indentation force increases nonlinearly and monotonically with increasing penetration depth from the contact zero point to a penetration depth of 6 Å. The loading segment of the indentation force–depth curve is almost consistent with the result predicted using the Hertz elastic contact theory. The simulation data in the elastic stage was fitted with the Hertz model to obtain the reduced modulus of the Ag film; then the reduced modulus was substituted into eqn (5) to calculate the Young's modulus, which was found to be 92.29 GPa for Ag. From Fig. 2(a) we also found that when the penetration depth reached 6 Å, the indentation force curve gradually deviated from the Hertz contact theory, and showed a noisy oscillation feature. When the indentation depth reached 12 Å, which is 15% of the total thickness, the Ag film was unloaded. As can be seen from Fig. 2(a), the indentation force gradually decreased with the withdrawal movement of the spherical indenter, and the fluctuation level of the indentation force also gradually decreased. When the indenter left the surface of the film, the residue penetration depth was 3.6 Å, leaving a permanent and irreversible plastic dent. Nanoindentation load–displacement curves of pure Ag films with different thicknesses are shown in Fig. 2(b). The nanoindentation Young's modulus E_N was obtained by nonlinear fitting of the curve in the 0–4 Å contact range to the Hertz contact model, and the calculated elastic moduli are listed in Table 2. As can be seen from Table 2, the elastic modulus E_N of Ag film under loading by nanoindentation shows little difference with variation in thickness and is close to the values reported in previous literature.

Fig. 2 Metallic Ag film: (a) simulation results of nanoindentation with a thickness of 8 nm and the fitted curve of the Hertzian elastic contact model; (b) nanoindentation load–displacement curves with different thicknesses.

Table 2 Elastic moduli E_N of pure Ag film and pure PMMA film obtained from nanoindentation simulation

	Thickness (nm)	E N(GPa)	References
Ag	4	86.93	63–83.8,^55 79–98,^53 85–112,^56 91 (ref. 57)
	6	89.59
	8	92.29
PMMA	4	4.99	5.52,^54 5.54,^58 4.4,^59 4.9,^60 5.07–5.2 (ref. 61)
	6	5.61
	8	4.9

---

Next, we performed similar nanoindentation simulations on PMMA film with a thickness of 8 nm, and the results are shown in Fig. 3(a). The Hertz model can still well describe the elastic contact response at the initial stage when indenting the PMMA polymer material. The Young's modulus of PMMA film with a thickness of 8 nm was 4.9 GPa calculated by the same procedure described above. Unlike metal materials, the elastoplastic transition of polymers is not obvious and is difficult to detect due to the absence of a dislocation phenomenon. Since in our following simulations PMMA is only located in the lower layer and does not directly come into contact with the spherical indenter, its surface will not undergo a large degree of distortion. Therefore, the plastic deformation of PMMA films barely occurs and will not be considered in this study. Nanoindentation load–displacement curves of pure PMMA films of different thicknesses are shown in Fig. 3(b), and the corresponding nanoindentation elastic moduli are given in Table 2. It was found that the nanoindentation Young's modulus E_N of PMMA film is nearly independent of film thickness.

Fig. 3 Polymer PMMA film: (a) simulation results of nanoindentation with a thickness of 8 nm and the fitted curve of the Hertzian elastic contact model; (b) nanoindentation load–displacement curves with different thicknesses.

Therefore, it can be inferred that using the simulation parameters adopted in this study accurately reproduces the mechanical properties of Ag and PMMA nanofilm, laying the foundation for the following simulations of composite films.

3.2 Mechanical properties of Ag/PMMA composite films

In this section, to investigate the relationship between mechanical properties and thickness of Ag/PMMA composite films, nanoindentation simulations were performed with the same thickness ratio of metal Ag and polymer PMMA, while the thicknesses of the two components of the film varied in an arithmetic sequence. To be specific, the thickness ratio was fixed at 1 : 1 for all four groups of models and the thickness combinations of these models included 2 + 2 nm (2 nm of Ag and 2 nm of PMMA), 4 + 4 nm (4 nm of Ag and 4 nm of PMMA), 6 + 6 nm (6 nm of Ag and 6 nm of PMMA), and 8 + 8 nm (8 nm of Ag and 8 nm of PMMA).

Fig. 4(a) shows the complete curves of nanoindentation force vs. penetration depth for different thickness combinations of bilayer composite films. As shown in Fig. 4(a), the indentation force of all composite films increased rapidly with increasing penetration depth in the initial elastic stage, except for the 2 + 2 model. For this case when the indenter moved down to a certain depth, the curve was no longer smooth and acute and irregular fluctuations occurred, and the curve was serrated until the end of the loading process. For the four models the corresponding indentation forces at the maximum penetration depths were recorded as 45 nN, 112 nN, 269 nN and 383 nN, respectively. Then the indenter was withdrawn at a speed of −5 m s⁻¹ and the nanoindentation simulation turned into the unloading phase. The unloading trend of each curve was also similar to the initial elastic loading process. Finally, after the indenter completely left the film surface, the indentation force decreased to zero. The residual penetration depths were 1.1 Å, 3.2 Å, 5.9 Å and 10.2 Å, respectively, which indicate that permanent plastic deformation had already occurred in the composite films and the thicker the composite film is, the greater the residual depth is.

Fig. 4 Ag/PMMA composite film with a thickness ratio of 1 : 1 where the unit of thickness is nm. (a) Nanoindentation load–displacement curves with different thicknesses. (b) An enlarged view of the initial elastic stage, where the simulated data and the Hertz fitting curve are represented by hollow scatter points and dotted lines with different colors, respectively.

Fig. 4(b) shows that, although the difference in the elastic modulus between metal Ag and polymer PMMA is considerable, the Hertz elastic contact theory can still well describe the force–displacement relationship of composite films in the initial elastic stage. At the same time, it is found that the curves of the initial elastic stage of each model become closer to each other with an increase in the film thickness. In particular, the indentation force–penetration depth curves of the 6 + 6 and 8 + 8 models coincide almost completely. In addition, with an increase in penetration depth, we also find that the 2 + 2 model first deviates from Hertz theory at 4.3 Å, but the other three curves all deviate from Hertz theory near an indentation depth of 6 Å. This is because the total thickness of the 2 + 2 model is only 4 nm, so the plastic deformation occurs earlier at a lower penetration depth. Based on the above results, we can infer that the maximum penetration depth that Hertz theory could predict seems to be independent of film thickness. The Hertz theory cannot accurately describe the whole process of elastic deformation of thin films and the deviation point cannot be used as a reliable marker of the elastic–plastic transition.

The rule of mixtures^62,63 has been adopted to provide a theoretical estimation, which is expressed as:

E_Mix = E_AgV_Ag + E_PMMAV_PMMA (6)

where E_Ag and E_PMMA are the elastic moduli of Ag and PMMA, respectively; V_Ag and V_PMMA are the volume fractions of Ag and PMMA, respectively; and E_Mix is the elastic modulus of the composite film from mixture theory. The values of E_Ag and E_PMMA are taken from the literature as 88 GPa⁵³ and 5.2 GPa,⁶¹ respectively.

Table 3 summarizes and compares the elastic moduli obtained by the nanoindentation method and by the rule of mixtures. By analyzing the data in Table 3, we found that as long as the thickness of Ag layer is larger than 4 nm, no matter what the thickness of PMMA layer is (6 + 2, 6 + 4, 6 + 6 or 8 + 8), the theory of mixtures could give an overall elastic modulus for the composite film close to that of the nanoindentation method with a discrepancy <10%. On the other hand, if the Ag film is too thin, for example 2 nm (2 + 2 and 2 + 6), the overall elastic moduli of the composite film obtained by theory of mixtures and the nanoindentation method have a much larger discrepancy of 36–40%. Therefore, it could be concluded that the Young's modulus of composite films can be obtained by the nanoindentation method using Hertz elastic contact theory, and the thickness of Ag films is a key factor in the accuracy of the results. The rule of mixtures given by eqn (6) follows the Voigt approximation (equal strain assumption). We further calculated the prediction results using the rule of mixtures based on the Reuss approximation (equal stress assumption), as shown in eqn (7):

E⁻¹_{Mix Reuss} = E⁻¹_AgV_Ag + E⁻¹_PMMAV_PMMA (7)

Table 3 Elastic moduli of Ag/PMMA composite films obtained from nanoindentation E_N and theory of mixtures E_Mix obtained by Voigt and Reuss approximation, respectively

h Ag + hPMMA (nm)	E N (GPa)	E Mix – Voigt (GPa)	Discrepancy – Voigt	E Mix – Reuss (GPa)	Discrepancy – Reuss
2 + 2	28.16	46.6	40%	9.82	−65%
4 + 4	39.83		15%		−75%
6 + 6	47.85		−3%		−79%
8 + 8	48.77		−5%		−80%
6 + 2	61.95	67.3	8%	17.94	−71%
6 + 4	50.06	54.88	9%	11.94	−76%
2 + 6	16.65	25.9	36%	6.80	−59%
4 + 6	33.48	38.32	13%	8.34	−75%

---

From Table 3 it is found that predictions given by the Reuss approximation have a larger discrepancy than that of the Voigt approximation. J. A. Ruud et al.⁶⁴ measured the elastic modulus of Ag–Ni multilayered film using nanoindentation, and they also found the rule of mixtures based on the Voigt approximation predicts well the elastic modulus of multilayered film.

The data in Table 3 also shows that the rule of mixtures based on the Voigt approximation (eqn (6)) could well predict the overall indentation response of the composite film as long as the Ag film thickness was larger than 4 nm, which could be used to design indentation experiments to extract the elastic modulus of the metal or polymer component from the indentation response of the composite film. For example, if the overall elastic modulus of a composite film has been obtained from a nanoindentation experiment and the elastic modulus of the polymer component is known, then the elastic modulus of the metal component could be calculated using the rule of mixtures as:

E_Metal = (E_Mix − E_PolyV_Poly)/V_Metal (8)

3.3 Influence of PMMA thickness on the mechanical properties of composite films

Based on the above research results, we know that the elastic modulus of Ag/PMMA composite film is significantly correlated with its thickness. In Section 3.2 the thickness of PMMA and Ag were changed simultaneously while maintaining a thickness ratio of 1 : 1. Next in order to separate the influence of PMMA film thickness on the Young's modulus of the composite films and the dislocation evolution of the Ag film with FCC structure, nanoindentation simulations were performed on composite films composed of Ag film with a thickness of 6 nm and PMMA films with thicknesses of 0 nm, 2 nm, 4 nm, and 6 nm.

Fig. 5(a) shows the curves of indentation force–penetration depth for composite films with thickness combinations of 6 + 0 (6 nm of Ag and no PMMA), 6 + 2 (6 nm of Ag and 2 nm of PMMA), 6 + 4 (6 nm of Ag and 4 nm of PMMA) and 6 + 6 (6 nm of Ag and 6 nm of PMMA). The calculation results of the elastic modulus are listed in Table 3. As can be seen in Fig. 5(b), the Hertz model can precisely describe the mechanical response of the initial elastic stage. With increasing PMMA thickness from 0 to 6 nm, the simulated data deviated from the Hertz contact model when the penetration depth reached 5 Å, 5.2 Å, 6 Å and 6 Å, respectively. This indicates that the PMMA film thickness has little effect on the prediction limit of the Hertz contact model.

Fig. 5 Ag/PMMA composite film with thickness combinations of h_Ag fixed at 6 nm and h_PMMA of respectively 2, 4, 6 nm. (a) Nanoindentation load–displacement curves with different thicknesses. (b) An enlarged view of the initial elastic stage, where the simulated data and the Hertz fitting curve are represented by hollow scatter points and dotted lines with different colors, respectively.

In addition, we also found that due to the introduction of PMMA, compared with the 6 + 0 model, the indentation force (the slope in Fig. 5(b)) of the other three models decreased to a large extent, and the former was almost twice as large as the latter three under the same penetration depth. As a consequence, compared with pure Ag film, the elastic modulus of the 6 + 2, 6 + 4, 6 + 6 composite films decreased by 30.9%, 44.1%, and 46.6%, respectively, indicating that the introduction of PMMA greatly reduces the stiffness of the composite film. Fig. 6 shows the distribution diagram of shear strain of 6 + 0, 6 + 2, 6 + 4 and 6 + 6 composite films at 6 Å indentation depth. When the penetration depth was 6 Å, the shear strain of 6 + 0 composite film (pure Ag film) near the contact area was significantly greater than that of 6 + 2, 6 + 4 or 6 + 6 composite films, which explains the higher elastic modulus of pure Ag film. For the three composite films with PMMA, their elastic responses showed little difference. As the indenter continued to load, the indentation force dropped and fluctuated to different degrees in all curves, and the curves of 6 + 4 and 6 + 6 composite films fluctuated significantly in the later stage of loading. Finally, with completion of the unloading process, we can see that the residual depth gradually increased with an increase in PMMA thickness, and were 1.8 Å, 2.3 Å, 4.2 Å and 5.9 Å, respectively, manifesting that the degree of plastic deformation is constantly increasing.

Fig. 6 Distribution diagrams of shear strain of 6 + 0, 6 + 2, 6 + 4 and 6 + 6 composite films at 6 Å indentation depth.

In order to explain the different mechanical responses of the combined films under nanoindentation loading conditions from the atomic level, we also extracted the volumetric strain distribution maps of composite films with different combinations under the same indentation depth in OVITO,⁶⁵ and the composite film model that had been well relaxed before indenter contact was selected as the reference configuration. As shown in Fig. 7, when the indentation depth is 8 Å, there is a large volumetric strain near the contact area of Ag film in the four groups of composite films, which is mainly caused by compression by the spherical indenter, and a small amount of accompanying tensile volumetric strain within the dislocation area. However, in the 6 + 2, 6 + 4 and 6 + 6 composite films, due to the addition of a PMMA film layer, a large amount of pressure is converted into elastic energy and stored in the underlying PMMA layer, which is manifested as obvious compressive volumetric strain between PMMA molecular chains, so the elastic modulus of PMMA shows a significant decrease compared with that of the 6 + 0 model. But at the same time, due to the staggered distribution between different molecular chains of polymer materials, some chain segments will also be subjected to tensile action, resulting in a staggered distribution mode in PMMA films, where a large amount of compressive strain is dominant and a small part of tensile strain is auxiliary.

Fig. 7 Distribution diagrams of volumetric strain of 6 + 0, 6 + 2, 6 + 4 and 6 + 6 composite films at 8 Å indentation depth. Blue and red in the color bars represent compressive and tensile strains, respectively.

Then, we identified the dislocation lines and defect structures near the contact area between the spherical indenter and Ag film with the help of the Dislocation Extraction Algorithm (DXA) in OVITO,⁶⁶ and further explored the elastic–plastic transition and plastic deformation mechanism of Ag in the composite films. For the sake of shortening the identification time and improving computational efficiency, we selected the cube region with a volume of 16 × 16 × 6 nm directly under the indenter as the plastic deformation zone, which contains the farthest distance that the dislocation line can extend during the nanoindentation process. This approach not only ensured the integrity of the dislocation structure during the entire deformation process but also reduced the computational burden.

Fig. 8(a)–(d) show the dislocation evolution and defect structure of four sets of models with thickness combinations of 6 + 0, 6 + 2, 6 + 4 and 6 + 6 at different penetration depths. The dislocation evolution undergoes three stages: the first dislocation nucleation, dislocation development towards the lower surface, dislocation explosion and dislocation configuration at maximum penetration depth. In each diagram of dislocation structure, the red atoms represent the HCP atoms on the stacking fault planes; the green and purple lines are the Shockley and stair-rod partial dislocation, respectively, while the red lines are the other types of dislocation. In order to clearly observe the dislocation evolution, we hide the perfect face-centered cubic (FCC) configuration and other atoms, including surface atoms and dislocation cores.

Fig. 8 Evolution of dislocation lines and defect structures in the plastic zone of metallic Ag film in composite films with thickness combinations of (a) 6 + 0, (b) 6 + 2, (c) 6 + 4 and (d) 6 + 6. The PMMA layer is hidden for clarity. Each group contains four snapshots of the dislocation structure: namely, the first dislocation nucleation, dislocation development towards the lower surface, the dislocation explosion, and dislocation configuration at maximum penetration depth.

In Fig. 8(a.1), (b.1), (c.1) and (d.1), we can observe that as the spherical indenter moves downward, Shockley dislocation with a Burgers vector of 1/6 〈112〉 is nucleated on the slip surface where the Ag atoms are most densely arranged. The dislocation lines of the four models all start from the upper surface and enter the Ag film at a depth of about 35 Å, and the stacking faults of the HCP structure are surrounded by these dislocation rings. However, the difference is that the initial penetration depths of the four models at which dislocation nucleation occurs are 4.5 Å, 5 Å, 7 Å and 8 Å, respectively. By comparing with the penetration depth at which the curve deviates from Hertz contact theory, we find that the penetration depth when the deviation occurs is quite close to the penetration depth of the dislocation nucleation in the 6 + 0 and 6 + 2 models, but with the increase in PMMA thickness, the gap between the two penetration depths is significantly widened in the 6 + 4 and 6 + 6 models. Many researchers believe that dislocation nucleation manifests as the first obvious force drop in the indentation force and displacement curve, which is also known as the pop-in phenomenon, and the appearance of dislocation nucleation usually means the beginning of plastic deformation.⁶⁷ But in this study, due to the addition of polymer PMMA, this criterion gradually fails with an increase in PMMA thickness. The reason is that during the initial loading process, partial elastic deformation is accumulated by the thicker PMMA film. The dislocation nucleation of metal Ag in the Ag/PMMA composite film can occur only under a larger indentation load or deeper penetration depth. In addition, the first decrease in force in these curves is fairly small, because in the process of nanoindentation by an indenter with a larger radius, the dislocation initially generated is discrete and unstable, and during unloading at this time the dislocation lines will disappear quickly, so it may not be enough to affect the variation tendency of indentation load and penetration depth curve. A similar phenomenon has also been reported in previous literature.^49,68,69 Therefore, our simulation results show that the bending of the indentation force–penetration depth curve for Ag/PMMA composite film is not a reliable indicator of dislocation nucleation. When the thickness of the PMMA film is close to that of the Ag film, the dislocation will be generated later than the drop point of the indentation force.

Then discrete dislocations develop towards the lower surface of the Ag film and cross each other, with a corresponding growth in defect atoms, as shown in Fig. 8(a.2), (b.2), (c.2) and (d.2). As the indenter continues to move down, the dislocation lines with different orientations in the plastic zone develop collectively and dramatically, with continuous emission and slip of Shockley partial dislocations until a dislocation explosion occurs, as shown in Fig. 8(a.3), (b.3), (c.3) and (d.3). After that, a stable dislocation network has formed in the plastic zone, and the dislocation develops steadily until the end of the loading process, as shown in Fig. 8(a.4), (b.4), (c.4) and (d.4).

Fig. 9(a) shows that the variation in the total length of the dislocation line in the plastic zone of Ag film varies with penetration depth during nanoindentation in four composite films with different thickness combinations. Firstly, it can be seen that the initial dislocation line length of all four composite films is about 140 Å, which illustrates that the increase of PMMA thickness does not change the incipient scale of dislocation nucleation. During the indentation, the dislocation lines of each model continued to develop with loading of the spherical indenter. However, we found that there are two rapid growth stages in each curve, which are marked I and II, respectively, in Fig. 9(a). For example, in stages a.I and a.II for the 6 + 0 model, the increasing slope of the dislocation curve during this period is significantly larger than that at other times. What is most noteworthy is that during the first rapid period of dislocation growth, the dislocation line penetrated the entire Ag film and developed towards the lower surface. Before that, the dislocation line only developed and interacted inside the Ag film. Subsequently, the development speed of the dislocation slowed down slightly, during which there also appeared a short plateau period. With further loading of the spherical indenter, a large amount of elastic energy is accumulated inside the composite films. Until the dislocation explosion occurred, most of the elastic energy in the Ag film overcame the previous dislocation structure and developed rapidly and abruptly. That is, the film entered the second rapid stage of dislocation growth, after which a stable and complex dislocation network had been formed. In addition, according to Fig. 5, it was also found that the indentation force–penetration depth curve shows large fluctuation during the rapid stage of dislocation development in each model, which indicates that only the collective and abundant dislocation motion in the Ag film can cause significant changes in the mechanical response of composite films. Finally, except for the rapid stage of dislocation development, it was also observed that the growth rates of the dislocation lines for the other three composite films (6 + 2, 6 + 4 and 6 + 6 models) are significantly lower than that of the pure Ag film (6 + 0 model), which also proves that the existence of PMMA can alleviate the rate of dislocation development within the Ag film and reduce the plastic deformation of Ag/PMMA composite film. By checking the residual strain when unloading the composite film before the occurrence of dislocations of the Ag layer (Fig. 9(b) and (c)), no obvious inelastic deformation was observed in either the Ag or PMMA layer, which indicates that the addition of polymer PMMA also improves the elastic limit of the composite film.

Fig. 9 (a) Variation of the total length of the dislocation lines in the plastic zone of Ag film with the penetration depth in composite films with different thickness combinations. (b) and (c) Atomic strain cloud image during loading and unloading of the 6 + 6 Ag/PMMA composite.

3.4 Influence of Ag thickness on mechanical properties of composite films

In this section, the effect of Ag layer thickness on the mechanical properties of Ag/PMMA composite films was investigated by tuning the thickness of the Ag films to 2 nm, 4 nm and 6 nm, while keeping the thickness of the PMMA film at 6 nm. The corresponding indentation force and penetration depth curves are shown in Fig. 10(a) and (b) shows the fitting result of the Hertz model in the elastic stage. Fig. 10(b) indicates that the relationship between indentation force and penetration depth is still in accordance with Hertz elastic contact theory in the initial elastic stage, and the maximum penetration depth applicable to the Hertz theory of the three composite films (2 + 6, 4 + 6 and 6 + 6) are 5.6 Å, 5.7 Å and 5.8 Å, respectively, indicating that the prediction limit of the Hertz model is almost independent of the thickness of the Ag film. Using eqn (5), the elastic moduli of the three composite films could be calculated as 16.65 GPa, 33.48 GPa and 47.85 GPa, respectively, as listed in Table 3. It could be concluded that the thickness of the Ag film will greatly affect the overall stiffness of Ag/PMMA composite films. With further loading of the spherical indenter, plastic deformation occurred in the composite films, and it was found that the 4 + 6 and 6 + 6 model curves show an obvious fluctuation phenomenon. As mentioned above, large-scale and collective dislocation activities may be responsible for the instability of the indentation force. Therefore, it is not difficult to speculate that the increase in Ag thickness leads to a sharp increase in dislocation level within its plastic zone, so the response curve shows a zigzag fluctuation. At the maximum indentation depth, the maximum indentation forces are 64.4 nN, 141 nN and 269 nN, respectively.

Fig. 10 Ag/PMMA composite film with thickness combinations of h_PMMA fixed at 6 nm and h_Ag being respectively 2, 4, 6 nm. (a). Nanoindentation load–displacement curves with different thicknesses. (b) The enlarged view of the initial elastic stage, where the simulated data and the Hertz fitting curve are represented by hollow scatter points and dotted lines with different colors, respectively.

Based on the above results, it is necessary to study and explain how the change in Ag thickness affects the mechanical properties of composite films from the microscopic deformation mechanism. Fig. 11 shows the atomic migrations of three composite films (2 + 6, 4 + 6 and 6 + 6) when the penetration depths were 5%, 10% and 15% of the total thickness. The atomic configuration before the indentation is used as the reference configuration. When the penetration depth was 5% of the total thickness of the composite films, the atomic displacement occurred only near the indentation center on the upper surface of the Ag film. However, when the indentation depth reached 10% of the total thickness of the composite film, it was found that there are two different deformation modes: in the 2 + 6 model, the larger atomic displacement is concentrated at the interface between the Ag film and the PMMA film. In addition, a larger bending deflection is also generated, leading to more serious stress concentration at the interface; while in the 6 + 6 model, the displacement of Ag film atoms is propagated over a wider range, and there is no bending deflection in the center of the Ag film. The overall PMMA surface is compressed by the Ag film; therefore, the homogeneous displacement of the upper surface of PMMA could be observed. As for the 4 + 6 model, its deformation mode is somewhere between the other two modes: namely, the coexistence of bending and indentation mechanisms. When the indentation depth reaches its maximum, the difference between these two deformation modes becomes more obvious. When the Ag film is relatively thin, the deformation mode of the Ag film in the composite film is bending. In the loading process of the spherical indenter, the central subsidence of the Ag film softens the overall stiffness of composite films, which contributes to the low hardness of the 2 + 6 composite film. When the Ag film is relatively thick, the Ag film in the composite film shows an indentation deformation mode. Only the upper surface of the Ag film is depressed under the action of the indenter, and at the same time, the pressure is evenly transferred to the upper surface of the PMMA, so the composite film shows greater hardness. The bending mode of thin Ag film on a PMMA substrate causes obstacles to determining the Young's modulus of Ag film directly from nanoindentation. A complex elasticity problem involving both the indentation and bending of the Ag film and the elastic deformation of the PMMA needs to be solved in order to extract the Young's modulus of Ag film from the overall indentation response. Niu et al.⁷⁰ selected polydimethylsiloxane (PDMS) as the substrate to measure the elastic modulus of graphene via AFM indentation. From the measured indentation response of graphene/PDMS, the elastic modulus of graphene can be determined via an inverse analysis using the finite element method.

Fig. 11 Atomic migration diagrams at different indentation depths for composite films with the thicknesses of PMMA layer fixed at 6 nm and the thicknesses of the Ag layer being 2 nm (2 + 6), 4 nm (4 + 6), or 6 nm (6 + 6).

Finally, the dislocation nucleation locations and defect evolution processes of two deformation modes in the Ag plastic zone were analyzed. As shown in Fig. 12, when the thickness of the Ag film was 2 nm or 4 nm, the dislocation nucleated from the interface between Ag and PMMA at penetration depths of 6.75 Å and 8.55 Å, respectively, and the nucleation scale was extremely small. However, with constant loading of the indenter, the dislocation line of the 2 + 6 model did not extend to the surface of the Ag film, and the dislocation movement only occurred in the lower part of the Ag film. While in the 4 + 6 model, the dislocation lines gradually crossed the Ag film from the bottom to the surface, forming a dislocation structure throughout the Ag film. For the 6 + 6 model, the dislocations initially nucleated from the upper surface of Ag film and then developed downward to form a dislocation structure crossing the Ag film thickness and forming a complex network.

Fig. 12 Initial dislocation nucleation locations and defect evolutions of composite films with the thickness of PMMA layer fixed at 6 nm and the thickness of the Ag layer being 2 nm (2 + 6), 4 nm (4 + 6) or 6 nm (6 + 6), respectively.

4. Conclusions

In this study, nanoindentation molecular dynamics simulation of Ag/PMMA composite films with different thickness combinations was carried out to investigate the influence of the thickness of each component on its mechanical properties. When the thickness ratio of the Ag layer to the PMMA layer was fixed at 1 : 1, the elastic modulus of the composite film increased with increased total thickness. And the elastic modulus of composite films tended to a constant value after the total thickness of the composite film reached 12 nm. The Hertz elastic contact theory can accurately describe the mechanical responses of the composite films at the initial elastic stage. The maximum penetration depth that the Hertz model could well describe is about 6 Å, and this limiting value is almost independent of film thickness. By comparing the elastic modulus obtained by nanoindentation simulation with the theoretical values given by the rule of mixtures, it is found that the rule of mixtures can well predict (discrepancy less than 10%) the overall indentation response of the composite film as long as the Ag film thickness reaches 6 nm. When the thickness of Ag layer is less than 4 nm, a large discrepancy (∼40%) between the elastic moduli obtained by the nanoindentation method and the rule of mixtures is observed.

The effects of changing Ag film thickness and PMMA film thickness on the mechanical properties of the composite films were studied by comparing and analyzing the composite films with different thickness ratios using strain distribution, dislocation evolution and atomic migration. The presence of the PMMA layer reduces the dislocation development speed of the Ag layer, but has less effect on the initial scale of dislocation nucleation. The indentation depth required for Ag dislocation nucleation increases with an increase in PMMA layer thickness, which indicates that the elastic limit of the composite film is improved. In addition, we also found that Ag layer in composite films experienced two stages of rapid dislocation development, during which the indentation force vs. penetration depth curve often fluctuated significantly; this phenomenon is particularly obvious when the thickness of PMMA is large.

Depending on the thickness of the Ag layer, there exist two deformation styles: bending-dominant mode and indentation-dominant mode. When the Ag layer is thinner, the center of Ag film undergoes a large deflection due to bending deformation, reducing the stiffness of the composite film. When the Ag layer is thick, the indentation deformation of the Ag surface dominates the deformation mode of the composite film, and the stiffness of the composite film increases. The points of deviation of the indentation force vs. penetration depth curve from Hertz elastic theory could not be used as a sign of dislocation nucleation of Ag/PMMA composite film, and the drastic fluctuation in the curve meant a large-scale collective movement of dislocation.

This study is expected to provide atomic-level insights into the mechanical contact properties of multi-layer heterogeneous nanocomposite films and to provide an important measurement technique of nanofilm properties based on the nanoindentation method for research into patterned nanoimprints, micro–nano processing technology, pressure-sensitive devices, flexible sensors, etc. The current work shows that when the Ag film is thick enough, the rule of mixtures matches well with the nanoindentation method. Thus, by rationally designing nanoindentation experiments, using the rule of mixtures and the indentation response of the composite film, it is possible to extract the elastic modulus of one component if the modulus of the other component is known. In future, more simulations could be conducted and the rule of mixtures could be modified to further extend its application range, including the parameter space where the metal film is much thinner.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

Thanks go to National Natural Science Foundation of China (NSFC-BRFFR grant no. 12111530281), BIT-BRFFR Joint Research Program (BITBLR2020018), State Key Laboratory of Explosion Science & Technology (QNKT20-01, YPJH20-6, JCRC18-01), and Beijing Institute of Technology Research Fund for financial support.

References

1. X. Liu, K. Cheng, P. Cui, H. Qi, H. Qin, G. Gu, W. Shang, S. Wang, G. Cheng and Z. Du, Nano Energy, 2019, 66, 104188.
2. Z. Dong, Q. He, D. Shen, Z. Gong, D. Zhang, W. Zhang, T. Ono and Y. Jiang, Microsyst. Nanoeng., 2023, 9, 31.
3. K. Wang, Y. Liu, H. Wang, Y. Liu, X. Yang and S. Sun, Biomater. Adv., 2022, 142, 213154.
4. D. T. Haynie, L. Zhang, W. Zhao and J. S. Rudra, Nanomed. Nanotechnol. Biol. Med., 2006, 2, 150–157.
5. H. Xu, M. K. Akbari, S. Kumar, F. Verpoort and S. Zhuiykov, Sens. Actuators, B, 2021, 331, 129403.
6. W. Ma and Z. Xing, Int. J. Heat Mass Transfer, 2013, 58, 639–651.
7. G. Jena and J. Philip, Prog. Org. Coat., 2022, 173, 107208.
8. X. Li, R. Liang, J. Tao, Z. Peng, Q. Xu, X. Han, X. Wang, C. Wang, J. Zhu and C. Pan, Acs Nano, 2017, 11, 3883.
9. Y. Zeng, P. Gordiichuk, T. Ichihara, G. Zhang, E. Sandoz-Rosado, E. D. Wetzel, J. Tresback, J. Yang, D. Kozawa and Z. Yang, Nature, 2022, 602, 91–95.
10. A. Shakil, S. Kim and A. A. Polycarpou, Adv. Mater. Interfaces, 2022, 9, 2101640.
11. M. Gu, W. J. Song, J. Hong, S. Y. Kim and B. S. Kim, Sci. Adv., 2019, 5, eaaw1879.
12. G. Bae, G. M. Choi, C. Ahn, S. M. Kim, W. Kim, Y. Choi, D. Park, D. Jang, J. W. Hong and S. M. Han, Adv. Funct. Mater., 2021, 31, 2010254.
13. S. Panda, S. Hajra, K. Mistewicz, P. In-na, M. Sahu, P. M. Rajaitha and H. J. Kim, Nano Energy, 2022, 107514.
14. A. R. Chowdhury, J. Jaksik, I. Hussain, R. Longoria, O. Faruque, F. Cesano, D. Scarano, J. Parsons and M. J. Uddin, Nano-Struct. Nano-Objects, 2019, 17, 148–184.
15. J. Guo, Y. Yu, L. Cai, Y. Wang and Y. Zhao, Mater. Today, 2021, 44.
16. C. Bavatharani, E. Muthusankar, S. M. Wabaidur, Z. A. Alothman, K. M. Alsheetan, M. mana AL-Anazy and D. Ragupathy, Synth. Met., 2021, 271, 116609.
17. S. Kang and T.-S. Kim, Appl. Surf. Sci., 2023, 608, 155113.
18. N. Gao, J. Yu, S. Chen, X. Xin and L. Zang, Synth. Met., 2021, 273, 116693.
19. Z. Lou, L. Wang, K. Jiang and G. Shen, Nano Today, 2019, 26, 176–198.
20. K. Ariga, M. Ito, T. Mori, S. Watanabe and J. Takeya, Nano Today, 2019, 28, 100762.
21. P. Liu, A. L. Cottrill, D. Kozawa, V. B. Koman, D. Parviz, A. T. Liu, J. Yang, T. Q. Tran, M. H. Wong and S. Wang, Nano Today, 2018, S1748013217305649.
22. G. Cao and H. Gao, Prog. Mater. Sci., 2019, 103, 558–595.
23. Y. Yang, X. Li, M. Wen, E. Hacopian, W. Chen, Y. Gong, J. Zhang, B. Li, W. Zhou and P. M. Ajayan, Adv. Mater., 2017, 29, 1604201.
24. X. Chen and Z. Yue, Int. J. Solids Struct., 2020, 191, 550–565.
25. G. Pharr and W. Oliver, MRS Bull., 1992, 17, 28–33.
26. G. Tang, D. Singh, Y.-L. Shen and N. Chawla, Mater. Sci. Eng., A, 2009, 502, 79–84.
27. J. Zhang, J. Niu, X. Zhang, P. Zhang, G. Liu, G. Zhang and J. Sun, Mater. Sci. Eng., A, 2012, 543, 139–144.
28. Z. Zeng, I. S. Flyagina and J.-C. Tan, Nanoscale Adv., 2020, 2, 5181–5191.
29. W. Tang, L. Shen and K. Xu, Thin Solid Films, 2005, 485, 72–76.
30. M. Ghidelli, M. Sebastiani, C. Collet and R. Guillemet, Mater. Des., 2016, 106, 436–445.
31. Q. Feng, X. Song, H. Xie, H. Wang, X. Liu and F. Yin, Mater. Des., 2017, 120, 193–203.
32. Y. Zhang, J. Xu, Y. Hu, S. Ding, W. Wu and R. Xia, Int. J. Mech. Sci., 2023, 249, 108254.
33. M. Humood, S. Chowdhury, Y. Song, P. Tzeng, J. C. Grunlan and A. A. Polycarpou, ACS Appl. Mater. Interfaces, 2016, 8, 11128–11138.
34. M. Hardiman, T. J. Vaughan and C. T. McCarthy, Compos. Struct., 2017, 180, 782–798.
35. S. Plimpton, J. Comput. Phys., 1995, 117, 1–19.
36. P. Williams, Y. Mishin and J. Hamilton, Modell. Simul. Mater. Sci. Eng., 2006, 14, 817.
37. P. Dauber-Osguthorpe, V. A. Roberts, D. J. Osguthorpe, J. Wolff, M. Genest and A. T. Hagler, Proteins: Struct., Funct., Bioinf., 1988, 4, 31–47.
38. T. Kumagai, S. Izumi, S. Hara and S. Sakai, Comput. Mater. Sci., 2007, 39, 457–464.
39. P. Peng, G. Liao, T. Shi, Z. Tang and Y. Gao, Appl. Surf. Sci., 2010, 256, 6284–6290.
40. K. Li, Z. Zhang, J. Yan, J. Yang and Z. Zhang, J. Mater. Sci. Technol., 2020, 57, 159–171.
41. J. I. Odujole and S. Desai, AIP Adv., 2020, 10, 095102.
42. S. Zhao, E. Hahn, B. Kad, B. A. Remington, C. Wehrenberg, E. M. Bringa and M. A. Meyers, Acta Mater., 2016, 103, 519–533.
43. H. Heinz, R. Vaia, B. Farmer and R. Naik, J. Phys. Chem. C, 2008, 112, 17281–17290.
44. M. Yaghoobi and G. Z. Voyiadjis, Comput. Mater. Sci., 2014, 95, 626–636.
45. L. Zhang and N. A. Fleck, Polymer, 2021, 226, 123748.
46. W. Xia, J. Song, D. D. Hsu and S. Keten, Macromolecules, 2016, 49, 3810–3817.
47. H.-N. Lee, R. A. Riggleman, J. J. de Pablo and M. Ediger, Macromolecules, 2009, 42, 4328–4336.
48. C. Peng, F. Zeng, B. Yuan and Y. Wang, Appl. Surf. Sci., 2019, 492, 579–590.
49. J. Knap and M. Ortiz, Phys. Rev. Lett., 2003, 90, 226102.
50. N. Tayebi, A. A. Polycarpou and T. F. Conry, J. Mater. Res., 2004, 19, 1791–1802.
51. J. Varillas, J. Očenášek, J. Torner and J. Alcalá, Acta Mater., 2021, 217, 117122.
52. R. Martinez and L. R. Xu, MRS Commun., 2014, 4, 89–93.
53. H. N. Pishkenari, F. S. Yousefi and A. Taghibakhshi, Mater. Res. Express, 2018, 6, 015020.
54. M. Muhammettursun, T. Bel, E. Kocacinar, E. Erman, F. B. Gul, A. Augousti and N. Baydogan, J. Appl. Polym. Sci., 2021, 138, 50689.
55. G. Cheng, X. Sun, Y. Wang, S. L. Tay and W. Gao, Surf. Coat. Technol., 2017, 310, 43–50.
56. Y. Cao, S. Allameh, D. Nankivil, S. Sethiaraj, T. Otiti and W. Soboyejo, Mater. Sci. Eng., A, 2006, 427, 232–240.
57. R. Zong, S. Wen, F. Zeng, Y. Gao, X. Li, B. He and F. Pan, Surf. Coat. Technol., 2007, 201, 7932–7938.
58. S.-P. Ju, H.-Y. Chen and C.-W. Shih, J. Nanopart. Res., 2018, 20, 1–17.
59. H. Li, J. Chen, Q. Chen and M. Liu, Mater. Des., 2021, 197, 109239.
60. L. Malekmotiei, A. Samadi-Dooki and G. Z. Voyiadjis, Macromolecules, 2015, 48, 5348–5357.
61. H. Lin, L. Lv and T. Jin, Appl. Sci., 2020, 10, 5454.
62. A. M. Díez-Pascual, M. A. Gómez-Fatou, F. Ania and A. Flores, Prog. Mater. Sci., 2015, 67, 1–94.
63. E. C. Romero, J. C. Osorio, R. T. Soto, A. H. Macías and M. G. Botero, Surf. Coat. Technol., 2019, 377, 124875.
64. J. A. Ruud, T. R. Jervis and F. Spaepen, J. Appl. Phys., 1994, 75, 4969–4974.
65. M. L. Falk and J. S. Langer, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1998, 57, 7192.
66. A. Stukowski, V. V. Bulatov and A. Arsenlis, Modell. Simul. Mater. Sci. Eng., 2012, 20, 085007.
67. C. L. Kelchner, S. Plimpton and J. Hamilton, Phys. Rev. B: Condens. Matter Mater. Phys., 1998, 58, 11085.
68. Y. Shao, X. Zhao, J. Li and S. Wang, Comput. Mater. Sci., 2013, 67, 346–352.
69. G. Ziegenhain, H. M. Urbassek and A. Hartmaier, J. Appl. Phys., 2010, 107, 061807.
70. T. Niu, G. Cao and C. Xiong, Int. J. Solids Struct., 2018, 132–133, 1–8.

---