Introduction

Since the early days of polymeric materials, filled systems remain essential for improving and/or inducing new functional properties otherwise unavailable in plastics. For a broad range of applications, fillers can improve the mechanical performance of the soft matrices, enhance the semi-conductive and heat transfer properties of otherwise insulating polymers, improve their gas barrier properties, etc. The development of novel fillers for polymers remains a perpetual challenge with e.g. new generations of 2D materials currently being developed with tailored functionalities (Gupta et al. 2015). These fillers could completely revolutionize the way to design polymer nanocomposites due to their extraordinary potential (Kim and Macosko 2008; Kim et al. 2010; Gkourmpis 2014; Paul and Robeson 2008; Potts et al. 2011). The main reason for this potential comes from the rich and unique heterostructures currently being investigated that could allow for the design of totally new classes of materials with improved multifunctional properties. Of particular importance are conductive anisotropic fillers (high aspect ratio), such as metal nanowires, graphene, and carbon nanotubes. Their steadily commercial introduction over the last 10–15 years is leading a revolution in polymer nanocomposites. Depending on the type of filler and the targeted application, the filler concentrations required can be fairly significant, a typical example being carbon black composites. This in turn can lead to an increase in the overall melt viscosity, which can lead to a number of processing difficulties (Shenoy 1999). Consequently, the potential of nanofillers incorporated in polymers to achieve high conductivity and superior mechanical properties due to their unique geometry at low or very low concentrations has attracted enormous scientific and commercial attention. Therefore, a detailed understanding of the structure-property relations of these nanocomposites is required. In this framework, a critical aspect for the performance of polymer nanocomposites, especially for their potential commercialization, is their melt flow properties. Particularly, from rheological point of view, the flow-field-filler interactions in the presence of small and large deformations/deformation rates are crucial for understanding the tailoring of filler concentration, dispersion, and morphology to achieve enhanced nanocomposite properties.

Polymer nanocomposites and percolation

The classification of nanofillers in this publication refers is based on their dimensionality as categorized by Ashby et al. (2009). Thus, 0-D nanofillers are defined as having all dimensions at nanoscale (nanoparticles, e.g., PCC, precipitated calcium carbonate; nG, gold nanoparticles), 1-D having one of the characteristic dimensions not at nanoscale (nanorods, nanowires, e.g., S/MWCNT, single and multi-walled carbon nanotubes; (MnO2)NW, manganese dioxide nanowires), 2-D having only one dimension at nanoscale (nanocoatings and nanofilms, e.g., G, graphene and its derivatives; incl. rGO, reduced graphene oxide; GnP, graphite nanoplatelets; nC, nanoclay), and 3-D with the overall filler dimensions exceed nanoscale, however, having a filler structure that involves the presence of features at nanoscale that are determinant for its reinforcement performance. It should be noted that for regulatory purposes, other dimensional categorization schemes can be used ISO 26824:2013, 1.1 (2013).

The introduction of fillers into a polymer matrix eventually leads to the creation of a network structure. Its connectivity and spatial arrangement and the resulting macroscopic effects can be described by the percolation theory. In all its variations, percolation theory focuses on critical phenomena that originate from the spatial arrangement of the given network and results in sharp transitions in the overall behavior of the system of interest Kirkpatrick (1973). Most commonly, percolation can be understood in the context of electrical conductivity (Gkourmpis 2016) but this is not restrictive. In such systems as the insulating polymer matrix is loaded with a conductive filler, the eventually formed network leads to a sharp insulation-conductor transition (Lux 1993) in accordance with percolation theory

$$ {\sigma}_{el}\approx {\sigma}_0{\left(\phi -{\phi}_c\right)}^t $$
(1)

with σ0 a pre-exponential factor that is dependent on the conductivity of the filler, the network topology, and the types of contact resistance. The terms ϕ and ϕc correspond to the filler concentration and the critical concentration at the transition (also known as percolation threshold) (Foygel et al. 2005). The critical exponent t is also conductivity-dependent with universal values of t ≈ 1.33 and t ≈ 2 for two and three dimensions respectively (Stauffer and Aharony 1987; Sahimi 1994). Despite the universality of the critical exponent, a wide range of values as high as t ≈ 10 have been reported (Gkourmpis 2014; Bauhofer and Kovacs 2009; Gkourmpis et al. 2019), but this divergence can mainly be associated with complex tunneling transport phenomena (Balberg 1987; Balberg 2009; Grimaldi and Balberg 2006; Johner et al. 2008) without disregarding more practical limitations such as filler variability in terms of production, entanglement, interconnections, and surface chemistry (Mutiso and Winey 2012).

The percolated network properties are expected to depend on the type of filler, e.g., some fillers are flexible (CNT) whereas certain morphological types favor inter-connectivity compared to others, e.g., 2D vs. 0D nanofillers (Hassanabadi and Rodrigue 2013). The mechanical response of the filler network depends also on the concentration, with a weakly connected network expected at and in the upper vicinity of the percolation threshold. With increasing filler concentration, the network connectivity increases thus resulting in a more structurally persistent network.

Thus far in our discussion, the role of the polymer matrix in the way the filler organizes has not been considered. Since percolation theory is in most of its variations an essentially statistical theory, spatial arrangement and dispersion are driven by the filler’s size and dimensionality. Obviously, such an approach might be successful on dilute systems but in the case of polymers, the existence of rich morphologies, extensive chain connectivity, and significantly higher molecular weights and viscosities provide significant limitations. The introduction of even simple fillers like carbon black in a polymer matrix has the potential to lead to significant agglomeration due to the inability of the filler to disperse in a totally random manner. The main mechanism behind the extensive agglomeration can be traced in the existence of weak attractive interactions between the fillers (Grimaldi and Balberg 2006), and the effect can be significantly augmented when using fillers that can reach dimensions up to several micrometers (Alig et al. 2012). In this picture, filler-filler, filler-polymer, and polymer-polymer interactions play a vital role, and polymer mobility is of paramount importance for the facilitation of the local and long range arrangements of the filler in the matrix (Pötschke et al. 2004). With the term polymer mobility, we refer to all the length scales associated with a polymer, ranging from the local motion and spatial segmental arrangements as seen by the chain conformation to the Rouse regime and chain reptation. Rheologically, network and dispersion properties have been typically assigned to the slope of the storage modulus G in the terminal regime at percolation for very low angular frequencies. Generally, it is considered that the stress response at small strains is the result of percolated particle network, whereas the large strain behavior is determined by the polymer chain dynamics (Senses and Akcora 2013).

Nonlinear oscillatory shear analysis

Nonlinear oscillatory shear analysis via Fourier-transform rheology (FT-Rheology) and large amplitude oscillatory shear (LAOS) is an emerging framework for the analysis of nonlinear effects in complex fluids. The main arguments for the analysis are the increased sensitivities compared to linear viscoelastic measurements (van Dusschoten and Wilhelm 2001; Hyun et al. 2011) and their potential to reveal material response features not observable in linear viscoelastic measurements (Ewoldt et al. 2008; Hyun et al. 2011). The basic nonlinear parameters and definitions necessary for the interpretation of experimental data are outlined here. Detailed descriptions of the nonlinear analyses employed can be found elsewhere, e.g., see Wilhelm (2002), Ewoldt et al. (2008), and Hyun et al. (2011).

In linear viscoelastic (strain-controlled) oscillatory shear tests, a sinusoidal strain input, γ(t) = γ0 sin (ωt), results in a sinusoidal shear stress output, σ12 = σ0 sin (ωt + δ), shifted with the phase angle δ, where γ0 and ω are the imposed strain amplitude and angular frequency. The shear rate can therefore be expressed as \( \dot{\gamma}(t)={\gamma}_0\omega \cos\ \left(\omega t\right) \), with \( {\dot{\gamma}}_0={\gamma}_0\omega \) as the shear rate amplitude. The Fourier spectrum of the shear stress, Σ12, therefore contains exclusively the fundamental intensity I1 of the imposed excitation, ω/ωi = 1, Fig. 1a. A nonlinear viscoelastic shear stress response will in contrast be non-sinusoidal/distorted and the corresponding Fourier spectrum will contain higher harmonics, In, in addition to I1, Fig. 1b. By comparing the shear stress, σ12, expressed from a Taylor series expansion of the shear viscosity, \( \eta \left(\dot{\gamma}\right) \), in dynamic oscillatory shear flow, to a corresponding Fourier series, combined with symmetry arguments, the time-dependent shear stress can be expressed as (Wilhelm et al. 1998; Hyun et al. 2011; Naue et al. 2018):

(2)

where I1, I3, I5... are the intensities of the fundamental applied angular frequency, ω0, and its higher (odd) harmonics. Thus, the higher harmonics can be used to quantify the nonlinear material response. The most basic parameter is the third relative higher harmonic, I3/I1 = I3/1, as it contains the dominant nonlinear contribution. The variation of I3/1 during a strain sweep test for a typical polymeric material exemplified by the matrix polymer (polypropylene) used in this study is presented in Fig. 1c. At small strain amplitudes (SAOS, small amplitude oscillatory shear), the measured signal corresponds to instrumentation noise and indicates the sensitivity limits of the torque sensor. At a critical shear strain amplitude, the nonlinearities become detectable and increase with \( {I}_{3/1}\propto {\gamma}_0^2 \), region called medium amplitude oscillatory shear (MAOS) or intrinsic large amplitude oscillatory shear ([LAOS]). Based on the quadratic I3/1 scaling in MAOS, the Q-parameter can be defined as nonlinear material parameter from FT-Rheology (Hyun and Wilhelm 2009):

$$ Q\equiv \frac{I_{3/1}}{\gamma_0^2} $$
(3)
Fig. 1
figure 1

Examples of a linear and b nonlinear shear stress response and corresponding Fourier spectra, and c dynamic strain sweep test comparing the dynamic moduli, G′ and G′′ and the third relative higher harmonic, I3/1

The large amplitude oscillatory shear regime (LAOS) is reached when the quadratic scaling with the strain amplitude is lost.

To attempt to gain a physical insight into material nonlinear behavior, the total stress can be decomposed into elastic and viscous components using Chebyshev polynomials as (Ewoldt et al. 2008):

$$ {\displaystyle \begin{array}{c}{\sigma}_{12}\left(t,\omega, {\gamma}_0\right)={\sigma}_e\left(\gamma (t)\right)+{\sigma}_v\left(\dot{\gamma}(t)\right)=\\ {}={\gamma}_0\sum {e}_i{T}_i\left(\gamma /{\gamma}_0\right)+{\gamma}_0\omega \sum {v}_i{T}_i\left(\dot{\gamma}/\dot{\gamma_0}\right)\end{array}} $$
(4)

where Tn is the nth-order Chebyshev polynomial of the first kind, and ei and vi are the elastic and respectively viscous Chebyshev coefficients. A visual inspection of the nonlinear response can be readily performed using elastic and viscous Lissajous-Bowditch diagrams, Fig. 2. In the viscous and elastic Lissajous-Bowditch (LB) representations, the outer closed trajectories represent σ(t)/σmax vs. γ(t)/γ0 and σ(t)/σmax vs. \( \dot{\gamma}(t)/{\dot{\gamma}}_0 \), respectively. The red/blue lines represent the total elastic/viscous stresses, i.e., σe/σmax vs. γ(t)/γ0 and σv(t)/σmax vs. \( \dot{\gamma}(t)/{\dot{\gamma}}_0 \), respectively. A linear viscoelastic response is thus represented by elliptic diagrams (the black lines in Fig. 2). Any deviations therefrom signify a nonlinear viscoelastic material response (light gray lines in Fig. 2). Based on the LB diagrams, two elastic and viscous material nonlinear parameters can be thus defined (Ewoldt et al. 2008) as

$$ {G}_M^{\prime}\equiv {\left.\frac{d\sigma}{d\gamma}\right|}_{\gamma =0}={e}_1-3{e}_3+\dots $$
(5)
$$ {G}_L^{\prime}\equiv {\left.\frac{\sigma }{\gamma}\right|}_{\gamma =\pm {\gamma}_0}={e}_1+{e}_3+\dots $$
(6)
$$ {\eta}_M^{\prime}\equiv {\left.\frac{d\sigma}{d\dot{\gamma}}\right|}_{\dot{\gamma}=0}={v}_1-3{v}_3+\dots $$
(7)
$$ {\eta}_L^{\prime}\equiv {\left.\frac{\sigma }{\dot{\gamma}}\right|}_{\dot{\gamma}=\pm \gamma \dot 0}={v}_1+{v}_3+\dots $$
(8)

representing the tangent slopes at zero strain/rate (•M) and secants at maximum strain/rate (•L), see Fig. 2. Thus, in the limit of linear viscoelasticity, the nonlinear elastic and viscous moduli are identical to the linear viscoelastic correspondents, i.e., \( {G}_M^{\prime }={G}_L^{\prime }={G}^{\prime}\left(\omega \right) \) and \( {\eta}_M^{\prime }={\eta}_L^{\prime }={\eta}^{\prime}\left(\omega \right) \), respectively. Using the relationship between the defined nonlinear viscous and elastic moduli in (5)–(Dealy and Wissbrun 1999), the distortions in the LB diagrams can be described using the sign of elastic and viscous Chebyshev coefficients as (Ewoldt and Bharadwaj 2013):

Fig. 2
figure 2

Examples of linear (black line) and nonlinear (gray line) a elastic and b viscous Lissajous-Bowditch diagrams and the defined associated moduli

e1 > 0 average elastic stiffening driven by

e3 > 0 large instantaneous strains

e3 < 0 large instantaneous rates-of-strain

e1 < 0 average elastic softening driven by

e3 > 0 large instantaneous rates-of-strain

e3 < 0 large instantaneous strains

v1 > 0 average viscous thickening driven by

v3 > 0 large instantaneous rates-of-strain

v3 < 0 large instantaneous strains

v1 < 0 average viscous thinning driven by

v3 > 0 large instantaneous strains

v3 < 0 large instantaneous rates-of-strain

A similar interpretation can be assigned to the strain-stiffening (S) and shear thickening (T) dimensionless ratios (Ewoldt et al. 2008):

$$ S\equiv \frac{G_L^{\prime }-{G}_M^{\prime }}{G_L^{\prime }} $$
(9)
$$ T\equiv \frac{\eta_L^{\prime }-{\eta}_M^{\prime }}{\eta_L^{\prime }} $$
(10)

In this case S > 0 indicates intracycle strain-stiffening, S < 0 intracycle strain-softening, T > 0 intracycle shear thickening, and T < 0 intracycle shear thinning.

Nonlinear rheology of polymer nanocomposites

Known nonlinear viscoelastic properties of nanofilled polymeric systems in the nonlinear oscillatory shear analysis framework outlined in the “Nonlinear oscillatory shear analysis” section are summarized in Table 1. The scientific studies are sorted by nanoparticle dimensionality and generic molecular topology. The literature survey is limited to thermoplastic nanocomposites tested in standard rotational rheometry. However, it should be noted that FT-Rheology has been applied to polymer composites and rubber-based composites and nanocomposites, e.g., see Leblanc (2006a, b) and Schwab et al. (2016). The range of polymer-nanofiller combinations cited in Table 1 reflects systems that generally provide good dispersion properties. Thus, there is a prevalence for 1D and 2D nanofillers combined mainly with linear/short chain branched (L/SCB) polymers. Particulate nanofillers, 0D, due to their morphological limitations in forming networks, are considered for comparison (Lim et al. 2013).

Table 1 Breviary of experimental findings from the scientific literature on polymer-based nanocomposites from nonlinear oscillatory shear experiments categorized by filler dimensionality and matrix molecular topology

Most studies focused on FT-Rheology and it is generally agreed that percolation thresholds can be detected with increased sensitivity in FT-Rheology compared to linear viscoelastic oscillatory tests. This has been shown for e.g. (L/SCB)/1D (Lim et al. 2013; Ahirwal et al. 2014; Hassanabadi et al. 2014) and (L/SCB)/2D (Hyun et al. 2012; Hassanabadi et al. 2014; Kádár et al. 2017; Gaska et al. 2017). The assessment is generally based on an increase in the magnitude of I3/1, often quantified via the Q-parameter. Exceptions are 0D nanofiller-based systems, e.g., (L/SCB)/0D (Lim et al. 2013), for which the percolation threshold was not reached due to the particle shape (Knauert et al. 2007; Hassanabadi and Rodrigue 2013). A relative decrease in Q0 with increasing filler content before the percolation threshold has also been reported (Ahirwal et al. 2014). The material nonlinear response as quantified by I3/1 has been shown to be sensitive to the filler dimensionality (Lim et al. 2013; Hassanabadi et al. 2014) and to the degree of orientation in anisotropic nanofillers as induced through AC/DC currents (Hyun et al. 2012). In addition, particularly in L/SCB polymers, a maximum in I3/1, dI3/1/0 = 0, was apparent in the measuring range for 1D (Lim et al. 2013) and 2D nanofillers (Lim et al. 2013; Du et al. 2018). Particular changes in I3/1 scaling have also been reported for 2D (Kádár et al. 2017; Gaska et al. 2017; Gaska and Kádár 2019) and 2D hybrid nanocomposites (Kádár et al. 2017; Gaska and Kádár 2019), i.e., \( {I}_{3/1}\propto {\gamma}_0^n \) with n ≠ 2 in long-chain branched (LCB) polymers. Overall, most nonlinear effects in the percolation region were attributed to dispersion properties (Lim et al. 2013; Hyun et al. 2012; Kamkar et al. 2018; Du et al. 2018). By comparing various filled systems, including nanofiller dimensionality, a new parameter was proposed (Lim et al. 2013) showing good correlation to how well a filler is dispersed in a polymer.

Chebyshev stress decomposition analysis has been applied to (L/SCB)/1D (Kamkar et al. 2018) and (L/SCB)/2D (Du et al. 2018) nanocomposites. Strain-stiffening, S, Eq. (9), and shear thickening, (T), Eq. (10), ratios have been thus shown to depend on composition and reinforcement type. S was positive in the nonlinear region, thus showing increasing strain-stiffening behavior with increasing filler concentration. T was shown to be negative and increased with the addition of fillers (Kamkar et al. 2018) but also strain amplitude dependent (strain sweep) (Du et al. 2018). Over a critical ϕ within the range of percolation thresholds estimated, at low strain amplitude T > 0 and at the high end of γ0 range T < 0. Thus, the behavior changes from shear thickening to shear thinning, with the shear thinning character increasing with ϕ in the nonlinear regime.

Present study

In this study, we present an industrially relevant melt mixed system where a novel 3D hierarchical reduced graphene oxide (HrGO) filler has been successfully dispersed in isotactic polypropylene. Consequently, the electrical percolation threshold was observed in the vicinity of ϕc ≈ 0.94 wt% and ranks lowest among similar systems in similar preparation conditions (Gkourmpis et al. 2019). Such percolation thresholds are comparable to graphene-based systems prepared using solution mixing techniques (Kim and Macosko 2008). The work concentrates on establishing a closer understanding of the particular specific interactions between the filler network and the polymer matrix in dynamic conditions, in relation to the filler network formation. This is pursued through the use of nonlinear rheological behavior in oscillatory shear. The nonlinear analysis is based on Fourier-transform rheology and Chebyshev polynomial stress decomposition. We show that the electrical percolation threshold corresponds to a unique nonlinear signature whose characteristics develop with increasing concentration which could potentially provide a clear distinction between the different network evolutionary mechanisms and the influence of shear thereon.

Materials and experiments

Materials

The polymer matrix used in this study was a commercial isotactic polypropylene (PP) supplied by Borealis AB (Stenungsund, Sweden). The PP used was highly isotactic (> 90%) having a molecular weight and dispersity index of Mw = 300 kg/mol and Mw/Mn = 8, respectively. As nanofiller, the present study features a novel type of 3D hierarchical structure that consists of stacks of 3–5 monolayers (as estimated by the manufacturer) of reduced graphene oxide (HrGO). Thus, while overall filler dimensions exceed nanoscale, the filler structure involves the presence of features at nanoscale level, i.e., the rGO sheets, and can be thus classified as a 3D nanofiller (Ashby et al. 2009). We note that the hierarchical structure is a direct result of the proprietary exfoliation process (Gkourmpis et al. 2019). The HrGO was produced via a modified Staudenmaier method (Staudenmaier 1898) and was provided by Cabot Corporation (Boston, MA, USA). A SEM image of the filler morphology in the powder state is shown in Fig. 3a.

Fig. 3
figure 3

a SEM micrograph of the hierarchical reduced graphene oxide (HrGO) in powder form and the associated structure. b Electrical conductivity as function of filler loading and the evaluation of the percolation threshold using Eq. (1). c SEM micrographs for PP. d SEM micrographs for PP-0.5% HrGO

Sample preparation and testing

The nanocomposites presented in this study were prepared by means of Brabender batch melt mixer Type W50 driven by a Brabender Plasticorder (Brabender GmbH, Duisburg, Germany). The compounding was performed at 210 °C in two steps (Gkourmpis et al. 2019). In the first step, the polymer was introduced in the mixing chamber and allowed to melt at 20 rpm for 15 min. Thereafter, the filler was added and the composite was then mixed for homogenization at 50 rpm for another 15 min. More details concerning the nanofiller structure, preparation, etching procedures, and the nanocomposite electrical, thermal, and mechanical performance can be found elsewhere (Gkourmpis et al. 2019). Electrical conductivity measurements evaluating the percolation threshold, ϕc ≈ 0.94 with a critical exponent t = 7.75, are shown in Fig. 3b (Gkourmpis et al. 2019).

Representative SEM micrographs for the unfilled PP and the PP/HrGO nanocomposites are shown in Fig. 3b, c. The microstructural observations of shown composites were carried out on cryo-fractured and chemically etched samples using a Leo Ultra 55 (Carl Zeiss AB, Oberkochen, Germany) Scanning Electron Microscope.

Linear and nonlinear oscillatory shear tests were performed on an Anton Paar MCR702 TwinDrive rheometer (Graz, Austria) in strain-controlled mode (separate motor-transducer). The frequency sweeps were performed in counter-oscillation. All tests were performed at 200 °C with a 25-mm plate-plate geometry and a gap of 1 mm. The nonlinear analysis of the shear stress output signal was performed in the framework of Fourier-transform analysis and Chebyshev polynomial decomposition, as summarized in the “Nonlinear oscillatory shear analysis” section.

Results and discussion

Linear viscoelastic frequency sweep results are presented for reference in Fig. 4. The complex viscosity functions, |η|, exhibit the expected dependence on concentration, with a significant increase for ϕ > 1 wt%, Fig. 4a. Within the applied angular frequency range, for ϕ ≤ 2.5 wt%, the viscosity functions can be readily fitted using the Carreau-Yasuda model, \( \left(\eta -{\eta}_{\infty}\right)\left({\eta}_0-{\eta}_{\infty}\right)={\left[1+{\left(\lambda \dot{\gamma}\right)}^a\right]}^{\left(n-1\right)/(a)} \), where η0, ∞ are the zero- and infinite-shear viscosities, n is the flow index, λ a characteristic relaxation time, and α the material exponent. A pronounced increase in viscosity at low shear rates is apparent at concentrations ϕ > 2.5 wt%, with the power law model, \( \eta =K{\dot{\gamma}}^{n-1} \), where K and n are the consistency and flow behavior indices, appropriately describing the shear thinning behavior for \( \dot{\gamma}>10 \) 1/s. This is considered to be indicative of a yield stress behavior, as strong attractive forces between nanofillers are known to lead to the formation of structures that could support a yield stress (Dealy and F. 1999). A similar qualitative concentration dependence applies to the dynamic moduli, Fig. 4b. Notably, a storage modulus, G, plateau in the limit of low angular frequencies, signaling the presence of an additional elastic contribution due to filler network formation is apparent only for concentrations ϕ > 2.5 wt%, in the range of frequencies investigated. This would suggest a percolation threshold significantly higher than the threshold determined from electrical measurements via percolation theory, ϕc ≈ 0.94 wt% (Gkourmpis et al. 2019). We note that evidence of percolation in the terminal region could be evidenced at angular frequencies below the lower limit investigated. However, measurements at very low angular frequencies are time consuming and an apparent plateau could also be reached due to e.g. oxidation effects.

Fig. 4
figure 4

Linear viscoelastic frequency sweeps of 3D hierarchical reduced graphene oxide-polypropylene (PP/HrGO) nanocomposites. a Magnitude of the complex viscosity, |η| and b storage modulus, G ′

The nonlinear material response as expressed by the third relative higher harmonic, I3/1, from dynamic strain sweep tests is presented in Fig. 5 for representative filler concentrations. In the absence of fillers (ϕc = 0), Fig. 5a, the linear-nonlinear transitions of the PP is as expected for a commercial polydisperse polymer (Hyun et al. 2011), with the SAOS region characterized by instrumentation noise, MAOS by \( {I}_{3/1}\propto {\gamma}_0^2 \), and the LAOS region beyond. The nonlinear material response for concentrations up to ϕ = 0.8 wt% is qualitatively similar to the unfilled polymer, see ϕ < 1 wt% in Fig. 5b. In contrast, for ϕ = 1 wt%, i.e., electrical percolation threshold, ϕ ≈ ϕ1, Gkourmpis et al. (2019), a significant change in I3/1 scaling is observed, Fig. 5c. Thus, the I3/1 scaling to the strain amplitude, γ0, becomes a function of the imposed angular frequency, ω, i.e., \( {I}_{3/1}\propto {\gamma}_0^n \), n = n(ω). Thus, the following scaling exponents can be approximated: \( {I}_{3/1}\propto {\gamma}_0^{0.7} \) for ω = 0.6 rad/s, \( {I}_{3/1}\propto {\gamma}_0^{1.1} \) for ω = 1 rad/s. For ω > 1 rad/s, it can be approximated that \( {I}_{3/1}\propto {\gamma}_0^2 \) similarly to concentrations below percolation, ϕ < ϕc. Above the electrical percolation threshold, ϕ = 1.5%, Fig. 5d, a MAOS scaling similar to the unfilled polymer is observed, independent of the applied ω. However, nonlinear material response features that are not present for ϕ < ϕc could be identified. Namely, at the SAOS-MAOS transition, a new region can be identified (γ0 ∈ (2, 10)%), where I3/1 ∝ γn, n ≈ 0 for ω < 2 rad/s. In addition, differences in the LAOS region could also be inferred.

Fig. 5
figure 5

Third relative higher harmonic, I3/1, from dynamic strain sweeps for the a unfilled polypropylene, PP (ϕc = 0), b before the electrical percolation threshold, ϕ < ϕc, c at the electrical percolation threshold (ϕ ≈ ϕc), and d directly above the percolation threshold. The shaded data points are estimated to be within the instrumentation noise region

The dependence on the applied angular frequency, ω, and at percolation, ϕ ≈ ϕc, could signify that the weakly formed network at percolation can be disrupted by the resulting shear rate. For high ω, I3/1 scaling is similar to the previous concentrations before percolation. The conjecture is further reinforced by the n ≈ 0 scaling observed for supercritical concentration, ϕ = 1.5 wt%, at low ω. The extrapolation to a n = 2 exponents, e.g., see Fig. 6, is based on solid mathematical foundations, i.e., see the ratio between I3 and I1 in Eq. 2. However, in this publication, we suggest that in certain conditions at and above the percolation threshold, the shear stress output, σ12, contains a contribution of the filler network, not observable in the linear viscoelastic data but rather amplified by the use of FT-Rheology. Thus, the scaling exponents of I3/1 should be perceived through the viewpoint of filler network properties.

Fig. 6
figure 6

aQ-parameter, Eq. (3), dependence on the strain amplitude and fittings of \( Q/{Q}_0={\left[1+{\left({C}_{1\gamma 0}\right)}^{C_2}\right]}^{\left({C}_3-1\right)/{c}_2} \). b Normalized zero-strain Q, Q0 = limγ → 0(Q), on the filler loading, ϕ, with the percolation threshold from the percolation theory determined from electrical conductivity, ϕc, also represented, Gkourmpis et al. (2019)

Changes in nonlinear I3/1 magnitude can be readily represented using the Q-parameter, Eq. 3, Hyun and Wilhelm (2009). This is exclusively based on the quadratic scaling of I3/1 with γ0 as deduced from Eq. 2 and thus disregards any possible anomalies at the linear-nonlinear transition. The Q-parameter variation with strain amplitude at ω = 1 rad/s for all filler concentrations is shown in Fig. 6a. The data points presented exclude the “instrumentation noise” as well as any significant non-quadratic scaling at the linear-nonlinear transition. Similarly to the procedure of Lim et al. (2013), the Q-parameter data was fitted using an equation similar to the Carreau-Yasuda model, i.e., \( Q/{Q}_0={\left[1+{\left({C}_1{\gamma}_0\right)}^{C_2}\right]}^{\left({C}_3-1\right)/{C}_2} \), with C1, 2, 3 being the corresponding fitting parameters. Therefrom, the zero-strain Q-parameter can be extrapolated, \( {Q}_0= li{m}_{\gamma_0\to 0}\left[Q\left({\gamma}_0\right)\right] \). The relative zero-strain Q-parameter, Q0(ϕ)/Q0(ϕ = 0), is shown in Fig. 6b while the errors in estimating the fit parameters are presented in Table S1 (Supplementary information). The electrical percolation threshold, ϕc, is overlaid on the graph. It can thus be stated that in quantitative terms, Q0(ϕ)/Q0(0) experiences a significant shift at ϕc = 0.8 wt%, below ϕc.

Based on Q-parameter data interpretation, it can be inferred that when considering the magnitude of the nonlinear material response, the relative increase in Q0(ϕ)/Q0(0) below the electrical percolation threshold could signify the detection of a rheological percolation threshold, meaning sufficient amount of filler is present to significantly alter molecular chain dynamics. This is due to the more gradual nature expected for the rheological threshold (Münstedt and Starý 2016). It should be noted that the high t-exponents, t ≈ 7.75, resulted from electrical conductivity measurements for the PP/HrGO nanocomposites (Gkourmpis et al. 2019) suggest that the filler network conducts predominantly via tunneling mechanisms (Balberg 1998; Gkourmpis 2016). Thus, it can be argued that the nanocomposite experiences an extended percolation threshold over a wider region, where at a rheological percolation, the chain mobility is significantly altered. This ultimately leads to the creation of conditions for tunneling that manifest themselves via the electrical percolation. It should be here noted that the nonlinear detection of the percolation threshold was performed at ω = 1 rad/s, with superior sensitivity compared to linear viscoelastic dynamic frequency sweeps, in accordance with similar studies.

To emphasize the shear rate amplitude dependence and the evolution of the nonlinear material response at the linear-nonlinear transition in strain sweep tests, the dependence of I3/1 on the shear stress amplitude \( {\dot{\gamma}}_0 \) at selected concentrations for ω = 0.6 rad/s is presented in Fig. 7a. The concentration dependence is illustrated for ω = 0.6 rad/s because it appears to be most susceptible to the presence of a filler network. Furthermore, the concentration range is further expanded beyond the data in Fig. 5 to capture the evolution of nonlinearities. The scaling indices, n, whereby \( {I}_{3/1}\propto {\gamma}_0^n \), identified beyond the instrumentation noise/SAOS in Fig. 7a, are summarized in Fig. 7b. As already mentioned in the description of Fig. 5, the quadratic scaling of I3/1 in the MAOS region is preserved for concentrations below percolation, i.e., n = 2. The percolation threshold marks a significant change in scaling, with n ≈ 0.7. Above the percolation point, n ≈ 0 at the linear-nonlinear transitions (SAOS-MAOS) followed by n ≈ 2 for higher strain/shear rate amplitudes. There is concurrently a decrease in the magnitude of I3/1 compared to ϕc, see also Fig. 6b. The same scaling is maintained until ϕ = 2.5 wt% together with an increase in I3/1 magnitude. In contrast, for ϕ = 3.5 wt%, three successive scaling regions could be identified, with n ≈ 1.1, n ≈ 0.2 and n ≈ 0.4. At this concentration, the nanocomposite is estimated to enter the conductor region of the percolation curve, see Fig. 3b. Similar studies have shown evidence of a percolated at the linear-nonlinear transition with similar scaling of nonlinearities (Gaska et al. 2017; Kádár et al. 2017; Gaska et al. 2019; Gaska and Kádár 2019). For a LCB/2D nanocomposite (LDPE/GnP)\( {I}_{3/1}\propto {\gamma}_0^{\approx 0} \) has been found after SAOS, followed by quadratic scaling (MAOS) (Gaska et al. 2017; Gaska et al. 2019; Gaska and Kádár 2019). Such scaling behavior could be comparable to PP/HrGO above ϕc, ϕ ∈ [1.5,2.5]. This could suggest the scaling is characteristic to the nonlinear dynamic response of a more robust filler network, i.e., having more connection points. The robust character of the network for LDPE/GnP could be due to the higher percolation threshold of the system (11 wt%) as and due to the morphology of the nanofiller (2D). In contrast, an important note here valid for the 3D hierarchical filler (HrGO) is that it can be considered as a “flexible” filler network due to the absence of physical connections between the individual graphene sheets, as described in the “Introduction” section. For a hybrid LCB/(2D,CB) nanocomposite (EBA/GnP,CB), it was shown that \( {I}_{3/1}\propto {\gamma}_0^{0.4} \), followed by a MAOS region (Kádár et al. 2017; Gaska and Kádár 2019). Due to the presence of CB in the composition, it could be inferred that the filler network could be more flexible and easy to distort. This could be a possible scenario for explaining the similarities in scaling between PP/HrGO and EBA/(GnP,CB), i.e., exponents 0.7 and 0.4, respectively, at percolation as well as why the EBA-based system reverts to a quadratic scaling at higher strain/strain rate amplitudes. However, for both the LDPE and EBA systems mentioned, the discussion disregards the role of the matrix topology that could influence the nonlinear materials response of the two systems. Furthermore, there could also be discrepancies between the error in the percolation thresholds compared, noting that ϕ = 1 wt% for PP/HrGO is very close to the threshold determined from the percolation theory.

Fig. 7
figure 7

a Third relative higher harmonic, I3/1, from dynamic strain sweeps for selected concentrations plotted as function of the strain rate amplitude, \( {\dot{\gamma}}_0 \), emphasizing the change in nonlinear material response at the linear-nonlinear transition at percolation and beyond. b Scaling indices, n, for the third relative higher harmonic, \( {I}_{3/1}\propto {\gamma}_0^n \) as function of filler concentration, ϕ, for regions identified beyond SAOS in (a)

The PP/HrGO nanocomposite elastic and viscous Lissajous-Bowditch (LB) diagrams, for a selected (ω, γ0) subset of the data in Fig. 5, are shown in Fig. 8. The material signatures are thus arranged in a Pipkin diagram for visual inspection. Linear viscoelastic material signatures (elliptic) are readily observable in both elastic and viscous LB diagrams, with nonlinear features apparent towards the higher end of (ω, γ0). The diagrams appear similar for concentrations below and above the percolation threshold, ϕ ≠ 1 wt% in Fig. 8. Markedly different material response could be observed at the percolation threshold, ϕ = 1 wt% in Fig. 8c, ϕ ≈ ϕc. Thus, the linear viscoelastic signature is already qualitatively different for ω < 1 rad/s, whereas for ω > 2 rad/s, both the linear and nonlinear signatures appear similar to results for ϕ ≠ 1 wt%.

Fig. 8
figure 8

Lissajous-Bowditch diagrams the a unfilled polypropylene, PP (ϕc = 0), b before the electrical percolation threshold, ϕ < ϕc, c at the electrical percolation threshold (ϕ = ϕc), and d directly above the percolation threshold. Top row, (a)–(d), represent the elastic Lissajous-Bowditch diagrams and bottom row, (e)–(h), represent the viscous Lissajous-Bowditch diagrams

Independently of the applied ω, v1 > 0, i.e., average viscous thickening nonlinear response. However, for ω < 1 rad/s, the nonlinear response is driven by large instantaneous strains, v3 < 0. In contrast, for ω > 1 rad/s, the nonlinear response is driven by large instantaneous rates-of-strain, v3 > 0. The nonlinear response at ω = 1 rad/s for ϕ = 1 wt% appears very similar ω = 0.6 rad/s for ϕ = 0.8 wt%. Interestingly, a similar behavior is observed above the percolation point at ϕ = 1.5 wt% (Fig. 9d). Thus, for ω ≤ 1 rad/s v3/1 < 0 and for ω > 1 rad/s v3/1 > 0. In addition, a γ0 dependence of v3/1 could be distinguished for ω ≤ 1 rad/s, with v3/1 > 0 at γ0 > 10% but changes sign at γ0 > 50%, i.e., intracycle strain-stiffening with increasing strains and shear thickening at intermediate strain (Du et al. 2018). For ϕ > 1.5 wt%, Fig. 10, a further increase in v3/1 takes place, with a strong relative relative increase for ϕ = 3.5 wt%. In addition, for ϕ = 3.5 wt%, v3/1(γ0) appears to show features in the linear-nonlinear transition similar to those observed from FT-Rheology, see Fig. 7.

Fig. 9
figure 9

Third relative elastic, e3/1, and viscous, υ3/1, Chebyshev coefficients for the a unfilled polypropylene, PP (ϕc = 0), b before the electrical percolation threshold, (ϕ < ϕc), c at the electrical percolation threshold (ϕ = ϕc), and d directly above the percolation threshold

Fig. 10
figure 10

Third relative elastic, e3/1, and viscous, υ3/1, Chebyshev coefficients for selected concentrations at ω = 0.6 rad/s

The nonlinear material response from Chebyshev polynomial decomposition could be discussed in terms of molecular chain conformational dynamics, confinement, and network orientation/fragmentations dynamics. Below percolation, ϕ < ϕc, the material nonlinear response is dominated by chain conformational dynamics and spatial arrangements. The effect of the shear flow in nonlinear oscillatory shear is expected to lead to a certain level of chain orientation, stretching, and disentanglement. The increased addition of fillers ultimately leads to strong confinement effects. Below percolation, there is an incomplete and unextended network. Therefore, it can be expected that the confinement effects would increase with ϕ. This appears to correspond to an increase in the strain-stiffening (e3/1 > 0, e3/1 ∝ ϕ) and a decrease in the shear thickening (v3/1 > 0, v3/1 ∝ 1/ϕ) nonlinear material response, as observed for ϕ < 0.8 wt%. This is in agreement with the gradual nature typically associated to the rheological percolation (Münstedt and Starý 2016). At ϕ = 0.8 wt%, the nonlinear signature appears sufficiently influenced by the filler concentration such that due to the altered chain conformation dynamics in the vicinity of the fillers that the shear thickening character is significantly decreased (v3/1 ≈ 0). Thus, there appears to be a similarity between v3/1 and the relative increase in Q0(ϕ)/Q0(0) at this filler concentration. While in this case elastic stiffening dominated the nonlinear material response, there could be a weak v3/1 = v3/1(γ0) inferred from the data, similar to the results of Du et al. (2018). Overall, the nonlinear response for ϕ = 0.8 wt% could be generically associated with the altered configurational characteristics of the polymer at the interface and in the interphase region, in comparison to the conformational characteristics of the polymer in the bulk. ϕ = 1 wt% marks the onset of a weakly connected network that could be disrupted easily with increasing ω, γ0. The ω-dependence is significantly featured in the viscous nonlinear response. For ω < 1, where theresponse is expected to be dominated by the weakly formed percolated filler network, the nonlinear signature is rate-of-strain dependent, v3/1 < 0. For ω > 1 as the weak network is expected to be discontinued and the chain confinement effects dominate, the nonlinear signature is strain dependent. Thus, the inferred nonlinear material signature at the low range of ω applied could be associated to a mechanical response of the network. Alternatively, this could be associated with parts of larger network formation structures that have been discontinued, as v3/1 is increasing with ω. For ϕ > ϕc in the semi-conductor region, the number of network connections is increased. This could be validated through the sign change in v3/1 that occurs at higher ω for ϕ = 1.5 wt% compared to ϕ = 1 wt%. In addition, the strain dependence observed for ϕ = 1.5 wt% could also be associated to a weakly formed network, similarly to the process assigned to the ω dependence reported in this publication. The behavior is qualitatively similar for all concentrations ϕ < 3.5 wt%. As the compositions approach the conductor region, a well-developed and strong network is expected which would result in a more complicated dynamic scenario. This fits within the framework of both the FT and Chebyshev decomposition results obtained for ϕ = 3.5 wt%. While, thus far the discussion has been overwhelmingly concentrated on conformational and interfacial aspects, the nonlinear dynamics could be dominated by the network orientation dynamics, e.g., see Natale et al. (2018). In a first approximation, it could be conjectured that due to the gradual nature of chain conformation and confinement, such effects could be perhaps associated dominantly to the elastic stiffening response. Meanwhile, network properties and their orientation/fragmentation dynamics to the viscous shear thickening nonlinear response.

Summary and conclusions

In this study, the percolation behavior of 3D hierarchical reduced graphene oxide (HrGO)-polypropylene (PP) nanocomposites has been evaluated from rheological point of view, with a focus on their nonlinear behavior. The electrical percolation threshold of the HrGO-PP nanocomposites (ϕc ≈ 0.94 wt%), Gkourmpis et al. (2019) is the lowest reported value in the scientific literature for similar systems prepared via melt mixing. An overview of the findings is presented in Table 2. Under the assumption that the nonlinear I3/1 scaling in MAOS should asymptotically be \( {I}_{3/1}\propto {\gamma}_0^2 \), a concurrent significant relative increase in Q0(ϕ)/Q0(0) was recorded at percolation. Thus, the percolation threshold was determined with increased sensitivity using nonlinear measurements compared to linear viscoelastic tests, in agreement with previous similar studies. However, in contrast to previous research, the electrical percolation threshold could be readily identified in the rheological data through a unique nonlinear signature. For ϕ = 1 wt%, ϕ ≈ ϕc, the nonlinear material response in dynamic strain sweep tests is dependent on the applied angular frequency, ω. This dependence was shown both in the third relative higher harmonic I3/1 scaling and in v3/1. Thus, at electrical percolation a \( {I}_{3/1}\propto {\gamma}_0^n \), n = n(ω) with n increasing with ω whereby for ω > 2 the typical MAOS scaling of n = 2 was observed, just as for concentrations below percolation. Correspondingly, a sign change in the third relative viscous Chebyshev coefficient, v3/1, i.e., from a intracycle viscous thickening behavior driven by large strains (v3/1 < 0) to being driven by large rates-of-strain (v3/1 > 0). As ω was increased I3/1 as well as v3/1 reverted to a behavior similar the concentrations before percolation. This could be potentially signaling the emergence of a weakly connected network that is being disrupted by the applied strains/strain rates. Thus, for supercritical concentrations, ϕ > ϕc, a scaling of \( {I}_{3/1}\propto {\gamma}_0^n \), with n = n(ω, γ0, ϕ) was observed. Similarly, v3/1 = v3/1(ω, γ0, ϕ). In particular, this behavior was inferred to be as a result of a more robust and developed filler network. Similar scaling in I3/1 has been reported for systems expected to have a stronger network at percolation due to the higher concentrations required to achieve it, e.g., see Gaska and Kádár (2019). Overall, the results could provide useful insights related to the detection of the percolation threshold in polymer nanocomposites, as well as into the filler network dynamics during flow.

Table 2 Breviary of experimental findings for the HrGO-PP nanocomposites investigated through nonlinear oscillatory shear experiments