Sage Journals: Discover world-class research

Abstract

The addition of fillers can significantly improve the mechanical behavior of polymers. The responsible mechanisms at the molecular level can be well assessed by particle-based simulation techniques, such as molecular dynamics. However, the high computational cost of these simulations prevents the study of macroscopic samples. Continuum-based approaches, particularly micromechanics, offer a more efficient alternative but require precise constitutive models for all constituents, which are usually unavailable at these small length scales. In this contribution, we derive a molecular-dynamics-informed constitutive law by employing a characterization strategy introduced in a previous publication. We choose silicon dioxide (silica) as an exemplary filler material used in polymer composites and perform uniaxial and shear deformation tests with molecular dynamics. The material exhibits elastoplastic behavior with a pronounced anisotropy. Based on the pseudo-experimental data, we calibrate an anisotropic elastic constitutive law and reproduce the material response for small strains accurately. The study validates the characterization strategy that facilitates the calibration of constitutive laws from molecular dynamics simulations. Furthermore, the obtained material model for coarse-grained silica forms the basis for future continuum-based investigations of polymer nanocomposites. In general, the presented transition from a fine-scale particle model to a coarse and computationally efficient continuum description adds to the body of knowledge of molecular science as well as the engineering community.

Keywords

Molecular dynamics continuum mechanics material characterization anisotropy elasticity

1. Introduction

The versatility of polymers can be further expanded by the addition of fillers. In particular, filler particles at the nanoscale have been found to significantly improve the polymers’ mechanical [1,2], thermal [3,4], electrical [5], and optical [6] behavior. Experimental studies show that the presence of filler particles leads to a significant modification of the properties of the surrounding polymer [7 –9]. This region of influence, referred to as the interphase, is considered to be responsible for the improved behavior of nanocomposites [10]. The experimental investigation of these nanocomposites is highly challenging due to the extremely small length scales and is further complicated by the large number of influencing factors [11] such as matrix-filler interaction [12] as well as filler size [13], shape [14], distribution [15], and content [16]. In order to identify promising polymer composite configurations more easily, empirical models and virtual testing strategies are increasingly resorted to.

1.1. Empirical models

Initially, attempts were made to estimate the properties of the polymer composites by empirically motivated, mostly analytical models, e.g., by Leidner and Woodhams [17] and Pukanszky [18]. These models, originally introduced for conventional composites, were extended to cover nanocomposites as well [19,20]. With these extensions, mechanical properties such as Young’s modulus [21] or tensile strength [22] can be derived from the material behavior of the individual constituents. However, the complex mechanisms within the interphases responsible for the observed reinforcements can only be represented in a highly simplified way in such empirical models, e.g., by the $B$ and $Z$ parameters in the work by Peng et al. [23].

1.2. Molecular dynamics studies

By way of contrast, particle-based numerical approaches such as molecular dynamics (MD) inherently capture the effects within the interphase since they are solely based on the particle interactions. Therefore, they are widely used to study nanocomposites: Lu et al. [24] investigate graphene-polymer nanocomposites by means of atomistic MD simulations and identify an interphase that exhibits 1.5 times the polymer matrix stiffness. In previous publications [25,26], we have studied the interphase effects in polystyrene-silica composites using an MD-continuum coupling method [27 –30]. This procedure enabled us to perform comparative pseudo-experiments on a large set of samples to identify elastoplastic effects occuring in the interphase [26]. Bocharova et al. [31] use coarse-grained (CG) MD simulations to complement their experimental findings and unravel the molecular mechanism responsible for observed reinforcements in poly(2-vinylpyridine)-silica nanocomposites. Despite these essential findings, MD simulations remain computationally too expensive to investigate engineering problems on macroscopic samples.

1.3. Micromechanics

To reduce the computational costs, micromechanics rely on efficient continuum material descriptions, while still resolving the composite at the level of the nanoparticles. To this end, representative volume elements (RVEs) that include all relevant geometrical information are considered and their mechanical behavior is homogenized to obtain the macroscopic properties of the composite [32]. However, the continuum constitutive laws of the constituents (filler and polymer matrix) necessary for the RVEs are usually not available at these extremely small length scales.

Therefore, the scope of this contribution is to present the calibration of a continuum mechanical constitutive law based on MD simulations of CG silica serving as an exemplary filler material (section 2.1). We employ a strategy introduced in our previous contribution [33] to characterize the material behavior using MD pseudo-experiments (section 2.2). To this end, we perform uniaxial tensile and simple shear tests (section 2.3) with time-proportional and time-periodic loadings. Within the limits of the small strains expected in a polymer nanocomposite, the examined silica exhibits linear elasticity with a pronounced anisotropy (sections 3.1 and 3.2). Based on these pseudo-experimental data, we calibrate linear elastic constitutive laws with different degrees of anisotropy, and thus accurately represent the material behavior (section 3.3).

This study is complementary to the characterization [33,34], and complex continuum modeling of polystyrene using a novel viscoelastic-viscoplastic constitutive law [35].

2. Model and simulation details

2.1. Coarse-grained silica

In this study, we rely on CG force fields introduced by Ghanbari et al. [36], who substitute each SiO₂ unit by one CG bead with an aggregated molar mass of 64.0 g mol−1. The CG interaction potentials were derived via iterative Boltzmann inversion [37] and comprise bond, angle, and nonbonded contributions. The obtained force fields were used and validated in the context of silica-enforced polystyrene in numerous studies [25,26,38 –40]. The CG beads form a crystalline structure due to the applied physically motivated force fields [36]. We obtain a triclinic unit cell with the lattice parameters $a = 7.8 Å$ , $b = 19.1 Å$ , and $c = 5.3 Å$ and the corresponding angles $α = 116.0 °$ , $β = 97.4 °$ , and $γ = 71.9 °$ . Since we are interested in the overall behavior of the CG silica, we are not considering finite spherical filler particles as in previous studies [26], but an infinite sample under periodic boundary conditions (PBC). Due to the triclinic structure, we set up a triclinic periodic simulation cell containing 1350 CG beads to ensure a box length of at least twice the cut-off radius in a single crystal arrangement, as shown in Figure 1.

Figure 1.

Simulation sample: three-dimensional periodic simulation cell with 1350 coarse-grained atoms in single crystal arrangement (triclinic lattice) with a mass density of 2.32 g cm⁻³.

Due to the regular structure of the single crystal, the resulting sample is not subjected to stochastic variations. Therefore, a statistical validation using different initial configurations, as carried out in comparable studies of polymers [33,34,41], is not necessary. After the initial placement of the beads, the sample undergoes an equilibration procedure (cf. section “Sample preparation details” in Appendix 1), resulting in a mass density of 2.32 g cm⁻³ which fits well with the literature [42].

2.2. Material classification scheme and characterization strategy

The mechanical response of a material to an external load can be assigned to one of the four paradigmatic material models: elasticity, plasticity, viscoelasticity, and viscoplasticity [43]. Merely, the presence of rate dependence and quasi-static hysteresis are decisive factors [44] as shown in Figure 2.

Figure 2.

Classification of material behavior based on rate dependence and quasi-static hysteresis into elasticity, plasticity, viscoelasticity, and viscoplasticity according to Haupt [44].

A practical approach to apply this classification scheme based on MD simulations and to quantify the individual contributions was presented as a characterization strategy in the work by Ries et al. [33]: First, time-proportional loading tests with different strain rates are performed to identify whether the material exhibits a rate dependence, i.e., viscous behavior. Second, time-periodic tests reveal whether the material is reversible by means of a quasi-static hysteresis. The application of quasi-static loads is not possible with MD due to the required computational capacities. However, we were able to show in the previous works [33,34] that the occurrence of a stable equilibrium hysteresis as a result of several loading–unloading cycles permits to evaluate the dissipation and thus the inelasticity of the material. Based on the obtained results, relaxation and creep tests may be carried out to decompose the inelasticity into its viscous and plastic components.

2.3. Applied deformation

In order to fully capture the constitutive behavior of the investigated material, both uniaxial and shear deformation, as depicted in Figure 3, need to be considered. In the case of an anisotropic, i.e., direction-dependent, material response, these principal load cases have to be applied for all three spatial directions. The silica considered here is typically embedded into a softer polymer matrix, e.g. in the form of spherical particles [11]. Due to the significant difference in the stiffness of the two phases, it can be assumed that the silica is only undergoing small deformations. In this case, linearized kinematics can be used, described by the Cauchy strain tensor

[ε]_{ij} = [\begin{matrix} ε_{xx} & ε_{xy} & ε_{xz} \\ ε_{yy} & ε_{yz} \\ sym . & ε_{zz} \end{matrix}] = [\begin{matrix} \frac{\partial u_{x}}{\partial x} & \frac{1}{2} [\frac{\partial u_{x}}{\partial y} + \frac{\partial u_{y}}{\partial x}] & \frac{1}{2} [\frac{\partial u_{x}}{\partial z} + \frac{\partial u_{z}}{\partial x}] \\ \frac{\partial u_{y}}{\partial y} & \frac{1}{2} [\frac{\partial u_{y}}{\partial z} + \frac{\partial u_{z}}{\partial y}] \\ sym . & \frac{\partial u_{z}}{\partial z} \end{matrix}]

(1)

relying on the assumptions of the infinitesimal strain theory [44].

Figure 3.

Applied deformation: (a) uniaxial tension with initial and deformed box length $L_{x}$ and $l_{x}$ , and atmospheric pressure $p_{atm}$ , and (b) simple shear with displacement $u_{x}$ , initial box length $L_{y}$ , and shear angle $γ_{xy}$ ; realized by canonical ensemble (NVT) and isothermal-isobaric ensemble (NPT).

2.3.1. Uniaxial tension

For uniaxial tension the tensile strain $ε_{ii}$ is given by

ε_{ii} = \frac{u_{i}}{L_{i}} = \frac{l_{i} - L_{i}}{L_{i}}

(2)

with displacement $u_{i}$ , initial and deformed box length in tensile direction $L_{i}$ and $l_{i}$ , respectively. In the lateral directions, the simulation cell is subjected to the atmospheric pressure $p_{atm}$ , as shown in Figure 3(a). This setup allows for unconstrained lateral contraction, comparable to real experiments. In MD, this is realized by combining the isothermal-isobaric ensemble (NPT) in lateral directions with a predefined deformation in tensile direction, which can be interpreted as a one-dimensional canonical ensemble (NVT) [33].

2.3.2. Simple shear

Simple shear, on the other hand, is induced by the parallel shifting of planes proportional to their position and is, therefore, by definition isochoric [46]. The associated shear strain $ε_{ij}$ is related to the displacement $u_{i}$ and the inital orthogonal box length $L_{j}$ , or the shear angle $γ_{ij}$ via

ε_{ij} = \frac{1}{2} \frac{u_{i}}{L_{j}} \approx \frac{1}{2} γ_{ij} .

(3)

The employed parallelepiped-shaped simulation cell necessitates the use of triclinic boxes (cf. Figure 3(b)). These feature additional angular degrees of freedom $γ_{xy}$ , $γ_{xz}$ , and $γ_{yz}$ , compared to conventional orthogonal cells. This simplifies the application of simple shear in MD: The requested shear angle is controlled directly according to equation (3), while all other degrees of freedom are fixed via an NVT ensemble.

2.3.3. Simulation details

We choose a Nosé-Hoover thermostat [47] and barostat [48] with a temperature and pressure coupling time of 500 fs and 500 fs, respectively. Contrary to our previous studies [33,34,49], the time step is kept constant at $Δ t = 5 fs$ . We apply the deformation continuously over the evaluated time $t$ either in time-proportional or time-periodic manner given by

ε (t) = \overset{\cdot}{ε} t and ε (t) = ε^{a} \sin (\frac{{\overset{\cdot}{ε}}^{\max}}{ε^{a}} t),

(4)

respectively. The former represents a linear ramp-up with a constant strain rate $\overset{\cdot}{ε}$ , while the latter follows a sinusoidal function with strain amplitude $ε^{a}$ and maximum strain rate ${\overset{\cdot}{ε}}^{\max}$ . Based on the findings of previous studies [33,34], we evaluate four load cycles for the time-periodic tests.

2.4. Stress measures

We evaluate the internal stress due to the external load using the virial stress tensor [50]

σ = - \frac{1}{V} \sum_{I = 1}^{n_{A}} \sum_{\binom{J = 1}{J > I}}^{n_{A}} f_{IJ} \otimes [r_{I} - r_{J}] - \frac{1}{V} \sum_{I = 1} m_{I} υ_{I} \otimes υ_{I}

(5)

with volume $V$ , number of atoms $n_{A}$ , atom position vectors $r_{I}$ and $r_{J}$ , atom velocity $υ_{I}$ , and atom mass $m_{I}$ . The pair-wise interaction force $f_{IJ}$ is obtained, as usual in MD, by derivation of the interaction potentials. Note that for the present systems, the stress is dominated by the virial components (minuend), while the kinetic contributions (subtrahend) only represent the samples’ thermal fluctuations, as demonstrated in section “Stress evaluation” in Appendix 1. According to Subramaniyan and Sun [51], the time-averaged virial stress coincides with the Cauchy stress commonly used in continuum mechanics. This hypothesis is, however, controversially discussed [52,53]. Since the stress obtained in MD simulations is usually subject to strong fluctuations, a Savitzky–Golay filter [54] is used for clearer visualization in all figures of section 3.

3. Results and discussions

3.1. Time-proportional loading

The approximately linear stress–strain curves due to uniaxial tension in $xx$ , $yy$ , and $zz$ directions are shown in Figure 4(a). The significantly softer material behavior in the $yy$ direction indicates anisotropy, as expected due to the crystalline microstructure. Furthermore, the investigated silica exhibits strain rate independence and thus does not feature viscous components.

Figure 4.

Time-proportional deformation up to 4% strain with strain rates of 0.1–10% ns⁻¹: (a) tensile stress over tensile strain for uniaxial tension simulations in $xx$ , $yy$ , and $zz$ and (b) shear stress over shear strain for simple shear simulations in $xy$ , $xz$ , and $yz$ ; standard deviation of stress fluctuations indicated by error bars.

Very similar results in terms of viscosity and anisotropy are observed for the simple shear deformation in Figure 4(b). Initially, we observe a slightly stiffer material response in the $yz$ direction. However, at strains larger than 1.8 %, substantial instantaneous stress drops occur that are caused by the sudden transition (“snapping”) of the crystal lattice to an energetically more favorable state. Since our sample comprises only a single crystal, these stress drops are particularly pronounced. Since we can already exclude viscosity at this point, the material is either purely elastic or elastoplastic (cf. Figure 2).

By evaluating the slope in the stress–strain diagram Figure 4(a), the effective moduli $E_{x} = 198 GPa$ , $E_{y} = 142 GPa$ , and $E_{z} = 200 GPa$ are identified. These quantities do not yet represent a comprehensive constitutive description as usually determined from experiments. In comparison with the reported range of 46–92 GPa [55,56], it is evident that our silica behaves distinctly stiffer. This is probably due to the shift in dynamics induced by coarse graining [57] and emphasizes the importance of a thorough investigation of the mechanical properties of MD models. A similarly high Young’s modulus was reported by Yu et al. [58] for silica nanobeams studied with MD simulations.

3.2. Time-periodic loading

The stress response induced by sinusoidal uniaxial tension–compression loading with a strain amplitude of 4% is plotted in Figure 5. For deformation in the $xx$ and $yy$ direction, we observe stress–strain hystereses, induced under compression and disappearing with an increasing tensile load. Interestingly, there are equivalent paths for this inelastic behavior evident in Figure 5(a): While in the first and fourth load cycle, the stress already returns to the reversible path in the compression range; in the second and third cycle, this occurs only in the tension range. This is observed analogously for deformation in $yy$ direction (cf. Figure 5(b)). The area enclosed by the stress–strain curves in cyclic tests is equivalent to the dissipated energy. Since dissipation only occurs for inelastic behavior and viscosity has already been excluded, the discussed hysteresis indicates plastic effects. In contrast, the $zz$ direction displays no hysteresis and thus no plasticity, evident in cf. Figure 5(c). As shown by the insets in Figure 5(a)–(c), plasticity does not yet occur at a strain amplitude of 1.5 %.

Figure 5.

Time-periodic uniaxial tension: tensile stress over tensile strain for sinusoidal loading in $xx$ (a), $yy$ (b), and $zz$ (c) direction, with four tension–compression cycles with strain amplitude $ε^{a} = 4 %$ ; time-progression of loading indicated by colorbar, later cycles overlay former ones. Insets: stress–strain relation for $ε^{a} = 1.5 %$ .

Analogously, the stress–strain curves for sinusoidal simple shear are combined in Figure 6. In contrast to uniaxial deformation, hystereses occur only in the tensile region for deformation in $xz$ and $yz$ direction. The underlying plasticity is attributed to the abovementioned snapping of the crystal lattice (cf. section 3.1). No inelastic effects occur for deformation in $xy$ direction (cf. Figure 6(a)) or for strain amplitudes smaller than $ε^{a} = 1.5 %$ , as shown in the insets in Figure 6(a)–(c).

Figure 6.

Time-periodic simple shear: shear stress over shear strain for sinusoidal loading in $xy$ (a), $xz$ (b), and $yz$ (c) direction, with four tension–compression cycles with strain amplitude $ε^{a} = 4 %$ ; time-progression of loading indicated by colorbar, later cycles overlay former ones. Insets: stress–strain relation for $ε^{a} = 1.5 %$ .

In summary, the investigated silica can be characterized as elastoplastic with pronounced anisotropy. However, there is a purely elastic range up to 1.5% strain, which is limited by the lattice snapping under shear deformation.

3.3. Calibration of the constitutive model

In the nanocomposite, the nanoparticles are randomly oriented, so the overall directional dependence is expected to disappear. For a precise continuum mechanical model of the neat filler, the observed anisotropy nevertheless needs to be accounted for. When used as a filler for nanocomposites, silica is embedded into a softer polymer matrix, e.g., polystyrene, and therefore experiences comparatively small deformations. Therefore, we aim to accurately model the observed material behavior in the elastic range up to 1.5 % strain by means of a constitutive law based on the generalized Hooke’s law, cf. equation (6) in section “Elastic material models” in Appendix 1. To this end, optimal elastic constants have to be identified for the elasticity matrix $C$ (7). For this optimization problem, the complete deformation state obtained in terms of $ε (t)$ from the MD simulations with time-proportional loading serves as input while the corresponding stress state $σ (t)$ , i.e., all stress entries are used as goal function. The deformation constraint thus does not coincide with the MD boundary conditions, which seems counterintuitive, but is admissible due to the stress–strain duality of pure elasticity. Taking into account only the tensile range effectively implies tension–compression symmetry, which is, however, supported by our pseudo-experimental data in the range of $\pm 1.5 %$ , cf. insets of Figures 5 and 6.

A gradient-based, steepest-descent optimization [59] (lsqcurvefit) is used to identify a set of parameters capable of representing all six load cases simultaneously. To facilitate the optimization, four support points at $ε = [0, 0.5, 1.0, 1.5]^{t}$ are considered per load case, which is sufficient due to the approximately linear curves.

We introduce $Δ$ as the root-mean-square deviation (RMSD) of the stress averaged over all support points. As expected, large deviations of $Δ^{iso} = 328 MPa$ occur for an isotropic constitutive law. A significant reduction is achieved with transverse isotropy in $xy$ , $xz$ , and $yz$ direction where we obtain $Δ^{tran, xy} = 168 MPa$ , $Δ^{tran, xz} = 166 MPa$ , and $Δ^{tran, yz} = 179 MPa$ , respectively. However, the fact that all transverse orientations yield similar results contradicts transversal isotropy.

Considering an orthothropic or a fully anisotropic constitutive law further decreases the deviation to $Δ^{orth} = 156 MPa$ and $Δ^{aniso} = 99 MPa$ , respectively. We compare their performance in detail by means of the stress–strain diagrams for uniaxial tension and simple shear in Figures 7 and 8. Both approaches are capable of predicting the important stress in the loading direction accurately. For uniaxial tension, deviations are observed mainly for shear stresses, which partly deviate distinctly from the goal of 0 $MPa$ . For simple shear, the fit for both constitutive laws performs very well. The only exception is $σ_{xz}$ for deformation in $xy$ direction (cf. Figure 8(a)), which cannot be captured at all in the orthotropic case and only partially in the fully anisotropic model. Evidently, the fully anisotropic elasticity matrix is better suited for modeling the observed stress response, but requires 21 elastic constants, and thus significantly more than the orthotropic model with only nine parameters. The resulting optimal elastic coefficients are given in equations (8) and (9), respectively.

Figure 7.

Uniaxial tension fit: comparison of orthotropic and fully anisotropic stress predictions against pseudo-experimental goal data for deformation in $xx$ , $yy$ , and $zz$ direction: (a) uniaxial tension xx, (b) uniaxial tension yy, and (c) uniaxial tension zz.

Figure 8.

Simple shear fit: comparison of orthotropic and fully anisotropic stress predictions against pseudo-experimental goal data for deformation in (a) $xy$ , (b) $xz$ , and (c) $yz$ direction.

In crystalline materials, the unit cell defines the structure of the crystal lattice and thus its symmetries. These symmetries are also present in the overall elastic properties for a monocrystalline sample, and thus the degree of anisotropy depends directly on the unit cell structure. In our case, the unit cell is triclinic, i.e., fully asymmetric, which is consistent with the sample’s fully anisotropic elasticity.

4. Conclusion and outlook

In this paper, we have approximated the mechanical behavior of silica at the nanoscale by a continuum mechanical constitutive law. For this purpose, a CG sample was examined using MD deformation tests under PBC. Based on uniaxial and shear tests, time-proportional as well as time-periodic, we were able to characterize the present material as elastoplastic with pronounced anisotropy. For the range of small strains up to 1.5%, which is relevant in the context of polymer nanocomposites, we have identified the elastic constants of the generalized Hooke model. These can accurately predict the material behavior in this range.

This study demonstrates that a constitutive law on the nanoscale, as it is necessary for the modeling of an RVE in micromechanics, can be derived on the basis of MD simulations.

This contribution is complementary to an earlier characterization and modeling of polystyrene as an exemplary matrix material [33 –35]. Together with our extensive efforts to characterize the interphase between polymer and filler particles [25,26], we can now provide continuum mechanical descriptions of all constituents of a nanocomposite. This enables us to set up RVEs to apply the molecular-level findings on the microscale.

Footnotes

Appendix 1 Acknowledgements

We would like to thank the Erlangen Regional Computing Center (RRZE) at the Friedrich-Alexander-Universität Erlangen-Nürnberg for providing the high-performance computing resources necessary for this contribution. Furthermore, we thank Mattia Livraghi from the PULS Group at Friedrich-Alexander-Universität Erlangen-Nürnberg for the very fruitful discussions.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: Sebastian Pfaller is funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)—396414850 (Individual Research Grant “Identifikation von Interphaseneigenschaften in Nanokompositen”). Maximilian Ries, Christof Bauer, Felix Weber, Paul Steinmann, and Sebastian Pfaller are funded by the DFG—377472739 (Research Training Group GRK 2423 “Fracture across Scales—FRASCAL”).

ORCID iDs

Maximilian Ries

Felix Weber

Data accessibility statement

We provide the data published in this contribution via Zenodo. .

References

Jancar

Douglas

Starr

, et al. Current issues in research on structure–property relationships in polymer nanocomposites. Polymer 2010; 51(15): 3321–3343.

Maillard

Kumar

Fragneaud

, et al. Mechanical properties of thin glassy polymer films filled with spherical polymer-grafted nanoparticles. Nano Lett 2012; 12(8): 3909–3914.

Abdolmohammadi

Siyamak

Ibrahim

, et al. Enhancement of mechanical and thermal properties of polycaprolactone/chitosan blend by calcium carbonate nanoparticles. Int J Mol Sci 2012; 13(4): 4508–4522.

Han

Fina

. Thermal conductivity of carbon nanotubes and their polymer nanocomposites: a review. Prog Polym Sci 2011; 36(7): 914–944.

Suematsu

Arimura

Uchiyama

, et al. Transparent BaTiO₃/PMMA nanocomposite films for display technologies: facile surface modification approach for BaTiO₃ nanoparticles. ACS Appl Nano Mater 2018; 1(5): 2430–2437.

Kim

Yang

Green

. Tailoring the refractive indices of thin film polymer metallic nanoparticle nanocomposites. Langmuir 2012; 28(25): 9735–9741.

Holt

Bocharova

Cheng

, et al. Controlling interfacial dynamics: covalent bonding versus physical adsorption in polymer nanocomposites. ACS Nano 2016; 10(7): 6843–6852.

Holt

Griffin

Bocharova

, et al. Dynamics at the polymer/nanoparticle interface in poly (2-vinylpyridine)/silica nanocomposites. Macromolecules 2014; 47(5): 1837–1843.

Maghami

Shahrooz

Mehrabani-Zeinabad

, et al. Characterization of the polymer/particle interphase in composite materials by molecular probing. Polymer 2020; 205: 122792.

10.

Cui

You

. Mechanism of mechanical reinforcement for weakly attractive nanocomposites in glassy and rubbery states. Macromolecules 2021; 54(2): 824–834.

11.

Koo

. Introduction to nanotechnology. In: Koo

(ed.) Fundamentals, properties, and applications of polymer nanocomposites. Cambridge: Cambridge University Press, 2016, pp. 3–21.

12.

Holt

Bocharova

Cheng

, et al. Interplay between local dynamics and mechanical reinforcement in glassy polymer nanocomposites. Phys Rev Mater 2017; 1(6): 062601.

13.

Blivi

Benhui

Bai

, et al. Experimental evidence of size effect in nano-reinforced polymers: case of silica reinforced PMMA. Polym Test 2016; 56: 337–343.

14.

Zhao

Senses

, et al. Role of filler shape and connectivity on the viscoelastic behavior in polymer nanocomposites. Macromolecules 2015; 48(15): 5433–5438.

15.

Hashemi

Jouault

Williams

, et al. Enhanced glassy state mechanical properties of polymer nanocomposites via supramolecular interactions. Nano Lett 2015; 15(8): 5465–5471.

16.

Alonso-Redondo

Belliard

Rolle

, et al. Robustness of elastic properties in polymer nanocomposite films examined over the full volume fraction range. Sci Rep 2018; 8(1): 1–8.

17.

Leidner

Woodhams

. The strength of polymeric composites containing spherical fillers. J Appl Polym Sci 1974; 18(6): 1639–1654.

18.

Pukanszky

. Influence of interface interaction on the ultimate tensile properties of polymer composites. Composites 1990; 21(3): 255–262.

19.

Ashraf

Peng

Zare

, et al. Effects of size and aggregation/agglomeration of nanoparticles on the interfacial/interphase properties and tensile strength of polymer nanocomposites. Nanoscale Res Lett 2018; 13(1): 1–7.

20.

Chopra

Deshmukh

Peshwe

. Theoretical prediction of interfacial properties of PBT/CNT nanocomposites and its experimental evaluation. Mech Mater 2017; 109: 11–17.

21.

Zare

Rhee

. A modeling approach for Young’s modulus of interphase layers in polymer nanocomposites. Phys Mesomech 2020; 23(2): 176–181.

22.

Zare

. Modeling approach for tensile strength of interphase layers in polymer nanocomposites. J Colloid Interf Sci 2016; 471: 89–93.

23.

Peng

Rhim

Zare

, et al. Effect of “z” factor for strength of interphase layers on the tensile strength of polymer nanocomposites. Polym Composite 2019; 40(3): 1117–1122.

24.

Detrez

Yvonnet

, et al. Identification of elastic properties of interphase and interface in graphene-polymer nanocomposites by atomistic simulations. Compos Sci Technol 2021; 213: 108943.

25.

Pfaller

Possart

Steinmann

, et al. Investigation of interphase effects in silica-polystyrene nanocomposites based on a hybrid molecular-dynamics–finite-element simulation framework. Phys Rev E 2016; 93(5): 052505.

26.

Ries

Possart

Steinmann

, et al. A coupled MD-FE methodology to characterize mechanical interphases in polymeric nanocomposites. Int J Mech Sci 2021; 204: 106564.

27.

Jain

Ries

Pfaller

, et al. Addressing surface effects at the particle-continuum interface in a molecular dynamics and finite elements coupled multiscale simulation technique. J Chem Theory Comput 2022; 18: 2375–2387.

28.

Pfaller

Rahimi

Possart

, et al. An Arlequin-based method to couple molecular dynamics and finite element simulations of amorphous polymers and nanocomposites. Comput Methods Appl Mech Eng 2013; 260: 109–129.

29.

Pfaller

Kergaßner

Steinmann

. Optimisation of the capriccio method to couple particle-and continuum-based simulations of polymers. Multiscale Sci Eng 2019; 1(4): 318–333.

30.

Zhao

Steinmann

Pfaller

. A particle-continuum coupling method for multiscale simulations of viscoelastic-viscoplastic amorphous glassy polymers. Int J Numer Methods Eng 2021; 122: 7431–7451.

31.

Bocharova

Genix

A-C

Carrillo

J-MY

, et al. Addition of short polymer chains mechanically reinforces glassy poly (2-vinylpyridine)–silica nanoparticle nanocomposites. ACS Appl Nano Mater 2020; 3(4): 3427–3438.

32.

Hill

. Elastic properties of reinforced solids: some theoretical principles. J Mech Phys Solids 1963; 11(5): 357–372.

33.

Ries

Possart

Steinmann

, et al. Extensive CGMD simulations of atactic PS providing pseudo experimental data to calibrate nonlinear inelastic continuum mechanical constitutive laws. Polymers 2019; 11(11): 1824.

34.

Ries

Steinmann

Pfaller

. Characterization of polystyrene under shear deformation using molecular dynamics. In: Abali

Giorgio

(eds) Developments and novel approaches in nonlinear solid body mechanics. Cham: Springer, 2020, pp. 219–229.

35.

Zhao

Ries

Steinmann

, et al. A viscoelastic-viscoplastic constitutive model for glassy polymers informed by molecular dynamics simulations. Int J Solids Struct 2021; 226: 111071.

36.

Ghanbari

Ndoro

TVM

Leroy

, et al. Interphase structure in silica-polystyrene nanocomposites: a coarse-grained molecular dynamics study. Macromolecules 2012; 45(1): 572–584.

37.

Reith

Pütz

Müller-Plathe

. Deriving effective mesoscale potentials from atomistic simulations. J Comput Chem 2003; 24(13): 1624–1636.

38.

Ghanbari

Rahimi

Dehghany

. Influence of surface grafted polymers on the polymer dynamics in a silica-polystyrene nanocomposite: a coarse-grained molecular dynamics investigation. J Phys Chem C 2013; 117(47): 25069–25076.

39.

Liu

Pfaller

Rahimi

, et al. Uniaxial deformation of polystyrene-silica nanocomposites studied by hybrid molecular dynamics-finite element simulations. Comp Mater Sci 2017; 129: 1–12.

40.

Munaò

Pizzirusso

Kalogirou

, et al. Molecular structure and multi-body potential of mean force in silica-polystyrene nanocomposites. Nanoscale 2018; 10(46): 21656–21670.

41.

Bauer

Ries

Pfaller

. Accelerating molecular dynamics simulations by a hybrid molecular dynamics-continuum mechanical approach. Soft Materials 2022; 1–16.

42.

Larrañaga

Lewis

. Hawley’s condensed chemical dictionary. 16th ed. New York: John Wiley & Sons, Inc., 2016.

43.

Steinmann

Runesson

. The catalogue of computational material models: basic geometrically linear models in 1D. Cham: Springer Nature, 2021.

44.

Haupt

. Continuum mechanics and theory of materials. Cham: Springer Science & Business Media, 2013.

45.

Pfaller

. Discrete and Continuous Methods for Modelling and Simulation of Polymer Materials. Friedrich-Alexander-Universität Erlangen-Nürnberg, Lehrstuhl für Technische Mechanik (LTM), 2021, 216.

46.

Treloar

LRG

. The physics of rubber elasticity. 3rd ed. Oxford: Clarendon Press, 1975.

47.

Evans

Holian

. The Nose–Hoover thermostat. J Chem Phys 1985; 83(8): 4069–4074.

48.

Martyna

Tobias

Klein

. Constant pressure molecular dynamics algorithms. J Chem Phys 1994; 101(5): 4177–4189.

49.

Ries

Weber

Striegel

, et al. Multiscale FE-MD coupling: influence of the chain length on the mechanical behavior of coarse-grained polystyrene. In: Proceedings of World Congress on Computational Mechanics, Paris, 11–15 January 2021.

50.

Thompson

Plimpton

Mattson

. General formulation of pressure and stress tensor for arbitrary many-body interaction potentials under periodic boundary conditions. J Chem Phys 2009; 131(15): 154107.

51.

Subramaniyan

Sun

. Continuum interpretation of virial stress in molecular simulations. Int J Solids Struct 2008; 45(14–15): 4340–4346.

52.

Zhou

. A new look at the atomic level virial stress: on continuum-molecular system equivalence. P Roy Soc Lond A Mat 2003; 459(2037): 2347–2392.

53.

Zimmerman

Webb

III Hoyt

, et al. Calculation of stress in atomistic simulation. Model Simul Mater Sc 2004; 12(4): S319.

54.

Savitzky

Golay

. Smoothing and differentiation of data by simplified least squares procedures. Anal Chem 1964; 36(8): 1627–1639.

55.

Carlotti

Doucet

Dupeux

. Elastic properties of silicon dioxide films deposited by chemical vapour deposition from tetraethylorthosilicate. Thin Solid Films 1997; 296(1–2): 102–105.

56.

Petersen

. Dynamic micromechanics on silicon: techniques and devices. IEEE T Electron Dev 1978; 25(10): 1241–1250.

57.

Karimi-Varzaneh

van der Vegt

Müller-Plathe

, et al. How good are coarse-grained polymer models? A comparison for atactic polystyrene. ChemPhysChem 2012; 13(15): 3428–3439.

58.

Sun

Zhang

, et al. Study on size-dependent Young’s modulus of a silicon nanobeam by molecular dynamics simulation. J Nanomater 2013; 2013: 319302.

59.

Coleman

. An interior trust region approach for nonlinear minimization subject to bounds. SIAM J Optimiz 1996; 6(2): 418–445.

60.

Ries

Bauer

Weber

, et al. Calibrating an anisotropic constitutive model based on molecular dynamics pseudo experiment. Geneva: Zenodo, 2021.

61.

Plimpton

. Fast parallel algorithms for short-range molecular dynamics. J Comput Phys 1995; 117(1): 1–19.

62.

Humphrey

Dalke

Schulten

. VMD—visual molecular dynamics. J Mol Graphics 1996; 14: 33–38.

63.

Kohlmeyer

Vermaas

Braun

. Topotools, 2020. https://doi.org/10.5281/zenodo.3845031.

64.

Tadmor

Miller

Elliott

. Continuum mechanics and thermodynamics: from fundamental concepts to governing equations. Cambridge: Cambridge University Press, 2012.

Characterization of the material behavior and identification of effective elastic moduli based on molecular dynamics simulations of coarse-grained silica

Abstract

Keywords

1. Introduction

1.1. Empirical models

1.2. Molecular dynamics studies

1.3. Micromechanics

2. Model and simulation details

2.1. Coarse-grained silica

2.2. Material classification scheme and characterization strategy

2.3. Applied deformation

2.3.1. Uniaxial tension

2.3.2. Simple shear

2.3.3. Simulation details

2.4. Stress measures

3. Results and discussions

3.1. Time-proportional loading

3.2. Time-periodic loading

3.3. Calibration of the constitutive model

4. Conclusion and outlook

Footnotes

Appendix 1

Acknowledgements

Funding

ORCID iDs

Data accessibility statement

References