Nonlinear “oddities” at the percolation of 3D hierarchical graphene polymer nanocomposites

The nonlinear rheology of a novel 3D hierarchical graphene polymer nanocomposites was investigated in this study. Based on an isotactic polypropylene, the nanocomposites were prepared using simple melt mixing, which is an industrially relevant and scalable technique. The novel nanocomposites stand out as having an electrical percolation threshold ( ≈ 0.94 wt%) comparable to solution mixing graphene-based polymer nanocomposites. Their nonlinear flow behavior was investigated in oscillatory shear via Fourier-transform (FT) rheology and Chebyshev polynomial decomposition. It was shown that in addition to an increase in the magnitude of nonlinearities with filler concentration, the electrical percolation threshold corresponds to a unique nonlinear rheological signature. Thus, in dynamic strain sweep tests, the nonlinearities are dependent on the applied angular frequency, potentially detecting the emergence of a weakly connected network that is being disrupted by the flow. This is valid for both the third relative higher harmonic from Fourier-transform rheology, I 3/1 , as well as the third relative viscous, v 3/1 , Chebyshev coefficient. The angular frequency dependency comprised non-quadratic scaling in I 3/1 with the applied strain amplitude and a sign change in v 3/1 . The development of the nonlinear signatures was monitored up to concentrations in the conductor region to reveal the influence of a more robust percolated network.


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 Electronic supplementary material The online version of this article (https://doi.org/10.1007/s00397-020-01204-w) contains supplementary material, which is available to authorized users. 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 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 γ̇t ð Þ ¼ γ 0 ωcos ωt ð Þ, with γ̇0 ¼ γ 0 ω as the shear rate amplitude. The Fourier spectrum of the shear stress, Σ 12 , therefore contains exclusively the fundamental intensity I 1 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, I n , in addition to I 1 , Fig. 1b. By comparing the shear stress, σ 12 , expressed from a Taylor series expansion of the shear viscosity, η γÞ ð , 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 where I 1 , I 3 , I 5 ... 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, I 3 /I 1 = I 3/1 , as it contains the dominant nonlinear contribution. The variation of I 3/1 during a strain sweep test for a typical polymeric material exemplified by t he 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 ∝γ 2 0 , region called medium amplitude oscillatory shear (MAOS) or intrinsic large amplitude oscillatory shear ([LAOS]). Based on the quadratic I 3/1 scaling in MAOS, the Q-parameter can be defined as nonlinear material parameter from FT-Rheology (Hyun and Wilhelm 2009) 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): where T n is the nth-order Chebyshev polynomial of the first kind, and e i and v i 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. γ̇t ð Þ=γ̇0, respectively. The red/blue lines represent the total elastic/viscous stresses, i.e., σ e /σ max vs. γ(t)/γ 0 and σ v (t)/σ max vs. γ̇t ð Þ=γ̇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 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 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): e 1 > 0 average elastic stiffening driven by e 3 > 0 large instantaneous strains e 3 < 0 large instantaneous rates-of-strain e 1 < 0 average elastic softening driven by e 3 > 0 large instantaneous rates-of-strain e 3 < 0 large instantaneous strains v 1 > 0 average viscous thickening driven by v 3 > 0 large instantaneous rates-of-strain v 3 < 0 large instantaneous strains v 1 < 0 average viscous thinning driven by v 3 > 0 large instantaneous strains v 3 < 0 large instantaneous rates-of-strain A similar interpretation can be assigned to the strainstiffening (S) and shear thickening (T) dimensionless ratios (Ewoldt et al. 2008): 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 Examples of linear (black line) and nonlinear (gray line) a elastic and b viscous Lissajous-Bowditch diagrams and the defined associated moduli rubber-based composites and nanocomposites, e.g., see Leblanc (2006a, b) and Schwab et al. (2016). The range of polymernanofiller 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 ).
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 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 I 3/1 , often quantified via the Q-parameter. Exceptions are 0D nanofillerbased systems, e.g., (L/SCB)/0D , 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 Q 0 with increasing filler content before the percolation threshold has also been reported (Ahirwal et al. 2014). The material nonlinear response as quantified by I 3/1 has been shown to be sensitive to the filler dimensionality 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 I 3/1 , dI 3/1 /dγ 0 = 0, was apparent in the measuring range for 1D  and 2D nanofillers Du et al. 2018). Particular changes in I 3/1 scaling have also been reported for 2D Gaska et al. 2017; and 2D hybrid nanocomposites , i.e., I 3=1 ∝γ n 0 with n ≠ 2 in longchain branched (LCB) polymers. Overall, most nonlinear effects in the percolation region were attributed to dispersion  ) ) (Hyun et al. 2012)  properties 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  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
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 M w = 300 kg/mol and M w /M n = 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.

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 plateplate 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 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 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, , 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, , where K and n are the consistency and flow behavior indices, appropriately describing the shear thinning behavior forγ > 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.
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 ω, I 3/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 I 3 and I 1 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 I 3/1 should be perceived through the viewpoint of filler network properties.
Changes in nonlinear I 3/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 I 3/1 with γ 0 as deduced from Eq. 2 and thus disregards any possible anomalies at the linear-nonlinear transition. The Qparameter 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 being the corresponding fitting parameters. Therefrom, the zero-strain Q-parameter can be extrapolated, . The relative zero-strain Q-parameter, Q 0 (ϕ)/Q 0 (ϕ = 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, Q 0 (ϕ)/Q 0 (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 Q 0 (ϕ)/Q 0 (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 . b Normalized zero-strain Q, Q 0 = 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) 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 linearnonlinear transition in strain sweep tests, the dependence of I 3/1 on the shear stress amplitudeγ 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 ∝γ n 0 , 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 I 3/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 I 3/1 compared to ϕ c , see also Fig. 6b. The same scaling is maintained until ϕ = 2.5 wt% together with an increase in I 3/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;. For a LCB/2D nanocomposite (LDPE/GnP)I 3=1 ∝γ ≈0 0 has been found after SAOS, followed by quadratic scaling (MAOS) (Gaska et al. 2017;. 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 ∝γ 0:4 0 , followed by a MAOS region . 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.
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 (b) Fig. 7 a Third relative higher harmonic, I 3/1 , from dynamic strain sweeps for selected concentrations plotted as function of the strain rate amplitude,γ 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 ∝γ n 0 as function of filler concentration, ϕ, for regions identified beyond SAOS in (a) 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%.
Independently of the applied ω, v 1 > 0, i.e., average viscous thickening nonlinear response. However, for ω < 1 rad/s, the nonlinear response is driven by large instantaneous strains, v 3 < 0. In contrast, for ω > 1 rad/s, the nonlinear response is driven by large instantaneous rates-of-strain, v 3 > 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 v 3/1 < 0 and for ω > 1 rad/s v 3/1 > 0. In addition, a γ 0 dependence of v 3/1 could be distinguished for ω ≤ 1 rad/s, with v 3/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 v 3/1 takes place, with a strong relative relative increase for ϕ = 3.5 wt%. In addition, for ϕ = 3.5 wt%, v 3/1 (γ 0 ) appears to show features in the linear-nonlinear transition similar to those observed from FT-Rheology, see Fig. 7. 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 strainstiffening (e 3/1 > 0, e 3/1 ∝ ϕ) and a decrease in the shear thickening (v 3/1 > 0, v 3/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 (v 3/1 ≈ 0). Thus, there appears to be a similarity between v 3/1 and the relative increase in Q 0 (ϕ)/Q 0 (0) at this filler concentration. While in this case elastic stiffening dominated the nonlinear material response, there could be a weak v 3/1 = v 3/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, v 3/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 v 3/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 v 3/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 e 3/1 > 0; v 3/1 > 0 ϕ = 0.8 wt% MAOS: I 3=1 ∝γ 2 0 weak relative increase in Q 0 (ϕ)/Q 0 (0) e 3/1 > 0; v 3/1 ≈ 0, especially for ω = 0.6 rad/s ϕ = 1 wt% MAOS: I 3=1 ∝γ n 0 , n = n(ω) (ϕ ≈ ϕ c ) n ≈ 0.7 for ω = 0.6 rad/s; n ≈ 1.1 for ω = 1 rad/s; n ≈ 2, for ω > 2 rad/s significant relative increase in Q 0 (ϕ)/Q 0 (0) e 3/1 > 0; v 3/1 = v 3/1 (ω) v 3/1 < 0 for ω < 1 rad/s (v 3 < 0) v 3/1 > 0 for ω > 2 rad/s (v 3 > 0) v 3/1 weak dependence on γ 0 1 < ϕ < 3.5 wt% SAOS-MAOS: I 3=1 ∝γ n 0 , n = n(ω, γ 0 , ϕ) e.g., ϕ = 1.5 wt%: for ω < 2 rad/s, γ 0 ∈ [1, 10] %, n ≈ 0; γ 0 ∈ [10, 40] %, n = 2 for ω ≥ 2, → n = 2 e 3/1 > 0; v 3/1 = v 3/1 (ω, γ 0 , ϕ): e.g., ϕ = 1.5 wt%: for ω < 2 rad/s, γ 0 ∈ [10, 40] %, v 3/1 < 0; γ 0 > 40%, v 3/1 > 0 for ω > 2, v 3/1 > 0 ϕ = 3.5 wt% I 3=1 ∝γ n 0 ; with increasing γ 0 : n ≈ 1.1, n ≈ 0.2 and n ≈ 0.4 significant relative increase in e 3/1 , v 3/1 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 I 3/1 scaling in MAOS should asymptotically be I 3=1 ∝γ 2 0 , a concurrent significant relative increase in Q 0 (ϕ)/Q 0 (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 I 3/1 scaling and in v 3/1 . Thus, at electrical percolation a I 3=1 ∝γ n 0 , 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, v 3/1 , i.e., from a intracycle viscous thickening behavior driven by large strains (v 3/1 < 0) to being driven by large rates-of-strain (v 3/1 > 0). As ω was increased I 3/1 as well as v 3/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 ∝γ n 0 , with n = n(ω, γ 0 , ϕ) was observed. Similarly, v 3/1 = v 3/1 (ω, γ 0 , ϕ). In particular, this behavior was inferred to be as a result of a more robust and developed filler network. Similar scaling in I 3/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 . 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.