Modelling electrical transport innano-particle composites using 3Dnetwork analysisbyKaan Owen WilliamsB.A.Sc., The University of British Columbia, 2009A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF APPLIED SCIENCEinThe Faculty of Graduate and Postdoctoral Studies(Electrical & Computer Engineering)THE UNIVERSITY OF BRITISH COLUMBIA(Vancouver)February 2014? Kaan Owen Williams, 2014AbstractConductive composites consist of a conductive filler dispersed within an in-sulating matrix. These composite materials have been known for many yearsand are regularly produced experimentally and commercially for a variety ofapplications. Novel techniques are now being found for creating compositesthat exhibit conductivity with less conductive filler material than classicalphysics suggests is sufficient if the particles are uniformly distributed. Sev-eral parties have offered physical explanations for the characteristics of theircomposites by incorporating a blend of classical and quantum physics butfew attempts have been made to compare explanations or develop any mech-anism to simulate the physics. The model presented in the present workincorporates first principles physics and semi-empirical theory to accountfor the distribution of particles within a composite and calculate resultantconductivity using three dimensional network analysis. Results from severalmodel iterations are presented and they are compared with published ex-perimental results. The model demonstrates that a random distribution ofspherical particles smaller than 200 nm at 3% loading, given realistic wavefunction decay rates and reasonable tunnelling barrier heights, cannot ex-plain experimentally observed conductivities in these composite materials.The final model, using a Voronoi tessellation approach, duplicates the be-haviour trend of the composites being simulated and illustrates some gapsin the present material science knowledge of conductive composites.iiPrefaceEarly work on the model presented herein was performed by Liam Russel in2010. His initial MATLAB embodiment was limited to the particle place-ment method described in Section 2.1.1 and was limited to simulations in-volving relatively few particles within relatively small volumes. Neverthelesshis framework established an approach that was substantially continued asthe present model was developed. In cases where the algorithms developedby Russel did not need logical adjustment (such as the particle connectiv-ity mapping described in Section 2.2) his code was re-written or adaptedto improve memory usage and computation time at larger scales. The onlycomponent that remained substantially unchanged is used to graphically dis-play model output. That component proved useful during the work discussedherein although it was not used to generate any of the figures presented.Figures 1.7 and 1.9 are used with permission from applicable sources.Figure 1.6 was generated by me to show a subset of the information presentedin a similar figure published by Kohjiya et al. and he is credited in thatfigure?s caption.A series of material analyses performed by Aaron Linklater cited in Fig-ure 1.8 and elsewhere are unpublished and were shared through an internalUBC report.iiiTable of contentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiTable of contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivList of tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiList of symbols . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiiList of abbreviations . . . . . . . . . . . . . . . . . . . . . . . . . . xivAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . xvDedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xvi1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Percolation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Tunnelling conduction . . . . . . . . . . . . . . . . . . . . . . 31.2.1 Sphere to sphere tunnelling . . . . . . . . . . . . . . . 31.2.2 Hard/soft shell approach . . . . . . . . . . . . . . . . 101.3 Particle agglomeration and filament formation . . . . . . . . 121.4 Carbon black as a conductive filler . . . . . . . . . . . . . . . 121.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 182 Experimental methods . . . . . . . . . . . . . . . . . . . . . . 192.1 Particle placement . . . . . . . . . . . . . . . . . . . . . . . . 212.1.1 Method 1: random placement of uniform particles . . 212.1.2 Method 2: clustering . . . . . . . . . . . . . . . . . . 212.1.3 Method 3: volume exclusion method . . . . . . . . . 222.1.4 Method 4: Voronoi wire frame method . . . . . . . . 23ivTable of contents2.1.5 Method 5: distributed radii . . . . . . . . . . . . . . . 252.1.6 Method 6: intersections allowed . . . . . . . . . . . . 262.2 Particle connectivity mapping . . . . . . . . . . . . . . . . . 282.3 Electrodes . . . . . . . . . . . . . . . . . . . . . . . . . . . . 292.4 Node trimming . . . . . . . . . . . . . . . . . . . . . . . . . . 293 Inter-particle conductivity . . . . . . . . . . . . . . . . . . . . 323.1 Absolute path length . . . . . . . . . . . . . . . . . . . . . . 323.2 Basic exponential . . . . . . . . . . . . . . . . . . . . . . . . 323.3 Conduction due to particle contact . . . . . . . . . . . . . . . 333.4 Conductivity summary . . . . . . . . . . . . . . . . . . . . . 344 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 374.1 Model size . . . . . . . . . . . . . . . . . . . . . . . . . . . . 374.2 Random placement of uniform particles . . . . . . . . . . . . 374.3 Contact conduction . . . . . . . . . . . . . . . . . . . . . . . 394.4 Clustering and filaments . . . . . . . . . . . . . . . . . . . . 434.5 Voronoi cell network . . . . . . . . . . . . . . . . . . . . . . . 505 Conclusions and future work . . . . . . . . . . . . . . . . . . . 565.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 565.2 MATLAB as a modelling tool . . . . . . . . . . . . . . . . . 575.3 Future work . . . . . . . . . . . . . . . . . . . . . . . . . . . 585.3.1 Greater complexity in modelling parameters . . . . . 585.3.2 Post-curing analysis of carbon black . . . . . . . . . . 585.3.3 Alternate techniques . . . . . . . . . . . . . . . . . . 59Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60AppendixA Simulation parameters . . . . . . . . . . . . . . . . . . . . . . . 66vList of tables1.1 Comparison of values for ? and ?. . . . . . . . . . . . . . . . 8viList of figures1.1 Figure 1.1A is an illustration showing the basic power lawrelationship described by combining equation (1.1) (solid line)and equation (1.2) (dashed line). In this illustration ?c =0.16 and t = 2.0 are used and ? is scaled for visual clarity.Figure 1.1B does not differentiate between equation (1.1) andequation (1.2) but shows the effect of scaling t. For the dottedline t = 1.5, for solid line t = 2, and for the dashed line t = 2.5. 21.2 Equation 1.3 is represented by 1R?as a function of ? withlines shown for four different values of ?. In this exampleA = ?(10 nm)2. . . . . . . . . . . . . . . . . . . . . . . . . . . 51.3 Comparison of potential barrier height and the constant ?according to equation (1.4). Recall that potential barriersin conductive composites are typically in the range 0.05 eV> ? > 1 eV [1]. The boundaries of that range are illustratedwith dotted lines for reference. . . . . . . . . . . . . . . . . . 61.4 Comparison of potential barrier height (?) and the inter-particle gap over which tunnelling must typically occur (?)using the modelling approach from Wang et al. [2] and Meieret al. [3]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91.5 The normalised critical distance, ?c2r , is shown with respect tovolume loading fraction of solid core spheres. The relationshipis described explicitly by equations (1.15) and (1.16) takenfrom [4]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111.6 Typical distance between nearest-neighbour carbon black ag-glomerates in natural rubber as a function of volume loadingfraction. Adapted from [5]. . . . . . . . . . . . . . . . . . . . 141.7 Scanning electron microscope micrograph of Vulcan XC72.Reprinted from [6] with permission from Elsevier. . . . . . . . 15viiList of figures1.8 Conductivity across a range of loadings of two different car-bon blacks in four different polymers. Taken from works pub-lished or shared by Niu et al. [7], Rwei et al. [8], Linklater [9],Li [10], Huang [11], and Huang et al. [12]. . . . . . . . . . . . 171.9 Size distribution of Vulcan XC72. Data represents about660 particles in about 200 images. Figure created by [13].Reprinted by permission of The Electrochemical Society. . . . 182.1 Processing the present model involves six steps that are de-scribed in the text. . . . . . . . . . . . . . . . . . . . . . . . . 202.2 A two dimensional arbitrary unit representation of the volumeexclusion method described in Section 2.1.3. Exclusion zonesare shown as circles with a dashed line (radius = 0.05). Par-ticles are shown as a circle with a solid line (radius = 0.01).In this example, exclusions zones are permitted to overlap,particles are not, and no particles may be in an exclusion zone. 232.3 A 2D Voronoi structure with cell edges shown as lines anddots at each seed point. . . . . . . . . . . . . . . . . . . . . . 242.4 Stereographic example of a 3D Voronoi structure with celledges shown as lines and dots at each vertex. Seed points areshown as a star. Hold image about 50 cm away and alloweyes to cross slightly. When three images appear, observe thecentral one. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 252.5 In this example set of particles the user set target mean ra-dius and standard deviation were 65.0 nm and 9.75 nm re-spectively. A log-normal distribution with those parametersis shown as a dashed line for reference. It has a mode radiusof about 62.8 nm. The resultant particles from four succes-sive models are combined and histogram data of the resultantradius distribution is shown as the solid line. Histogram binsspan about 0.32 nm. Particle radii vary from 31.2 nm to 124.5nm with a mean and standard deviation of 64.3 nm and 9.55nm respectively. . . . . . . . . . . . . . . . . . . . . . . . . . . 27viiiList of figures2.6 In this two dimensional arbitrary unit example conductiveparticles are shown as ? and electrodes are shown as *. Con-nections between sites (particles or electrodes) are shown asa line. Diagram A shows a connected network that has fewelectrode connection sites. Diagram B shows the same net-work but with more electrode connection sites. Observe thedifference between the two diagrams in particular near x-axislocations 0.2, 0.5, and 0.8. . . . . . . . . . . . . . . . . . . . . 302.7 This two dimensional arbitrary unit example contains thesame connected network as Figure 2.6B. Conductive parti-cles are shown as ? and electrodes are shown as *. Connec-tions between sites (particles or electrodes) are shown as aline. Three particle sites that do not contribute useful elec-trical nodes have been identified (with circles) as suitable forexclusion from the conductivity calculation. . . . . . . . . . . 313.1 Geometry between overlapping spheres, ? < 0. When ? > 0,rc = 0, and the spheres do not overlap. . . . . . . . . . . . . . 343.2 Conductance between particles arranged according to increas-ing gap size. Particle size varies as illustrated in Figure 2.5and the effect of this variation can be seen as a small amountof noise in the line. A negative gap indicates particles overlapand the particles are therefore in physical contact. . . . . . . 364.1 Conductivity is shown to increase asymptotically to a max-imum value as the edges length of a cube being simulatedincreases. Individual results are shown as a circle. The long-dash line indicates mean value and the dotted lines are onestandard deviation above and below the mean. The solid lineshows the exponential increase in computation time as theedge size increases. . . . . . . . . . . . . . . . . . . . . . . . . 384.2 Particles that are deemed ?connected? will conduct chargeaccording to equation (3.2). The amount of conduction thatoccurs depends on the gap between particles. The right axis(and thin line) show histogram counts for gaps of a particularlength from a typical execution of the model. Histogram binsspan about 5.5 nm. For reference, the left axis (and thickdashed line) show the results of the tunnelling equation (3.2)for a gap of that size. . . . . . . . . . . . . . . . . . . . . . . . 40ixList of figures4.3 Five unique models were created with the same initial pa-rameters (described in the text). Once each model was cre-ated, the radius of interaction was varied to create maps ofconnected particles. If a model exhibited end-to-end connec-tivity of its particles at a particular radius of interaction, aconductivity value for that combination is shown. Data tothe left of the vertical dashed line is intermittent and data tothe right of the vertical dashed line begins to converge. . . . . 414.4 Particles with a negative inter-particle distance are said to bein contact and conduct according to equation (3.11). Parti-cles not in physical contact that are deemed ?connected? willconduct charge according to equation (3.2). The right axis(and thin line) show histogram counts for gaps of a particu-lar length from a typical execution of the model. Histogrambins span about 5.8 nm. The left axis (and thick dashed line)show the results of the contact equation (3.11) and tunnellingequation (3.2) for a gap of that size. . . . . . . . . . . . . . . 424.5 Five unique models were created with nearly all the sameinitial parameters to generate each set of results (each set isshown with a different symbol). The only difference betweeneach set of model parameters is the numerical value for ?(indirectly defined by setting the barrier height to three dif-ferent values). For ?, ? = 0.1 meV. For ?, ? = 1 ?eV. For?, ? = 10 neV. Once each model was created, the radius ofinteraction was varied to create maps of connected particles.If a model exhibited end-to-end connectivity of its particlesat a particular radius of interaction, a conductivity value forthat combination is shown. Data to the left of the verticaldashed line is noisy and intermittent and data to the right ofthe vertical dashed line begins to converge. . . . . . . . . . . 444.6 Three hundred particles are distributed in two dimensionalspace (arbitrary units). Seventy five of them are placed ran-domly and the remainder are placed in succession where eachmust be within 0.5 units of a previously placed particle. Eachparticle is shown as ? and its radius of interaction is shownas ? (ri ? 0.015 units). For the purposes of this illustration,consider particles connected if their ??s touch. The result issimilar to a random distribution. Some connections exist butrarely involving many particles. . . . . . . . . . . . . . . . . . 46xList of figures4.7 Three hundred particles are distributed in two dimensionalspace (arbitrary units). Seventy five of them are placed ran-domly and the remainder are placed in succession where eachmust be within 0.05 units of a previously placed particle.Each particle is shown as ? and its radius of interaction isshown as ? (ri ? 0.015 units). For the purposes of this illus-tration, consider particles connected if their ??s touch. Severalgroupings of many particles are apparent but the groupingsdo not often interconnect. . . . . . . . . . . . . . . . . . . . . 474.8 Conductivity values as a function of radius of interactionwhen an enforced proximity of 25 nm is implemented. Eachdata point represents the mean result of five simulations andeach line is associated with a different enforcement probabil-ity as shown in the legend. Overall volume loading of particlesis held constant at ? = 0.03. As described in Section 2, con-ductivities near 10?9 S/m are of unreliable precision due toa near-singularity during matrix inversion. This accounts forthe erratic low mean conductivity values. . . . . . . . . . . . 484.9 Conductivity values as a function of radius of interactionwhen an enforced proximity of 100 nm is implemented. Eachdata point represents the mean result of five simulations andeach line is associated with a different enforcement probabil-ity as shown in the legend. Overall volume loading of particlesis held constant at ? = 0.03. . . . . . . . . . . . . . . . . . . . 494.10 Histogram data of typical mean edge lengths for a model witha target density of vdens = 5 ? 106 mm?3. Each histogrambin spans about 0.95 ?m. . . . . . . . . . . . . . . . . . . . . 504.11 Conductivity as a function of volume loading for a range offilament radius values. . . . . . . . . . . . . . . . . . . . . . . 524.12 Model results are shown here with published experimentaldata [8, 9]. (Individual attributions are shown in Figure 1.8.)Model results are those previously shown in Figures 4.3, 4.5,4.9, and 4.11. For reference, the upper limit of conductivityis also shown. . . . . . . . . . . . . . . . . . . . . . . . . . . . 534.13 Voronoi model results are compared to a line fit for the powerlaw behaviour discussed in Section 1.1. . . . . . . . . . . . . . 544.14 Conductivity as a function of radius of interaction for fourdifferent volume loadings. Data to the lower-left of the dashedline is noisy and intermittent in end-to-end connectivity anddata to the upper-right of the dashed line is convergent. . . . 55xiList of symbolsSymbol Units DescriptionA m2 effective cross-sectional area where tunnelling occursAc m2 cross-sectional area of particle contacta m distance between sphere-centresa0 fitting constanta1 fitting constanta2 fitting constanta3 fitting constantC? F capacitance between two sites separated by ?f denotes a function defined elsewheref0 s?1 resonant frequencyG S conductanceGt S inter-particle tunnelling conductanceG0 S base conductanceG1 S filler conductanceGc S contact conductanceme kg electron mass ? 9.11 ? 10?31n number of particlesqe C elementary charge ? 1.60 ? 10?19R? ? resistance between two sites separated by ?R0 ? characteristic resistancer m radius of a spherical conductive particler1 m radius of sphere oner2 m radius of sphere tworc m radius of the contact circle made by two overlapping spheresri m radius of interactiont critical exponent of conductanceV m3 volumevdens m?3 quantity of Voronoi cells per unit volumevrMax m maximum thickness of filaments in Voronoi cells? m distance between conductive particlesxiiList of symbols?c m critical distance? m?1 wave function decay rate? relative permittivity?0 A2 S4 kg?1 m?3 vacuum permittivity ? 8.85? 10?12? rad angle considered in Figure 3.1? m characteristic tunnelling length? ratio of a circle?s circumference to its diameter ? 3.14? S m?1 conductivity?c S m?1 internal conductivity of conductive particles?? S m conductivity of a fixed cross-sectional area? volume fraction of conductive particles within an insulator?c critical volume fraction of conductive particles within an insulator? eV potential barrier height~ m2 kg s?1 reduced Planck constant ? 1.05 ? 10?34% parts per hundredxiiiList of abbreviationsAbbreviation Description2D two dimensional3D three dimensionalA ampC coulombcm centimetrecos cosineeV electron-voltF faradGB gigabyteGHz gigahertzHz hertzkg kilogramk? kiloohmm metremeV millielectron-voltmm millimetreneV nanoelectron-voltnm nanometreRAM random access memoryrad radianS siemenss secondsin sine?eV microelectron-volt?m micrometre?S microsiemens? ohm? degree CelciusxivAcknowledgementsI would like to express special appreciation and thanks to my supervisors DrKonrad Walus and Dr Boris Stoeber. Your support, guidance, and patiencehelped me tremendously at all stages of my post-graduate work. I extendmy sincere thanks also to Dr John Madden for providing feedback during theearly stages of this work and serving on my defence committee. All three ofyou provided a tremendous amount of support for my education and researchduring both my undergraduate and post-graduate studies. You have grantedme opportunities and insights that have been both unique and valuable.I would also like to thank Henryk Modzelewski for his assistance withthe computing cluster used to perform most of the simulations that drovethis research and for reminding me that just because I sometimes get theresults I expect does not mean they are necessarily correct.Thank you to Shabnam Shambayati for her five rules to surviving gradschool and to Sheng Pan for ensuring my workspace was always a good placeto share adventure stories. May you live in interesting times.Without the support of Arwen, Theo, Tess, Dylan, and Sara this mightnever have happened. Thank you.xvDedicationEngineering is only as useful as those it serves. While this endeavour wasintellectually intriguing in a direct way, it also positioned me to build manyof the new relationships that now fuel me emotionally. This work is dedi-cated to those who have fuelled me and helped me better serve.Foremost: Tess Baker, Sara Pour, Mark Abbott, and Sean McHugh.What next?xviChapter 1BackgroundElectrically conductive composites can be created through the dispersion ofconductive filler particles throughout an insulating matrix. Potential fillerparticles include carbon allotropes such as carbon black or nanotubes; or-ganic conductors such as polypyrrole; and metals. Charge flow within theseconductive composites has been attributed to a variety of transport mecha-nisms. A definitive trend shared by all the composites is that a compositewith very little conductive filler has a conductivity near that of the insulat-ing material and the conductivity increases toward that of the conductivefiller as the filler is added. The current theory describing the physics thatdrive that conductivity change are explored further in the following sections.1.1 PercolationIn the simplest case, conduction can be said to follow basic percolationtheory whereby particles must touch in order for charge transport to oc-cur [14, 15]. Net conductivity is then dependent on the total number ofphysically connected pathways between test points. At low loading volumesof conductive particles the mean distance between filler particles is large andthe conductance of the composite is near that of the insulating matrix. Asthe volume loading of conductive particles, ?, is increased, the mean dis-tance between filler particles decreases and they begin to form linkages. Aslinkages are sufficiently formed to create a conductive pathway through thematrix the conductivity of the particles begins to contribute to the conduc-tivity of the bulk. At a critical loading fraction, ?c, there is a sharp increasein conductivity, ?, as the conductive pathways reach sufficient abundanceas to be common-place [16].The critical loading fraction in a simple model such as this can be cal-culated statistically. A random filling routine is used to place particles inthree-dimensional space. The critical loading fraction is taken to be themean loading fraction at which 50% of these models exhibit a percolatingcluster [4]. At this critical loading fraction, the material is said to be atits percolation threshold. At loadings near the critical loading fraction, the11.1. Percolationconductivity can change drastically (over several order of magnitude) withsmall changes in the loading fraction [14, 16]. This simple percolation-styleconductivity model follows a power-law behaviour of the form? ? (?c ? ?)?t, (1.1)where t is the critical exponent and ? < ?c. A continued increase in ? beyond?c results in a gradual plateau of ? that approaches the bulk conductivityof the filler particle. This gradual plateau also continues the power-lawbehaviour but with some signs reversed.? ? (?? ?c)t, (1.2)where ? > ?c. The shape of this generalised relationship is shown in Fig-ure 1.1. Simple percolation theory predicts a volume loading fraction of?c = 0.16 and t = 2.0 for spheres packed cubically in 3D space and isnot dependent on the size of the spheres or the size of the volume beingconsidered [14, 17].0 0.1 0.2 0.310?2010?1510?1010?5100??(S/m)A0 0.1 0.2 0.310?2010?1510?1010?5100??(S/m)BFigure 1.1: Figure 1.1A is an illustration showing the basic power lawrelationship described by combining equation (1.1) (solid line) and equa-tion (1.2) (dashed line). In this illustration ?c = 0.16 and t = 2.0 are usedand ? is scaled for visual clarity. Figure 1.1B does not differentiate betweenequation (1.1) and equation (1.2) but shows the effect of scaling t. For thedotted line t = 1.5, for solid line t = 2, and for the dashed line t = 2.5.When experimental data is reviewed the simple percolation theory isquickly shown to be insufficient. Empirical analysis gives a critical expo-nent in the range of t = 1.3 to t = 20 [18?22]. There are conflicting views on21.2. Tunnelling conductionwhich physical properties influence which variables. Balberg, working withcarbon black in polyethylene, claims ?c and t are systematically related tothe structure of the carbon black under test, as well as the moulding con-ditions [17]. Lower structure (more spherical) results in higher t values. Asstructure becomes higher (less spherical), t approaches a value of 2.0. Thisis attributed to the particle contact requirement of the simple percolationtheory [17]. Others [14, 16] claim ?c of the composite depends on the shapeof the particles and the homogeneity of the dispersion while the critical ex-ponent t depends on the material combination. Low values of ?c have alsobeen attributed to the formation of long conductive filaments and other in-homogeneities in the filler distribution [18]. Some [16, 23, 24] have made thespecific claim that t values are not at all correlated to particle shape or ?c.Clearly there is a lack of consensus as one moves from the purely theo-retical model involving perfect spheres in stochastically predictable arrange-ments to laboratory devices involving irregular shapes, sizes, and distribu-tions. Much of this discord is rooted in the range of conduction influencersin situations where the particles seem unlikely to always be in contact. Inthese cases we turn to quantum tunnelling for additional insight.1.2 Tunnelling conduction1.2.1 Sphere to sphere tunnellingWhen the dispersed particles are very small or have high aspect ratios, signif-icant conduction is observed at volume loadings lower than the percolationthreshold predicted by basic percolation theory. This is often attributed toa quantum tunnelling conduction mechanism whereby particles need not bein physical contact for charge transport to occur. The amount of chargetransported between conductors that are not in contact falls off exponen-tially as the distance, ?, between them increases and is typically considerednegligible between conductors more than ? = 10 nm apart [16]. Despitethis general consensus, there are published findings of tunnelling conductionoccurring over gaps as large as 2.5 ?m [25].The amount of charge transported is also heavily dependent on the heightof the potential barrier that occurs at the insulator-conductor junction. Thepotential barrier height, ?, between adjacent particles in composites com-prised of a metallic conductor and a polymeric insulator is typically in therange 0.05 eV > ? > 1 eV [1].31.2. Tunnelling conductionThe resistance, R?, between two sites is often given asR? =16?2~?3A?q2ee?? , (1.3)where? = 2~?2me?, (1.4)~ is the reduced Planck constant, A is the effective cross-sectional area wheretunnelling occurs (between spherical particles), and qe and me are the chargeand mass of an electron respectively [2, 26]. In these analyses ? is gener-alised as the typical distance between particles in a large system with manyparallel branches. A plot of equation 1.3 is shown in Figure 1.2. It is ap-parent that for the ? values being considered particles must be very closein order to register a meaningful conductance between them. Although thedevelopment of this equation describes A as the effective cross-sectional areawhere tunnelling occurs between spherical particles, a rigorous definition of?effective cross-sectional area? is not provided. The tunnelling conductancebetween spheres is going to be different than the tunnelling conductance be-tween parallel plates of the same surface area (or projection area). A 1963analysis of inter-particle conduction [27?30] that leads to equation (1.3) istypically at the root of modern analyses of this topic in the literature. Theequation is of the same form as seen in analyses of scanning tunnelling mi-croscopy although in those analyses a parallel plate assumption is typicallyeither explicit or implied, even when tunnelling is occurring from a mono-atomic needle tip. Point-to-point tunnelling approximations result in a ne-glecting of the terms ?A? (as is done later in equation (3.3)) [31?33]. Thesediscrepancies in geometry assumptions nevertheless remain unresolved.? is defined from first principles as a description of wave function decaywithin the tunnelling barrier. The tunnelling barrier is dependent on theeffective local work function [33]. An intuitive way to consider ??1 is totreat it as a number representative of the distance over which tunnellingoccurs (equivalent to exploring the case where ?? = 1). The location of ? inthe exponential (1.3) gives it the units m?1 in agreement with the descrip-tion in (1.4). Figure 1.3 shows the relationship described by equation (1.4)but inverts ? to better serve its consideration as an indicator of tunnellingdistance.As can be extracted from the line in Figure 1.3, the typical boundaryvalues for ? stated earlier, 0.05 eV and 1 eV, give values for ??1 of 0.437 nmand 0.0976 nm respectively. Despite the agreement of units, the valuesobtained seem too small to be a literal representation of typical tunnelling41.2. Tunnelling conduction1 2 3 4 5 610?2010?1510?1010?5100? (nm)1 /R ?(S) ? = 0.05 ev? = 0.14 ev? = 0.37 ev? = 1.00 evFigure 1.2: Equation 1.3 is represented by 1R?as a function of ? with linesshown for four different values of ?. In this example A = ?(10 nm)2.51.2. Tunnelling conduction0.0001 0.001 0.01 0.1 1 100.010.1110? (eV)1 /?(nm)Figure 1.3: Comparison of potential barrier height and the constant ? ac-cording to equation (1.4). Recall that potential barriers in conductive com-posites are typically in the range 0.05 eV > ? > 1 eV [1]. The boundariesof that range are illustrated with dotted lines for reference.61.2. Tunnelling conductiondistances. For reference, scanning tunnelling microscopy typically operatesover a distance of 0.7 nm to 0.3 nm, equivalent in this analysis to a barrierheight range of 0.02 eV to 0.1 eV [33, 34]. There is only a small amountof overlap in those ranges suggesting a direct comparison between ? andtunnelling distance may not be always fair. Nevertheless, this definitiondoes offer a starting point for a model that requires a value for ? whenmapping connections between particles.Numerical values for ? and ? are not typically available for a preparedsample of conductive composite. In some cases they are treated as constantsuseful for fitting a trend line to experimental results. This approach is usedby Wang et al. [2] in a composite of polydimethylsiloxane and carbon blackwhere a unitless combined value of ?? = 1.786 is used to fit a trend line.Combining (1.4) with this definition gives the equality?? = 2~?2me?? = 1.786. (1.5)Solving for distance as a function of barrier height we find?(?) = 1.786~2?2me??0.1743?? , (1.6)for ? in nm and ? in eV. Presumably the barrier height is relatively fixedwithin any single device but this relationship gives a range of possible valueswithin which Wang claims his device operates.An alternative way to experimentally determine values for ? and ? is bymeasuring the frequency response of a tunnelling device. Treating the deviceas a simple Resistor-Capacitor circuit, the dominant resonant frequency canbe used to extract a mathematical description of the typical distance betweenthe particles that are transporting charge. Resistance of the circuit is takento be the one described in equation (1.3). The capacitance between any twoparticles, C?, is described byC? = ??0A? , (1.7)where ?0 is the permittivity of free space and ? is the relative permittivity ofthe insulator. A capacitance equation of this form describes the capabilityof parallel plates to store charge. In the present examination the equation isassumed to approximately hold for semi-spherical plates in order to considerprevious work in this field [3]. Capacitance and resistance can be combinedto describe the resonant frequency,71.2. Tunnelling conductionf0 =12?R?C?, (1.8)of the device in Hz. Substitute (1.3) and (1.7) into (1.8).f0 =3?q2e32?3~??0e???, (1.9)where ? is defined by (1.4). Solve for ? as a function of ?.? =ln( 3?q2e32?3f0~??0)? . (1.10)This analysis method is used by Meier et al. [3] to consider a compositeof polystyrenebutadiene and carbon black and results in the relationship?(?) ? ln? + 43.520.5?? , (1.11)for ? in nm and ? in eV.Whereas [2] fixed the relationship ?? = 1.786 in equation (1.6), ?? ef-fectively ranges from 19.46 to 21.8 in equation (1.11) using data from [3].In Figure 1.4 equation (1.11) is compared to equation (1.6) to gain furtherperspective on the barrier height vs tunnelling gap relationship extractedfrom the work of [2] and [3].The results from [3] suggest typical distances between particles that aresubstantially further apart than the results from [2]. These extremes ofthe data in Figure 1.4 are summarised in Table 1.1. The ? ranges differby a factor of 10. This serves to further illustrate the disagreement amongresearchers regarding typical properties of conductive composites.? (eV) 1 0.01 Data source? (nm) 0.17 1.7 Wang et al. (2010) [2]? (nm) 2.1 19.0 Meier et al. (2007) [3]Table 1.1: Comparison of values for ? and ?.In both sets of results, there must be very small potential barriers (? <0.01 eV) for devices where particles are typically more than 20 nm apart.81.2. Tunnelling conduction0 5 10 15 2000.10.20.30.40.50.60.70.80.91? (nm)?(eV) model by Wang et al. (2010)model by Meier et al. (2007)Figure 1.4: Comparison of potential barrier height (?) and the inter-particlegap over which tunnelling must typically occur (?) using the modelling ap-proach from Wang et al. [2] and Meier et al. [3].91.2. Tunnelling conduction1.2.2 Hard/soft shell approachIn the case of small particles [15] or generally spherical particles with highaspect ratio tendrils [35], it can be instructive to consider the particles ashard solid spheres that have a soft outer shell. The hard solid cores behave asin the typical percolation model. They are incompressible and have a fixedability to conduct when in contact with another particle. The soft outer shelldefines an effective radius of interaction. If two particles under considerationhave overlapping soft shells, they are considered capable of transportingcharge. In a simplified model, this conduction can be considered to behavein an on/off fashion where only spheres with an overlapping soft shell need beconsidered. In this case the critical distance, ?c, is the largest distance (coresurface to core surface) that must be considered in order for a conductingcluster to span the sample 50% of the time [36].Ambrosetti et al. [16] use this model at length and offer the resistanceequationR? = R0e2?c? , (1.12)where R0 is a constant and ? is a characteristic tunnelling length. Comparingequations (1.3) and (1.12) another relationship is apparent.? = 2? , (1.13)or, substituting (1.4),?(?) = ~?2me?. (1.14)This gives two opportunities to interpret physical meaning when ?, ?, or? is used as a fitting constant during the analysis of experimental data.A value for any of them suggests both a characteristic tunnelling lengthand a characteristic potential barrier height. Using (1.14) with the previ-ous employed reference barrier heights yields ?(0.05 eV) ? 0.873 nm and?(1 eV) ? 0.1952 nm. The expression (1.13) results in exactly twice the val-ues shown in Figure 1.3 and implies the relationship ?? = 2. As before thereference barrier heights are not in agreement with expected inter-particledistances.The work done by Heyes et al. employs a series of Monte Carlo andmolecular dynamics simulations to produce a set of equations that describepercolation events when randomly packing spheres using the soft shell andsolid core approach to modelling [4]. Ultimately this results in a high orderequation for critical distance with extremely good fit over the range 0.25 ?101.2. Tunnelling conduction10?7 < ? < 0.62 and normalised to particle solid core diameter (expressedhere as 2r). The critical distance is defined as?c2r = 1 +1f(?) , (1.15)and illustrated in Figure 1.5 usingf(?) = (a0?)13 + (a0?)23 + (a0?) +43 (a0?)43 + 53 (a0?)53+2 (a0?)2 + a1?3 + a2?4 + a3?6,(1.16)where a0 = 2.8902, a1 = 129.644, a2 = ?98.334, and a3 = 1302.40 [4].0.05 0.1 0.15 0.211.522.533.54?normalizedcriticaldistanceFigure 1.5: The normalised critical distance, ?c2r , is shown with respect tovolume loading fraction of solid core spheres. The relationship is describedexplicitly by equations (1.15) and (1.16) taken from [4].The random packing approach to the development of this model is in-structive as the numerical outcomes reinforce the notion that random pack-ing is unlikely to be taking place during the mixing and curing of conduc-tive particles in an insulating polymer matrix. Solving equation (1.15) for111.3. Particle agglomeration and filament formation? = 0.05 and r = 25 nm yields a critical distance of about 90 nm or, us-ing equations (1.12) and (1.14), the tunnelling barrier is on the order of20 ?eV. If the typical radius is increased to 65 nm, the critical distance andtunnelling barrier change to about 240 nm and 4 ?eV respectively. It isnot typically reasonable to expect tunnelling to occur over distances thatlarge nevertheless composites prepared at those loadings and particle sizeshave been shown to exhibit conductivity [8, 37]. An apparent conclusion isthat this simple model makes assumptions that are either not valid or failto account for something in those experiments. The model requires furtherrefinement.1.3 Particle agglomeration and filamentformationAn analysis of carbon black in poly(methyl methacrylate) has revealed thatthe particles formed wire-like structures with wire diameters as small as24 to 31 nm and edge lengths over 1 ?m [38]. This work reveals a distinctlyVoronoi-like structure. Samples were prepared using carbon black with amean diameter of 21 nm. Final preparation employed a press sinter methodthat squeezes the sample during curing. The cured material is then thinlysliced and imaged. A figure of merit in this work is the ratio of filler particlesize to wire diameter. This suggests it is reasonable for wire-like filamentsto form that may be only slightly thicker than the mean diameter of theirconstituent particles [38]. Scanning electron microscopy of carbon black witha mode diameter near 50 nm in polyethylene and polypropylene is in generalagreement, showing filaments and Voronoi cell-like structures with cell edgeson the order of several ?m long and filaments with variable diameters fromless than 0.2 ?m to greater than 1 ?m [39].1.4 Carbon black as a conductive fillerCarbon black is most commonly manufactured by combustion of heavy oils.The carbon precipitates out as small particles. The small particles fuse toform aggregates. Large aggregates tend to have more structure (be lessspherical) with a substantial number of branches. The high surface areapromotes conduction and promotes situations where aggregates, influencedby Van der Waals forces, cluster together into larger groups called agglom-erates [8]. Carbon black can also begin as carbon pellets that are shearedand crushed to form smaller particles. During the shearing and crushing121.4. Carbon black as a conductive fillerprocess carbon black again aggregates and agglomerates because of interac-tions between particle surfaces [40, 41]. In both cases, this aggregated statetypically results during the manufacturing process. Aggregates in mixturewith a rubbery matrix tend to further coalesce to form larger agglomer-ates. During mixing it is believed that agglomerates often break apart andre-form but aggregates remain substantially unchanged [5]. A scatteringprofile study of carbon black before and after combination with polymersshows that the smallest size of carbon black agglomerates grows during themixing and curing process suggesting that the carbon black has bonded intofilaments bounded by polymer chains [42]. The size and shapes of the newcarbon black units depends on the polymer, the mixing conditions, and thecuring conditions [5, 8, 42?46].It is common to describe the behaviour of carbon black conductive com-posites using a percolation model, a tunnelling model, or a combination ofthe two. The two models are limited in their compatibility and still failto explain the myriad differences in experimental results. The most sub-stantive incompatibilities stem from the requirement that in a tunnellingmodel particles must be within a few nanometres of each other yet mustnot touch while they must form a geometric continuum for percolation tooccur. While there is consensus that the conductivity of a particular carbonblack-polymer composite will vary as the volume loading of carbon black,it has also been noted that changing the type, and therefore the carbonblack particle structure and size, may also yield very different conductivitydespite a constant volume loading. Those same geometric differences alsoaffect tunnelling performance as high aspect ratio geometries can lead tolarge electric fields that distort the tunnelling barrier potential profile andresult in augmented quantum tunnelling. Work done by Verhelst reveals thisvariation in the extreme. Five different carbon blacks (with mean diametersof 30 nm, 10 nm, 29 nm, 42 nm, and 29 nm) are mixed with polystyrene-butadiene. The test samples exhibit very different percolation curves andrequire different volume loadings (10.4%, 19.1%, 28.7%, 43.3%, and 53.1%respectively) to achieve the same resistivity (about 7.4 k?m) [46]. In gen-eral qualitative terms this means that particle structure and size must becaptured in any analytical model either directly or, more often, indirectlywithin fitting constants such as a unique percolation threshold, percolationcritical exponent, and potential barrier height for each combination of a typeof carbon black, a polymer, and possibly also each set of mixing conditions.There is agreement among scientists that, for a fixed volume loading,the presence of smaller agglomerates enables a more favourable distributionof particles with smaller inter-particle distances [1, 8, 45]. This increases131.4. Carbon black as a conductive fillerthe likelihood of charge transfer through the insulating barriers betweenparticles.Absolute minimal inter-particle distances are suggested by some of thechemical properties of carbon. The length of a carbon to carbon bondis 0.14 nm and the inter-planar spacing of graphene sheets in graphite is0.34 nm [47]. It has been shown that as loading increases, the distancebetween particle cluster centres decreases but the distance between particlecluster edges approaches a threshold value near 3 nm [5]. This result isshown in Figure 1.6 [5]. An implication of this and the acknowledgementthat agglomerates both break apart and form during the mixing process isthat aggregates closer than 3 nm tend to fuse into agglomerates and thatfreshly broken agglomerates tend to separate by at about 3 nm.0 0.1 0.2 0.30246810?distancebetweenneighbouringagglomerates(nm)Figure 1.6: Typical distance between nearest-neighbour carbon black ag-glomerates in natural rubber as a function of volume loading fraction.Adapted from [5].One of the most popular forms of carbon black used for creating con-ductive composites is Vulcan XC72 made by Cabot Inc [48]. A scanningelectron microscope micrograph of Vulcan XC72 is shown in Figure 1.7.141.4. Carbon black as a conductive fillerFigure 1.7: Scanning electron microscope micrograph of Vulcan XC72.Reprinted from [6] with permission from Elsevier.151.4. Carbon black as a conductive fillerFigure 1.8 contains data points from experimental works done across arange of loadings of two different carbon blacks in four different polymers.These results will serve as a reference for results produced during simulation.The second carbon black, Black Pearls 2000, is also produced by Cabot Inc.Size profiles of Vulcan XC72 done by the author, a colleague [9], andothers suggests a log-normal distribution of particle diameters [37, 46, 49,50]. Spectroscopy measurements and scanning electron microscope mea-surements of particles tend to differ in conclusion regarding the most abun-dant unit size. Spectroscopy measurements tend to give a diameter rang-ing from 100 nm (in water with a surfactant) [9] to 400 nm (in oil) [50]whereas scanning electron microscope analysis of pure carbon black sug-gests a diameter range from 30 nm to 125 nm with a mode diameter near55 nm [11, 37, 46, 49]. A spectroscopic analysis of the material shown inFigure 1.7 (and published with that image) suggests a mode diameter inexcess of 300 nm in clear contradiction with the image in Figure 1.7 [6]. Liuet al. examined transmission electron microscope images of carbon blackthat had been applied in suspension to a substrate and then allowed to dry.The resultant size profile is shown in Figure 1.9. This supports the notionthat particles in suspension are forming agglomerates much larger than thesize of any single particle or aggregate [13].Despite this potential confusion there is consensus in the assertion thatthe particles have a high degree of structure and are prone to aggregation andagglomeration so it may be reasonable to attribute some of this differencein size measurement results to agglomeration within the dispersal mediumused for spectroscopy. Studies of Vulcan XC72 cluster size done after mixingand curing support claims that within a rubbery matrix aggregates on theorder of 100 nm to 150 nm are typical and they readily join to form largeragglomerates [5, 39].For the purpose of modelling, Vulcan XC72 at a loading of 3% by volumeis used as a starting reference as it is a volume loading that shows measur-able conductivity in published works and the particle size distributions forVulcan XC72 are reasonably well described. A typical diameter near 125nm is assumed for small aggregates: some aggregation is assumed to havetaken place with the expectation that the algorithms used to model particlearrangement will capture agglomeration phenomenon.161.4. Carbon black as a conductive filler0 0.05 0.1 0.15 0.2 0.25 0.310?810?610?410?2100102??(S/m) Vulcan XC72 in polydimethylsiloxane (Niu et al.)Vulcan XC72 in polydimethylsiloxane (Rwei)Vulcan XC72 in polydimethylsiloxane (Rwei)Vulcan XC72 in polydimethylsiloxane (Linklater)Vulcan XC72 in polydimethylsiloxane (Li)Blackpearls 2000 in polypropylene (Huang)Blackpearls 2000 in ethylene?octene (Huang)Blackpearls 2000 in ethylene?vinyl acetate (Huang et al.)Figure 1.8: Conductivity across a range of loadings of two different carbonblacks in four different polymers. Taken from works published or shared byNiu et al. [7], Rwei et al. [8], Linklater [9], Li [10], Huang [11], and Huanget al. [12].171.5. SummaryFigure 1.9: Size distribution of Vulcan XC72. Data represents about660 particles in about 200 images. Figure created by [13]. Reprinted bypermission of The Electrochemical Society.1.5 SummaryPrior work has shown that a combination of percolation and quantum tun-nelling effects must be considered in order to successfully model the internalbehaviour of conductive composites. While the physics of individual com-ponents seem reasonably well understood in isolation there remains someambiguity with respect to the combined effects as well as their interactioneffects. Carbon black has been demonstrated as a useful conductive filler,but characterisations of the filler are inconsistent.In general, a successful model must allow for variability in particle size,incorporate agglomeration effects and filament formation, and account forconduction due to both contact and tunnelling. The following sections out-line the modelling approach used and review the mathematics and physicsthat drive the model. Modelling results are then presented and discussed.Finally some thoughts are presented for future work that could serve tomove this model forward or otherwise improve the work described herein.18Chapter 2Experimental methodsModelling experiments were conducted using software written in MATLABand executed on a dedicated high performance computing cluster. The clus-ter consists of 20 Dell PowerEdge 1950 servers each with a pair of 2.66 GHzquad core Xeon processors and between 8 and 32 GB of RAM.In order to approximate the random nature of physical mixing and varia-tion in particle size, several randomization algorithms were created with pa-rameters tuned semi-empirically. A pass in the simulator typically involvedsix steps: load or generate all needed model parameters and establish modelboundaries, place particles and electrodes according to a set of predefinedrules, create a map describing the connections to and from all particleswithin one radius of interaction of each other, calculate the inter-particleconductivity for each connection, calculate the bulk equivalent conductivityof the model, and finally store or present the results in data files or figures.This is diagrammed in Figure 2.1.The first and final steps are reasonably simple and computationally triv-ial. Particle placement methods are described in Section 2.1. Connectionsbetween particles are mapped by examining the space one radius of inter-action from the surface of each particle. Using the language introduced inSection 1.2.2, it is a check to see if the soft shell of one particle reaches outas far as the hard shell of any other particles. The computation of inter-particle conductivity is relatively simple though the underlying physics ofthe calculation warrants analysis. This analysis is done in Section 3. Thebulk conductivity calculation treats the system as a large resistor array andsolves the system using nodal analysis. Computing the solution involvesthe inversion of a large matrix. If the resistor network does not complete apath from end-to-end or has only a weak end-to-end connection a conduc-tivity of or near zero would be expected. This situation creates a singularitywhose computation produces numerical values that are very small and maybe positive or negative. Negative results are easy to identify as being errantbut small positive values may be the result of either no connection or avery weak connection. The lower numerical limit of reliable simulation datapresented herein is somewhere between 10?8 and 10?9 S/m. Alternately19Chapter 2. Experimental methodsFigure 2.1: Processing the present model involves six steps that are describedin the text.202.1. Particle placementput, an estimated numerical error of ?10?8 S/m should be applied to allcomputational results presented herein.2.1 Particle placementParticles are placed in a semi-random fashion throughout a specified vol-ume. Several different distribution algorithms were explored. Several of themethods are described below. All the methods approximate particle shapeas a sphere of fixed radius. In some cases the spheres are of uniform sizeand in others cases there is a distribution of sizes.2.1.1 Method 1: random placement of uniform particlesThis is the most basic placement algorithm and is the foundation upon whichthe others are built. The number of particles needed, n, is calculated asn = V ?43?r3, (2.1)where r is the particle radius and V is the total volume being simulated.Particles are seeded by randomly assigning particle centre coordinates. Avalidity check then occurs to ensure a newly placed particle does not intersectan existing particle. Intersecting particles are deemed to be invalidly placedand are then removed. The process loops until all particles have been placed.2.1.2 Method 2: clusteringParticles are generated as in placement Method 1 but the validity check hasan extra component. In addition to avoiding intersections particles must,with a user defined probability, be within a user defined distance of any ex-isting particle. This allows the user to approximate empirical observationsof particle agglomeration [44] and tendency to form filaments [51]. In eachsuccessive modelling pass the algorithm is repeated for any particles thatfailed intersection related validity checks. This can result in slightly morerandomly placed particles than defined by the user but allows a model whoseother parameters are overly constrained to run to completion with the userthen being informed of any resultant deviations from input parameters. Aweakness of this approach is that computation time quickly becomes pro-hibitively long for small particles in large systems. This weakness was oneof the motivations to explore other approaches that could produce similar212.1. Particle placementor better results. Further insight into this modelling methodology as well assome of the results it produced are presented in Section 4.4.2.1.3 Method 3: volume exclusion methodBefore seeding particles, spherical exclusion zones are created based on auser defined exclusion fraction and exclusion size. These two tunable pa-rameters allow the user to approximate an empirical phenomenon wherebyparticles tend to create shell-like structures around areas that contain fewor no particles [39]. This is generally attributed to surface effects [52, 53].This methodology has been used by others to model the formation ofcarbon black filaments in rubbers [54]. The model created by Jean et al.was formulated and tested using imaged slices of composite material wherethe particles have a mean radius of 20 nm and the slices have a thicknessof 40 nm (? 10 nm). During the extensive analysis of particle aggregationand agglomeration, conclusions are drawn regarding the structure of theparticle distribution. It is noted that a Voronoi tessellation is well suitedto this type of modelling although a volume exclusion methodology is thenpursued to approximate the structure. The authors show that while thecarbon black particles are spherical and take on a log-normal distribution ofsizes, they form aggregates that can be considered spheroids and are aboutten times as large as their constituent particles. The authors use thesespheroidal approximations of aggregates to create a suitable set of sphericalexclusion zones within which no particles may exist. Spherical particles arethen randomly placed within all space except that defined by the exclusionzones. Parameters such as the size and density of these zones are tuned byan automated process that compares slices of the resultant 3D space withimages of the material slices until the particle and aggregate shapes anddensities match to within a pre-set confidence interval [54]. Despite thiseffort, the presented model still fails in two notable ways. It does not offermeaningful insight into the volume resistivity of composites modelled andthe model does not survive scaling (it percolates less as the model scalesup). The former may be an oversight on the part of the authors as they dopresent resistivity data for their reference material but the latter is clearlyan area for improvement as the usefulness of their model is surely tied to itsability to predict the behaviour of increasingly larger volumes of material.Within the present work a volume exclusion modelling method was im-plemented as per the following. Uniform particles were generated as inplacement Method 1, however, particles were allowed to fall within an ex-cluded zone with a user defined probability. The probability of being allowed222.1. Particle placementwithin an exclusion zone or not was included as a user tunable variable inorder to more closely approximate empirical results. This does not affect vol-ume loading fractions but increases the likelihood of particles being nearerto one another. A two dimensional representation of this approach is shownin Figure 2.2. This method yielded results in some cases that were similarto the shell and filament effects described in published works [5, 39, 44, 51?53], however the computation time for the simulations became prohibitivelylong.0 0.2 0.4 0.6 0.8 100.20.40.60.81Figure 2.2: A two dimensional arbitrary unit representation of the volumeexclusion method described in Section 2.1.3. Exclusion zones are shown ascircles with a dashed line (radius = 0.05). Particles are shown as a circle witha solid line (radius = 0.01). In this example, exclusions zones are permittedto overlap, particles are not, and no particles may be in an exclusion zone.2.1.4 Method 4: Voronoi wire frame methodThe wire frame method seeks to improve on the filamentation results ofplacement Method 3 while simultaneously reducing the computation timeinvolved in creating and analysing the model. Rather than creating exclusion232.1. Particle placementzones as a sphere using a centre point and radius, which must then becompared to every particle placed to check for intersection, exclusion zonesare only defined as a point. The array of exclusion points is used to generatea three dimensional Voronoi structure wire-frame in the modelling spacewhere the particles are placed on the ?wires?.A Voronoi diagram divides space using boundaries drawn between seedpoints such that only one seed point is enclosed within each region. All spacewithin each region is closer to its own enclosed seed point than it is to anyother seed point. In two dimensions this is relatively easy to visualize and isshown in Figure 2.3. This style of geometry is consistent with observationsof natural phenomenon such as polycrystalline micro-structures in metallicalloys [55?57]. It is also used in a variety of mapping applications [57, 58].In three dimensional space the principle is the same but the product isvisually more complex. Rather than a map of convex polygons that encloseseed points in a two dimensional plane, the map is instead divided intoconvex polyhedra that enclose seed points in three dimensional space [57].To aid in visualization an example 3D Voronoi structure is presented inFigure 2.4 as a stereographic image set.0 0.2 0.4 0.6 0.8 100.10.20.30.40.50.60.70.80.91Figure 2.3: A 2D Voronoi structure with cell edges shown as lines and dotsat each seed point.Analyses of composite materials have shown Voronoi-like structures ofcarbon black throughout [5, 38, 51, 54]. Work by Jean et al. using carbonblack in natural rubber revealed a characteristic Voronoi edge length of96 nm to 128 nm, 5 - 6 times the size of the constituent particles (whosemean radius is given as 20 nm) [54]. Work by Levine et al. using carbon242.1. Particle placement23456782345678912345678123456782345678912345678Figure 2.4: Stereographic example of a 3D Voronoi structure with cell edgesshown as lines and dots at each vertex. Seed points are shown as a star.Hold image about 50 cm away and allow eyes to cross slightly. When threeimages appear, observe the central one.black in poly(methyl methacrylate) measured unbroken wire segments upto 100 ?m long and mean wire thickness less than ten times the size of theconstituent particles (whose mean radius is given as 21 nm) [38]. Thesenumbers serve as references when adjusting model parameters.In the present model, seed points are randomly assigned. The volumeand perimeter of each Voronoi cell is calculated. During particle placement,each particle is randomly assigned to either be placed anywhere in the modelvolume, or be placed within a user defined range of one of the Voronoi cellsegments. Assignment to a particular segment is weighted to ensure allsegment sections have equal probability of being selected for a particle, inother words, an equal probability per unit length. Finally, particles arechecked for intersection.2.1.5 Method 5: distributed radiiParticles were sized according to a log-normal distribution centred aboutpublished mode particle size values as described in the literature [6, 50, 54].The mode size of the particles is preferred to the mean as the mode sizeidentifies the most commonly occurring (or the typical) particle size but thetwo values differ for a log-normal distribution. During discussions of inputparameters and output results a combination of mode and mean must how-252.1. Particle placementever be used. This is because while ?mode? is more useful when definingthe type of distribution the randomization algorithm is to produce, it is notpossible to directly calculate the mode of a discrete number of points eachof which is drawn from a continuous distribution. In cases where a pub-lished size range was used as a reference the range values were treated asthe 95% confidence intervals for the log-normal distribution. This combinedmode and variance were set at execution time but they were permitted todrift. The drift would occur as a result of particles being rejected due tointersection. During each subsequent placement pass, new particles weregenerated using the same probability distribution function. Because therewas a greater likelihood of larger particles experiencing intersections, thelikelihood of having an abundance of smaller particles increased with eachpass. Placement passes continued until the total required volume of parti-cles was reached. An example of the output from this algorithm is shownin Figure 2.5. Using a log-normal distribution allowed particles to take on arandom radius while maintaining non-negative values and favouring smallerparticles over larger ones. This was deemed appropriate as mechanical mix-ing has been shown to break apart large agglomerates [5]. All simulationsresults presented herein include algorithms that create particles whose sizesare distributed in this way.2.1.6 Method 6: intersections allowedThe placement methods described above all served to better approximatethe formation of agglomerates, clusters, and filaments of particles describedin literature [5, 44, 53]. They also improved the bulk connectivity of themodel. The connectivity improved because the mean distance between par-ticles decreased and particle networks became more likely to span the extentof the model. This accounts for connectivity between particles where a tun-nelling or hopping type conduction mechanism is dominant but does notreadily account for connectivity that might result if two particles were tobe in physical contact. The placement algorithm, despite all the additionalconstraints introduced, still assigns a bounded random coordinate to eachparticle, a bounded random radius, and then ensures it does not intersectother particles. This results in a situation where the likelihood of two par-ticles being in contact is approximately zero and, even if it were to occur,it would be a point-contact with no cross-sectional area and thus be unsuit-able for conventional conductivity consideration. To approximate situationswhere particles may come into contact and stay in contact with each other,slight intersection of particles was allowed and the amount permissible was262.1. Particle placement20 40 60 80 100 120 14005001000150020002500particle radius (nm)numberofparticlesFigure 2.5: In this example set of particles the user set target mean radiusand standard deviation were 65.0 nm and 9.75 nm respectively. A log-normaldistribution with those parameters is shown as a dashed line for reference.It has a mode radius of about 62.8 nm. The resultant particles from foursuccessive models are combined and histogram data of the resultant radiusdistribution is shown as the solid line. Histogram bins span about 0.32 nm.Particle radii vary from 31.2 nm to 124.5 nm with a mean and standarddeviation of 64.3 nm and 9.55 nm respectively.272.2. Particle connectivity mappingsubject to user adjustment.2.2 Particle connectivity mappingHaving placed all the particles using one of the methods described above,they must now be evaluated to determine their connectivity to each otherparticle being modelled. The degree of effective connectivity between par-ticles is dependent on their proximity to one another. Two connectivitysituations are considered. When particles are in direct contact with one an-other the notion that they are connected is obvious. A case requiring greaterreflection occurs when particles are near one another but do not touch. Aspreviously discussed in Section 1.2.2, one way to check for and model thistype of connectivity is by defining particles as having a soft, permeable shellthat surrounds a solid core [4, 59, 60]. If a soft shell touches a neighbouringsolid core, particles are said to connect. Solid cores are not typically permit-ted to overlap. In a basic mathematical model this type of connectivity hasa binary state: particles either connect or they do not; they either form anetwork or they do not. An alternate interpretation is to use this approachto check for connectivity but define the connection in quantum mechanicalterms as a gap over which tunnelling occurs. The soft shells then define aradius of interaction where tunnelling conduction will be computed. Therule prohibiting hard shells from overlapping is relaxed in order to allowrandom placements that result in direct, classically conductive connectionsto occur within agglomerates.Because tunnelling conduction falls off exponentially as tunnelling dis-tance increases, there will be a threshold particle gap or soft shell size abovewhich the amount of charge transport is negligible and particle connectionscan be said to be absent. Computationally it is practical to assign an upperbound on the distance worthy of consideration. Beyond that upper bound,called the radius of interaction, ri, the amount of charge being transportedis taken to be sufficiently negligible that it need not be considered.In a model filled with randomly placed particles all potentially contribut-ing to charge tunnelling or hopping it is conceivable that each particle trans-mits charge to or receives charge from every other particle. This is capturedin a map that describes the distance between all particles. Particles thatintersect and are therefore in direct contact will have a negative gap equalto the amount of linear overlap. Non-overlapping particles are considered assites where tunnelling is occurring. As previously discussed, some of thosecharge flows will be negligible. The negligible flows can be identified by282.3. Electrodestheir correspondingly large distance between particles allowing the map tobe pruned to only include the particle connections spanning gaps less thanthe radius of interaction.2.3 ElectrodesAs the ultimate goal of the model is to represent a bulk material and ap-proximate its end-to-end conductivity the model requires that considerationbe given to electrodes. In physical embodiments the electrodes are typicallya good conductor with physical properties that suit the associated experi-ment or application. Examples include conductive grease, evaporated met-als, foil, and poly3,4-ethylenedioxythiophene/polystyrenesulfonate (a con-ductive polymer) [15, 22, 35].In our model, electrodes are represented fundamentally the same as otherconductive points but they have some special properties: they are geomet-rically confined to the planes of the electrodes, they have no thickness, andthey are sufficiently abundant as to approximate a continuum. The first twocriteria are absolute and the third is subjective. When modelling electrodesin this way there are two possible outcomes if the abundance of electrodesites is ill defined. If there are too few electrode sites it is likely that de-spite a possible internal network of particles there will be a lack of electrodeconnection sites within the radius of interaction of the particles in the bulkdespite many particles being near the electrode plane. An example of thisis shown in Figure 2.6A. This is analogous to an overly high contact resis-tance. If there are too many electrode sites, connections from the nearbyinternal particles will be very likely to occur and potentially over representthe connectivity. Figure 2.6B shows the same internal network as Figure2.6A, but with an increased number of electrode sites. This is analogousto an overly low contact resistance. In this model the latter situation isfavourable. An artificially low contact resistance ensures that the resistanceof the bulk dominates.In order to ensure good connectivity without creating an overwhelmingnumber of electrode sites, electrode sites that do not connect with any bulkparticles are removed before performing conduction calculations.2.4 Node trimmingAs the model requires computationally intense mathematics and the inver-sion of large matrices there is incentive to use the fewest number of electrical292.4. Node trimming0 0.2 0.4 0.6 0.8 100.20.40.60.81A0 0.2 0.4 0.6 0.8 100.20.40.60.81BFigure 2.6: In this two dimensional arbitrary unit example conductive parti-cles are shown as ? and electrodes are shown as *. Connections between sites(particles or electrodes) are shown as a line. Diagram A shows a connectednetwork that has few electrode connection sites. Diagram B shows the samenetwork but with more electrode connection sites. Observe the differencebetween the two diagrams in particular near x-axis locations 0.2, 0.5, and0.8.302.4. Node trimmingnodes possible. Each particle within the model will develop its own potentialand participate in charge transport as long as it is either connected to atleast two other particles or connected to one other particle and an electrode.Each electrode is composed of many particles without a radius and all itsconstituent particles are considered to be equipotential. One electrode isgrounded and the other is held at a fixed voltage. These details allow theimplementation of an algorithm for determining if any electrode sites can beignored when tabulating the electrical nodes. Any particle without a con-nection to another particle is not considered an electrical node. Any particlewith only one connection to another particle must also have a connectionto an electrode or else it is not considered an electrical node. These par-ticles that are marked as not being electrical nodes will be ignored duringconductivity calculations. An example of this is shown in Figure 2.7.0 0.2 0.4 0.6 0.8 100.20.40.60.81Figure 2.7: This two dimensional arbitrary unit example contains the sameconnected network as Figure 2.6B. Conductive particles are shown as ? andelectrodes are shown as *. Connections between sites (particles or electrodes)are shown as a line. Three particle sites that do not contribute useful elec-trical nodes have been identified (with circles) as suitable for exclusion fromthe conductivity calculation.31Chapter 3Inter-particle conductivityOnce the mapping of connected particles is complete and trimmed each con-nected pair is considered individually to determine the conductance betweenthem. Particles may be in contact or have some distance between them. Thepresent work developed a variety of methods used to model the charge trans-port occurring in conductive composites. The first method presented is verybasic and the subsequent methods are more sophisticated.3.1 Absolute path lengthIn the simplest of modelling situations a conductance, G, is assigned pro-portional to the linear distance between the particles.G = ??? , (3.1)where ?? is a conductivity parameter that describes the conductivity of afixed cross section, and ? is the distance between the two particle surfaces.This method accounts solely for the empirical observation that as the meandistance between particles decreases, conductance increases. It does notconsider any of the more complicated physics that contribute to conductivity.The computation time is very short and this method can be instructivewhen looking for general changes in the model due to alterations in inputparameters such as geometry, particle size, and volume loading. It providesvery little physical meaning but does offer insight into whether the numberof end-to-end connections has increased or decreased due to a parameterchange in the simulation.3.2 Basic exponentialUsing an exponential algorithm to calculate inter-particle conductance, Gt,approximates a tunnelling type of charge flow. The equation has the formGt = G0 e??? + G1 (3.2)323.3. Conduction due to particle contactwhere G0 is a base conductance and G1 is the conductance of the insulatingfiller. The physical meaning of ? is explored in the definitions of equa-tions (1.4) and (1.13). This relationship between conductance and distanceis applied to all gaps greater than zero and less than the upper limit imposedby the radius of interaction and captures point-to-point charge tunnellingacross the distance between particle edges without considering the cross sec-tion of the particles. Using this information G0 is defined using (1.3) as areference.G0 =3q2e16?2~ . (3.3)This gives G0 ? 4.6 ?S.Models that describe imaging systems such as scanning tunnelling mi-croscopes successfully use this analytical approach for scanning tips with aradius of curvature less than 200 nm suggesting that this is a reasonableapproach for models with particles sizes less than that limit [33].3.3 Conduction due to particle contactConduction arising from particle contact is assumed to follow classical con-duction theory. The amount of contact occurring is related to the amountthe particles overlap in the model. The conductance in this area, Gc, isapproximated asGc = ?cAca , (3.4)where ?c is the internal conductivity of the conductive particles, Ac is thecross-section of the contact area, and a is the distance between sphere-centres.Ac = ?r2c , (3.5)where rc is the radius of the contact circle made by the two overlappingspheres. From Figure 3.1 the distance between sphere-centres isa = r1 + r2 + ?, (3.6)where r1 and r2 are the respective radii for the two spheres being considered.The angle between the point of overlap and the straight line betweencentres is shown as angle ?. The law of cosines is used to describe therelationship between these lines and ?,r22 = r21 + a2 ? 2 r1 a cos(?). (3.7)333.4. Conductivity summaryr1 r2rcr1r2-?a?Figure 3.1: Geometry between overlapping spheres, ? < 0. When ? > 0,rc = 0, and the spheres do not overlap.This can be rearranged to find? = cos?1(r21 ? r22 + a22r1a). (3.8)The contact radius,rc = r1 sin(?), (3.9)is found as a function of ? through trigonometry. Substituting (3.6) and (3.8)into (3.9) givesrc = r1 sin(cos?1(r21 ? r22 + (r1 + r2 + ?)22r1(r1 + r2 + ?))), (3.10)for ? < 0.Substituting (3.5) and (3.6) into (3.4) producesGc =?c?r2cr1 + r2 + ?. (3.11)Equation (3.11) is computed numerically using rc as defined by (3.10)above.3.4 Conductivity summaryFor a model that considers conduction over a distance, as well as contactconduction we can combine equations (3.2) and (3.11) to give a piecewise343.4. Conductivity summarydefinition for G that depends on ?.G(?) ={?cpir2cr1+r2+? if ? < 0G0 e??? + G1 if ? ? 0.(3.12)Figure 3.2 shows a combined example of the results from using (3.2) tofind tunnelling conduction and solving (3.11) to find the contact conduction.This plot was generated from a set of simulations with a target mean particleradius of 65 nm and target standard deviation of 9.75 nm (actual radiusdistribution shown in Figure 2.5). As the gap increases beyond 15 nm (notshown) the conductance continues to decay in a log-linear fashion. As thegap approaches 0 nm, the tunnelling conduction begins to increase sharply.Once contact has been made, classical charge transport dominates quickly.A minimum contact conduction is set at G0. In this example, this minimumis surpassed when the overlap exceeds 0.01 nm resulting in a contact radiusapproaching 2.5 nm. While the relationship between gap and conductivityis relatively smooth, a small amount of noise is visible and is due to therandom matching of particles of differing radius.353.4. Conductivity summary?15 ?10 ?5 0 5 10 1510?710?610?510?410?310?2? (nm)conductance(S) tunnellingcontactFigure 3.2: Conductance between particles arranged according to increasinggap size. Particle size varies as illustrated in Figure 2.5 and the effect ofthis variation can be seen as a small amount of noise in the line. A negativegap indicates particles overlap and the particles are therefore in physicalcontact.36Chapter 4Results4.1 Model sizeModelling an abundance of particles with a radius less than 100 nm in a cubicspace whose edges are on the order of millimetres or larger is not computa-tionally trivial. Models were computed on a cluster with nodes containing8 cores and memory availability ranging from 8 to 32 GB as described inSection 2. Using this equipment to run simulations, a cubic space with edgelengths on the order of tens of microns may be simulated. It is reasonableto consider these sizes are representative of bulk conductivity provided thatan edge length is chosen such that variations in edge length do not causesubstantial changes in simulation results. In order to determine where thatthreshold lays, a series of simulations were performed with increasing edgelengths of cubic space. The time and memory required to run a simulationincreases substantially with an increase in cube dimensions. For instance,to simulate a cube with an edge length of 18 ?m takes about 24 hours and acube with an edge length of 20 ?m takes about 48 hours. It was found thatthe conductivity did not change substantially for edge lengths above 18 ?mwhen simulating particles with a mean radius less than 100 nm at a load-ing of 3% by volume using the most complex modelling scenarios discussedbelow. This edge length value is used for all further simulations describedin this work as a compromise between asymptotic saturation and modelcomputation time. These results are shown in Figure 4.1.4.2 Random placement of uniform particlesUniform particles are placed randomly in 3D space as per the method de-scribed in Section 2.1.1. Particles are checked to ensure they do not overlap.Initial tests were done at 3% volume loading with particles of diameter125 nm. The basic exponential (3.2) is used to compute conduction usingG0 as defined by (3.3).In order to promote conduction in the model a rather low barrier height is374.2. Random placement of uniform particles10 12 14 16 18 200.010.020.030.040.050.06?(S/m)cube edge length (?m)08001600240032004000time(minutes)Figure 4.1: Conductivity is shown to increase asymptotically to a maximumvalue as the edges length of a cube being simulated increases. Individualresults are shown as a circle. The long-dash line indicates mean value andthe dotted lines are one standard deviation above and below the mean. Thesolid line shows the exponential increase in computation time as the edgesize increases.384.3. Contact conductioninitially used, ? = 0.1 meV. This is substantially lower than barrier heightsthat are typical in these materials but it serves as a useful starting place asthe model is developed. Using (1.4), this gives ??1 = 9.8 nm. The radiusof interaction is swept over a range of values to determine typical distancesat which key inter-particle connections are being made. Figure 4.2 showshistogram data of connection lengths along with the associated numericalconductance per connection of that length. Figure 4.3 shows the net con-ductivity of the model as the radius of interaction is changed. To create thisdata set five simulations were run through all the steps required to place par-ticles. Then the conductivity of each model was evaluated where the radiusof interaction was swept from 100 nm to 300 nm. This effectively permitsthe simulation to consider, at first, only nearby neighbours, and, in the end,more distant neighbours. In each case, if the network of conductive parti-cles does not create a complete path through the model then the resultantconductivity is that of the insulating matrix and is not shown. At a radiusof interaction greater than 200 nm, all five models exhibit conductivity near0.22 ?S/m. When only considering particles less than 180 nm apart, onlythree data points are generated; two of the models failed to exhibit end-to-end connectivity of conductive particles. With the radius of interactionbelow 150 nm only one of the five models continued to exhibit connectivity.It can be concluded therefore that in order to repeatedly achieve a modelwith end-to-end connectivity using this particle placement method, the sim-ulation must allow tunnelling events to occur over distances greater than200 nm. This minimum value for the radius of interaction already uses overlyoptimistic potential barrier values suggesting that the placement methodused fails to result in an appropriate layout of conductive particles.4.3 Contact conductionAs described previously, the initial model calculates the tunnelling conduc-tion between particles but as the particle generation algorithm does notpermit particle edges to overlap there will not ever be a case where contactconduction occurs. By allowing a small overlap during the particle place-ment phase of the simulation (as described in Section 1.4), situations willarise that are analogous to the formation of agglomerates during a mixingprocess.Figure 4.4 illustrates the jump in conductivity that occurs at contact.For gaps greater than zero this result is identical to that shown in Figure4.2.394.3. Contact conduction0 50 100 150 200 250 30010?1910?1610?1310?1010?7? (nm)tunnelling conductance (S) 0123x 104number of inter?particle gapstunnelling conductancenumber of gapsFigure 4.2: Particles that are deemed ?connected? will conduct charge ac-cording to equation (3.2). The amount of conduction that occurs dependson the gap between particles. The right axis (and thin line) show histogramcounts for gaps of a particular length from a typical execution of the model.Histogram bins span about 5.5 nm. For reference, the left axis (and thickdashed line) show the results of the tunnelling equation (3.2) for a gap ofthat size.404.3. Contact conduction100 150 200 250 30010?1010?910?810?710?6ri (nm)?(S/m)Figure 4.3: Five unique models were created with the same initial param-eters (described in the text). Once each model was created, the radius ofinteraction was varied to create maps of connected particles. If a modelexhibited end-to-end connectivity of its particles at a particular radius ofinteraction, a conductivity value for that combination is shown. Data to theleft of the vertical dashed line is intermittent and data to the right of thevertical dashed line begins to converge.414.3. Contact conduction0 50 100 150 200 250 30010?1910?1710?1510?1310?1110?910?710?510?3? (nm)tunnelling conductance (S) 0123x 104number of inter?particle gapstunnelling conductancenumber of gapsFigure 4.4: Particles with a negative inter-particle distance are said to be incontact and conduct according to equation (3.11). Particles not in physicalcontact that are deemed ?connected? will conduct charge according to equa-tion (3.2). The right axis (and thin line) show histogram counts for gaps ofa particular length from a typical execution of the model. Histogram binsspan about 5.8 nm. The left axis (and thick dashed line) show the results ofthe contact equation (3.11) and tunnelling equation (3.2) for a gap of thatsize.424.4. Clustering and filamentsAs per equation (3.2), conductance between two particles is exponen-tially dependent on the height of the barrier between them. In Figure 4.5results are shown for three different values of ? in order to see the effecta changing barrier height has on conductivity. Barrier heights 0.1 meV,1 ?eV, and 10 neV are equivalent to values of ??1 of 9.76 nm, 97.6 nm,and 976 nm respectively. While the units of these values and their discus-sion in Section 1.2 suggest a possible connection to physical properties, thevalues do not make sense when compared to reasonable potential barriersand reasonable tunnelling distances. To the right of the dashed vertical line,conductivities are seen to converge on a conductivity near 0.1935 ?S/m,0.887 S/m, and 4.45 S/m respectively. It is noteworthy that although thenet conductivity rises dramatically with an increasing ??1 due to its placein the exponent of equation (3.2), the inter-particle distance at which themodel begins to conduct remains relatively constant.This set of results illustrates that although a model seeking to representsome of the physical realities of agglomerate formation must allow particlesto physically connect in some way, the simple introduction of an algorithmto this effect has not solved problems previously described with respectto particle proximity and unrealistic tunnelling distances. The varying of ?showed an increase in the converged conductivity but did not alter the radiusof interaction at which convergence began. The simulation still requires aradius of interaction above 200 nm to reliably produce results.4.4 Clustering and filamentsIn order to approximate the natural formation of filaments during the mixingand curing process an algorithm is used whereby, with a specified probability,particles must be within a specified range of a particle that has already beenplaced (as described in Section 2.1.3). This results in localised variabilityin volume loading although the bulk volume loading remains constant. Asexpected, specifying relatively large distances in the enforcement algorithmhas no effect. Specifying progressively smaller enforced distances increasesthe likelihood of particles being within a short radius of interaction of eachother. Specifying relatively small distances in the enforcement algorithmhas the effect of many particles being within range of each other but has theside-effect of creating particle clusters that do not always connect to otherclusters. To aid in the visualization of the problems described, two dimen-sional examples of them are shown in Figure 4.6 and Figure 4.7. Figure 4.6shows a set of particles where the enforced proximity value is relatively large.434.4. Clustering and filaments50 100 150 20010?1010?5100ri (nm)?(S/m)Figure 4.5: Five unique models were created with nearly all the same initialparameters to generate each set of results (each set is shown with a differentsymbol). The only difference between each set of model parameters is thenumerical value for ? (indirectly defined by setting the barrier height to threedifferent values). For ?, ? = 0.1 meV. For ?, ? = 1 ?eV. For ?, ? = 10 neV.Once each model was created, the radius of interaction was varied to createmaps of connected particles. If a model exhibited end-to-end connectivityof its particles at a particular radius of interaction, a conductivity valuefor that combination is shown. Data to the left of the vertical dashed lineis noisy and intermittent and data to the right of the vertical dashed linebegins to converge.444.4. Clustering and filamentsIt is functionally equivalent to a completely random placement and no im-provements are noted. Some particles connect but rarely are more than fouror five particles connected in succession. Figure 4.7 shows a set of the samenumber of particles where a relatively short enforced proximity has resultedin localised accumulations of particles. Where connections occur there isoften a large number of particles interconnected but the structures might bebetter described as globs than filaments. Each of these challenges is morepronounced when considering the problem in three dimensions.Although the enforced proximity did not immediately increase end-to-end connectivity further study was performed by sweeping through a varietyof enforcement probabilities and enforcement proximities. The radius of in-teraction required for bulk conductivity to occur is actually larger than with-out the enforced proximity (an increase to about 200 nm from the 150 nmshown in Figure 4.5). The effect of this on the bulk conductivity seems ap-parent in Figure 4.8 where an enforced proximity of 25 nm produces a pro-gressively lower mean conductivity as enforcement probability is increasedat the upper radii of interaction. While this trend does seem to appear,overall conduction is weak and near the lower computational limits of thesimulator (as described in Section 2) limiting the reliability of conclusionsdrawn from this data set. An optimum enforced maximum range value wasfound near 100 nm. This maximum range value resulted in increased for-mation of filaments that span the entire bulk. These results are shown inFigure 4.9. In both figures the probability of proximity enforcement is sweptfrom 0.50 to 0.95. The effect of increased enforcement probability is shownto amplify the effects discussed above; for short range enforcement (25 nm)the bulk conductivity drops and for medium range enforcement (100 nm)the conductivity increases.Ultimately this methodology resulted in zones of relatively higher par-ticle density rather than forming filaments. When the higher density zonesoverlapped increased conductivity was observed and large portions of themodel space were left without any meaningful number of conductive par-ticles. The intention of this modelling approach was to more closely ap-proximate experimentally observed particle distributions, including filamentformation, and, as a result, exhibit better conductivity with a shorter ra-dius of interaction. Although a set of parameters was found that improvedperformance at a shorter radius of interaction, the placement algorithm didnot generate filament-like formations. In addition, it is difficult to reason-ably defend the choice of modelling parameters used to create the optimizedresult.454.4. Clustering and filaments0 0.2 0.4 0.6 0.8 100.20.40.60.81Figure 4.6: Three hundred particles are distributed in two dimensional space(arbitrary units). Seventy five of them are placed randomly and the re-mainder are placed in succession where each must be within 0.5 units ofa previously placed particle. Each particle is shown as ? and its radius ofinteraction is shown as ? (ri ? 0.015 units). For the purposes of this illus-tration, consider particles connected if their ??s touch. The result is similarto a random distribution. Some connections exist but rarely involving manyparticles.464.4. Clustering and filaments0 0.2 0.4 0.6 0.8 100.20.40.60.81Figure 4.7: Three hundred particles are distributed in two dimensional space(arbitrary units). Seventy five of them are placed randomly and the re-mainder are placed in succession where each must be within 0.05 units ofa previously placed particle. Each particle is shown as ? and its radius ofinteraction is shown as ? (ri ? 0.015 units). For the purposes of this illus-tration, consider particles connected if their ??s touch. Several groupings ofmany particles are apparent but the groupings do not often interconnect.474.4. Clustering and filaments0 50 100 150 200 25010?1010?910?810?7ri (nm)?(S/m) proximity probability = 0.50proximity probability = 0.75proximity probability = 0.85proximity probability = 0.95Figure 4.8: Conductivity values as a function of radius of interaction whenan enforced proximity of 25 nm is implemented. Each data point representsthe mean result of five simulations and each line is associated with a differ-ent enforcement probability as shown in the legend. Overall volume loadingof particles is held constant at ? = 0.03. As described in Section 2, conduc-tivities near 10?9 S/m are of unreliable precision due to a near-singularityduring matrix inversion. This accounts for the erratic low mean conductivityvalues.484.4. Clustering and filaments0 50 100 150 200 25010?1010?910?810?710?610?5ri (nm)?(S/m) proximity probability = 0.50proximity probability = 0.75proximity probability = 0.85proximity probability = 0.95Figure 4.9: Conductivity values as a function of radius of interaction whenan enforced proximity of 100 nm is implemented. Each data point representsthe mean result of five simulations and each line is associated with a differentenforcement probability as shown in the legend. Overall volume loading ofparticles is held constant at ? = 0.03.494.5. Voronoi cell network4.5 Voronoi cell networkImplementing a node networking scheme that more closely follows a Voronoinetwork as described in Section 2.1.4 resulted in a substantial jump in modelconductivity. Key model parameters are the quantity of Voronoi cells oc-curring per unit volume, vdens, and the maximum thickness of the filamentsformed, vrMax. Values for vdens were determined empirically by runningsuccessive simulations and comparing the resultant cell edge lengths to typ-ical values as described in Section 2.1.4. Cell edge lengths were targetedto be in the range from several nm?s to 10?s ?m. This was found to oc-cur in models with 106 mm?3 < vdens < 107 mm?3. Figure 4.10 showsan example histogram of edge lengths for a model with a target density ofvdens = 5? 106 mm?3.0 2 4 6 8 10 120102030405060edge length (?m)mean number of edgesFigure 4.10: Histogram data of typical mean edge lengths for a model witha target density of vdens = 5? 106 mm?3. Each histogram bin spans about0.95 ?m.Filament maximum thickness values were swept from 10 to 30 timesthe mean particle radius. In this case that equates to approximately 1 ?m504.5. Voronoi cell network< vrMax < 2 ?m. The de-facto minimum thickness is the diameter oneparticle (or possibly no particles as the case may be). This range of modelinputs is consistent with the works described in Section 1.3. The filamentsare not forced to be consistently of the thickness defined by vrMax, ratherthey are permitted to contain bulges and narrows provided none of the bulgesexceed the defined maximum thickness. Figure 4.11 shows the results of thisvrMax sweep as a function of volume loading and gives the resultant modelconductivity. Each parameter set shows a conductivity evolution with a sim-ilar trend as loading increases. The thicker filaments likely exhibit slightlylower conductivity because they result in conductive particles spread moreradially along the axis of the filament whereas a narrower filament restrictsthe particles to more confined and regular filaments. For reference, Fig-ure 4.12 shows the results of the vrMax = 1.06 ?m sweep within the contextof the experimental results previously shown and individually attributed inFigure 1.8. Also shown for reference are converged values from results pre-viously discussed and shown in Figures 4.3, 4.5, and 4.9. The Voronoi modelproduces a trend line that runs parallel to those measured in experimentswith carbon black and polydimethylsiloxane. Despite the shape of the re-sults being similar, it is apparent that the magnitude of the Voronoi modeldata (about 10?3 S/m at 3% volume loading) is substantially smaller thanthe magnitude of the experimental data (about 10?1 S/m at 3% volumeloading). The previously used modelling methods all generated results evensmaller (about 10?6 S/m at 3% volume loading). Also included for referenceis an upper limit to carbon conductivity. The limit is computed by consid-ering all carbon matter present at a particular volume loading arranged intoa column that perfectly and uniformly connects through the model space.The Voronoi model data takes on a shape that conforms to the powerlaw behaviour discussed in Section 1.1. Figure 4.13 compares that data toequation(1.2) with ?c = 0.05 and t = 5. A constant of proportionality,105, is applied to overlap the power law shape with the model data. Whenperforming the fit an appropriate value for ?c is not immediately clear. Itrepresents the inflection point of the power law curve. In the case of thesemodel results that point is below the minimum computable conductivity.The parameters chosen give a reasonable fit across the data available. Bal-berg?s assertions suggest this value of t indicates the particles are ?morespherical? though more quantitative implications remain unclear [17].In order to ensure end-to-end model connectivity during model devel-opment, the radius of interaction has been maintained relatively high suchthat conduction decay is dependant on the value of ?. The latest modelswere repeated while evaluating the variation in conductivity as a function514.5. Voronoi cell network0 0.01 0.02 0.03 0.04 0.0510?1510?1010?5100??(S/m) vrMax = 0.75 ?mvrMax = 1.06 ?mvrMax = 1.38 ?mvrMax = 1.69 ?mvrMax = 2.00 ?mFigure 4.11: Conductivity as a function of volume loading for a range offilament radius values.524.5. Voronoi cell network0.01 0.02 0.03 0.04 0.0510?1010?810?610?410?2100102104??(S/m) experimentalrandom placementcontact allowedclusteringvoronoi cellsupper limit of ?Figure 4.12: Model results are shown here with published experimental data[8, 9]. (Individual attributions are shown in Figure 1.8.) Model results arethose previously shown in Figures 4.3, 4.5, 4.9, and 4.11. For reference, theupper limit of conductivity is also shown.534.5. Voronoi cell network0 0.01 0.02 0.03 0.04 0.0510?1010?810?610?410?2100??(S/m) Voronoi model datapower law (?c = 0.05, t = 5)Figure 4.13: Voronoi model results are compared to a line fit for the powerlaw behaviour discussed in Section 1.1.544.5. Voronoi cell networkof radius of interaction. This is illustrated in Figure 4.14 for four differentvolume loadings.50 100 150 20010?1010?5ri (nm)?(S/m) ? = 0.01? = 0.02? = 0.03? = 0.04Figure 4.14: Conductivity as a function of radius of interaction for fourdifferent volume loadings. Data to the lower-left of the dashed line is noisyand intermittent in end-to-end connectivity and data to the upper-right ofthe dashed line is convergent.In contrast to prior plots of conductivity vs radius of interaction, in Fig-ure 4.14 it is clear that end-to-end connectivity occurs when consideringneighbours less than 100 nm apart, an improvement over prior iterations ofthe model. The saturation of conductivity that occurs shortly after connec-tivity first occurs is illustrative of the exponential dependence on distance intunnelling scenarios. This overall improvement is attributed to the filamentsrepresented by the Voronoi network. Although the minimum radius of inter-action needed for each of the simulation sets shown in Figure 4.14 is greaterthan might be expected for a tunnelling scenario, it is an improvement overprior models.55Chapter 5Conclusions and future work5.1 ConclusionsThe model created and evaluated as part of this work synthesizes theorydrawn from first principles, theoretical explanations of observed phenomena,and basic analysis of experimental data. The semi-empirical approach usedcombines many parameters, most of which can be said to have a physicalinterpretation. Many previous works have neglected these physical interpre-tations in favour of good agreement between analytical expressions and thedata being published concurrently with those expressions. This model dif-ferentiates itself by identifying obvious or potential physical meaning behindconstants that could otherwise be treated as simple fitting parameters. Incases where parameters are permitted to move outside apparently reasonablebounds a discussion of assumptions and motivations is provided.A great deal of uncertainty can be found in the parameter ? and, throughthe relationship defined in (1.4), ? in the tunnelling conduction exponential.Although the physics describing the height of the barrier being tunnelledthrough holds in simple situations such as with a scanning tunnelling mi-croscope, it seems apparent that in the case of many conductive compositesthere are more factors that warrant consideration. While not clearly under-stood they either modify or act in addition to barrier height effects. Ratherthan introduce additional parameters, these effects are left bundled in thebarrier height parameter to illustrate this gap in understanding.The initial model shows that a random distribution of spherical particles,given realistic ? and ?, cannot explain experimentally observed conductiv-ities in these composite materials. The final model shows conductivity insimulations loaded less than ? = 5% by volume. The conductivity is severalorders of magnitude lower than experimental results that use carbon parti-cles of similar dimension and upon whose structure the model was based.This discrepancy further illustrates that the model still lacks a comprehen-sive set of physics and input parameters.Ultimately, it shows that bulk conductivity at loadings well below theclassical percolation threshold can be simulated by combining known physics565.2. MATLAB as a modelling toolwith empirical observations of micro-structure in conductive composites cre-ated using particles on the size order of 200 nm and smaller. The trend insimulation results agrees with the trends observed in published experimentaldata. The differences between the simulation results and the experimentalresults serve to motivate future exploration of these materials in pursuit ofa more complete understanding of the phenomena at work.5.2 MATLAB as a modelling toolAll simulations were conducted in MATLAB Version 7.9.0.529 (R2009b).Most algorithms were tested on a standard laptop computer before exe-cution on the high performance computing cluster described in Section 2.Compute time for simulations was a key concern during model development.The core algorithm for serially generating lists of locations in space for par-ticle placement and the earliest version of the user interface were created byLiam Russel in 2010. A serial approach to subsequent modelling algorithmsfollowed naturally from this start. As model complexity grew the computetime for modelling also grew rapidly, in part due to its serial nature. Mostmodel components were later rewritten to take advantage of memory allo-cation and parallel processing functionality in MATLAB. A model that wasdesigned with the intention of performing parallel operations at all possibletimes could serve to improve simulation durations on a platform such as ahigh performance computing cluster.Long compute times (during development they reached over seven daysfor a single simulation), although undesirable, can be overcome with pa-tience. The most substantial obstacle to model computation was the cubicgrowth in memory requirements as the model size increased with one par-ticular matrix inversion making the greatest draw on resources. On severaloccasions memory usage exceeded the capability of the computing clusterand resulted in a failure of the simulation after many hours or days of com-putation.MATLAB provides a good set of tools for developing parallel process-ing algorithms and for memory management. Ultimately the limitationsdescribed above can be attributed somewhat to choice of modelling algo-rithms and predominantly to hardware limitations.575.3. Future work5.3 Future work5.3.1 Greater complexity in modelling parametersOne approach to increasing the ?realism? of this model would be to incor-porate more first principles physics and more empirical observations. Forinstance, incorporating a set of modelling routines that move particles ac-cording to physical processes at work during mixing could offer additionalinsight. This approach could take advantage of known surface tension ef-fects that affect agglomeration and filament formation. The model couldalso be modified to incorporate the assertion by Kohijiya that aggregatescloser than 3 nm tend to fuse into agglomerates and that freshly brokenagglomerates tend to separate by 3 nm [5].A series of experiments could be done to determine whether hopping ortunnelling is dominant. The processes could be differentiated by observingthe temperature variance of conduction. The temperature dependence ofhopping may reveal whether or not this phenomenon warrants inclusion inmodels of this sort.One caution though is that the incorporation of any poorly understoodparameters can lead to greater uncertainty in the model and may ultimatelyover-tune it in the same way that adding additional high order componentsto a line fit may seem at first to give better results even if it is not reallyrevealing more about the governing mechanisms.5.3.2 Post-curing analysis of carbon blackFurther analysis of the size and shape of carbon aggregates and agglomeratesthat have undergone a mixing and curing process is warranted.Despite difficulties slicing thin samples of many conductive compositematerials, it is ultimately possible to cool a sample sufficiently for thin slicesto be made. Combining a cryostat microtome with high resolution imagingtools could shed substantial light on the structures formed by carbon blackin the scenarios being considered by this model. A challenge in doing thisis that many cryotomes have a lowest operating temperature near -40 ? asthey are designed for biological systems. In those systems softer rubberssuch polydimethylsiloxane may not be sufficiently frozen to produce a cleancut. Once that obstacle is overcome and the corresponding digital images areanalysed, the resultant numerical data describing a range of micro or nanoscale physical parameters for a suite of composites could prove instrumentalin the further tuning of this model or the creation of a new one.585.3. Future workTo gain further insight into the carbon black it may be valuable to gen-tly dissolve the polymer component of a cured sample and then image theprecipitated carbon particles with as little disturbance as possible. A com-parison with imaging and size profiling of the same carbon done pre-curingmay prove insightful.5.3.3 Alternate techniquesAt the outset it was hoped this model could be used as a predictive tool toanticipate the electrical characteristics of conductive composites and thatit could then be further refined to explore the importance of deformationof irregular geometries in future applications. The approach taken was onethat captured and described the inter-particle charge transport phenomenonat work within the composite. As a result of this choice it does offer insightinto those physics but it falls short of being a true predictive model. Whilea truly predictive model that relies solely on governing physics may be themost satisfying solution, the present gaps in understanding suggest there iswork to be done before that model exists. In the interim it may be usefulto divide the approach in two. The present model can be expanded andrefined to become a better representation of the science behind conductivecomposites while a separate model may be useful as a purely predictivetool. The purely predictive tool could potentially take more of a ?black-box? approach that did not consider physics but rather simply consideredinput and output parameters and is trained against known data using thewide range of adaptive machine learning algorithms available in the field ofcomputer science.59Bibliography[1] X.-W. Zhang, Y. Pan, Q. Zheng, and X.-S. Yi, ?Time dependence ofpiezoresistance for the conductor-filled polymer composites,? Journalof Polymer Science Part B: Polymer Physics, vol. 38, pp. 2739?2749,Nov. 2000.[2] S. Wang, P. Wang, and T. Ding, ?Resistive viscoelasticity of siliconerubber/carbon black composite,? Polymer Composites, vol. 32, pp. 29?35, Jan. 2011.[3] J. Meier, J. Mani, and M. Klu?ppel, ?Analysis of carbon black network-ing in elastomers by dielectric spectroscopy,? Physical Review B, vol. 75,p. 054202, Feb. 2007.[4] D. M. Heyes, M. Cass, and A. C. Branka, ?Percolation threshold ofhard-sphere fluids in between the soft-core and hard-core limits,? Molec-ular Physics, vol. 104, pp. 3137?3146, Oct. 2006.[5] S. Kohjiya, A. Kato, and Y. Ikeda, ?Visualization of nanostructure ofsoft matter by 3D-TEM: Nanoparticles in a natural rubber matrix,?Progress in Polymer Science, vol. 33, pp. 979?997, Oct. 2008.[6] J. Bang, K. Han, S. Skrabalak, H. Kim, and K. Suslick, ?PorousCarbon Supports Prepared by Ultrasonic Spray Pyrolysis for DirectMethanol Fuel Cell Electrodes,? Journal of Physical Chemistry C,vol. 111, pp. 10959?10964, July 2007.[7] X. Z. Niu, S. L. Peng, L. Y. Liu, W. J. Wen, and P. Sheng, ?Char-acterizing and Patterning of PDMS-Based Conducting Composites,?Advanced Materials, vol. 19, pp. 2682?2686, Sept. 2007.[8] R. S.-P., K. F.-H., and C. K.-C., ?Dispersion of carbon black in acontinuous phase: Electrical, rheological, and morphological studies,?Colloid & Polymer Science, vol. 280, pp. 1110?1115, Dec. 2002.[9] A. Linklater, ?Characterisation of Vulcan XC72 in PDMS composites,?UBC Internal report, 2012.60Bibliography[10] Z. H. Li, ?Effects of carbon blacks with various structures on vulcan-ization and reinforcement of filled ethylene-propylene-diene rubber,?eXPRESS Polymer Letters, vol. 2, pp. 695?704, Sept. 2008.[11] J.-C. Huang, ?Carbon black filled conducting polymers and polymerblends,? Advances in Polymer Technology, vol. 21, pp. 299?313, Jan.2002.[12] J.-C. Huang and C.-L. Wu, ?Processability, mechanical properties, andelectrical conductivities of carbon black-filled ethylene-vinyl acetatecopolymers,? Advances in Polymer Technology, vol. 19, pp. 132?139,Jan. 2000.[13] Z. Y. Liu, J. L. Zhang, P. T. Yu, J. X. Zhang, R. Makharia, K. L.More, and E. a. Stach, ?Transmission Electron Microscopy Observationof Corrosion Behaviors of Platinized Carbon Blacks under Thermal andElectrochemical Conditions,? Journal of The Electrochemical Society,vol. 157, no. 6, p. B906, 2010.[14] J. Strumpler, R. and Glatz-Reichenbach, ?Conducting Polymer Com-posites,? Journal of Electroceramics, vol. 3, no. 4, pp. 329?346, 1999.[15] M. Knite, V. Teteris, A. Kiploka, and J. Kaupuzs, ?Polyisoprene-carbonblack nanocomposites as tensile strain and pressure sensor materials,?Sensors and Actuators A: Physical, vol. 110, pp. 142?149, Feb. 2004.[16] G. Ambrosetti, C. Grimaldi, I. Balberg, T. Maeder, A. Danani,and P. Ryser, ?Solution of the tunneling-percolation problem in thenanocomposite regime,? Physical Review B, vol. 81, p. 155434, Apr.2010.[17] I. Balberg, ?A comprehensive picture of the electrical phenomena incarbon blackpolymer composites,? Carbon, vol. 40, pp. 139?143, Feb.2002.[18] J. Kuba?t, R. Kuz?el, I. Kivka, P. Bengtsson, J. Prokes?, and O. Ste-fan, ?New conductive polymeric systems,? Synthetic Metals, vol. 54,pp. 187?194, Mar. 1993.[19] I. Balberg, ?Tunneling and nonuniversal conductivity in composite ma-terials,? Physical Review Letters, vol. 59, pp. 1305?1308, Sept. 1987.61Bibliography[20] R. Newnham, D. Skinner, and L. Cross, ?Connectivity andpiezoelectric-pyroelectric composites,? Materials Research Bulletin,vol. 13, pp. 525?536, May 1978.[21] S. Kirkpatrick, ?Percolation and Conduction,? Reviews of ModernPhysics, vol. 45, pp. 574?588, Oct. 1973.[22] W. Y. Hsu, W. G. Holtje, and J. R. Barkley, ?Percolation phenomenain polymer/carbon composites,? Journal of Materials Science Letters,vol. 7, pp. 459?462, May 1988.[23] S. Vionnet-Menot, C. Grimaldi, T. Maeder, S. Stra?ssler, and P. Ryser,?Tunneling-percolation origin of nonuniversality: Theory and experi-ments,? Physical Review B, vol. 71, p. 064201, Feb. 2005.[24] W. Bauhofer and J. Z. Kovacs, ?A review and analysis of electrical per-colation in carbon nanotube polymer composites,? Composites Scienceand Technology, vol. 69, pp. 1486?1498, Aug. 2009.[25] W. Zhang, A. a. Dehghani-Sanij, and R. S. Blackburn, ?Carbon basedconductive polymer composites,? Journal of Materials Science, vol. 42,pp. 3408?3418, Mar. 2007.[26] X. Zhang, Y. Pan, Q. Zheng, and X. Yi, ?Piezoresistance of conductorfilled insulator composites,? Polymer International, vol. 50, pp. 229?236, Feb. 2001.[27] J. G. Simmons, ?Low-Voltage Current-Voltage Relationship of TunnelJunctions,? Journal of Applied Physics, vol. 34, no. 1, p. 238, 1963.[28] J. G. Simmons, ?Electric Tunnel Effect between Dissimilar ElectrodesSeparated by a Thin Insulating Film,? Journal of Applied Physics,vol. 34, no. 9, p. 2581, 1963.[29] J. G. Simmons, ?Generalized Formula for the Electric Tunnel Effect be-tween Similar Electrodes Separated by a Thin Insulating Film,? Journalof Applied Physics, vol. 34, no. 6, p. 1793, 1963.[30] J. G. Simmons and G. J. Unterkofler, ?Potential Barrier Shape De-termination in Tunnel Junctions,? Journal of Applied Physics, vol. 34,no. 6, p. 1828, 1963.[31] P. Bagwell and T. Orlando, ?Landauers conductance formula and itsgeneralization to finite voltages,? Physical Review B, vol. 40, pp. 1456?1464, July 1989.62Bibliography[32] F. Bassani, G. L. Liedl, P. Wyder, J. Gomez-Herrero, and R. Reifen-berger, ?Scanning Probe Microscopy,? Encyclopedia of Condensed Mat-ter Physics, pp. 172?182, 2005.[33] P. K. Hansma and J. Tersoff, ?Scanning tunneling microscopy,? Journalof Applied Physics, vol. 61, no. 2, p. R1, 1987.[34] C. J. Chen, Introduction to Scanning Tunneling Microscopy. OxfordUniversity Press, 2 ed., 2008.[35] D. Bloor, K. Donnelly, P. J. Hands, P. Laughlin, and D. Lussey, ?Ametalpolymer composite with unusual properties,? Journal of PhysicsD: Applied Physics, vol. 38, pp. 2851?2860, Aug. 2005.[36] G. Ambrosetti, N. Johner, C. Grimaldi, T. Maeder, P. Ryser, andA. Danani, ?Electron tunneling in conductor-insulator composites withspherical fillers,? Journal of Applied Physics, vol. 106, no. 1, p. 016103,2009.[37] S. Kato, K. Hashimoto, and K. Watanabe, ?Microbial interspecies elec-tron transfer via electric currents through conductive minerals.,? Pro-ceedings of the National Academy of Sciences of the United States ofAmerica, vol. 109, pp. 10042?6, June 2012.[38] L. E. Levine, G. G. Long, J. Ilavsky, R. a. Gerhardt, R. Ou, andC. a. Parker, ?Self-assembly of carbon black into nanowires that form aconductive three dimensional micronetwork,? Applied Physics Letters,vol. 90, no. 1, p. 014101, 2007.[39] H.-P. Xu, Z.-M. Dang, D.-H. Shi, and J.-B. Bai, ?Remarkable selectivelocalization of modified nanoscaled carbon black and positive temper-ature coefficient effect in binary-polymer matrix composites,? Journalof Materials Chemistry, vol. 18, no. 23, p. 2685, 2008.[40] J. Fro?hlich, W. Niedermeier, and H.-D. Luginsland, ?The effect of filler-filler and fillerelastomer interaction on rubber reinforcement,? Compos-ites Part A: Applied Science and Manufacturing, vol. 36, pp. 449?460,Apr. 2005.[41] O. A. Al-Hartomy, F. Al-Solamy, A. Al-Ghamdi, N. Dishovsky,M. Ivanov, M. Mihaylov, and F. El-Tantawy, ?Influence of CarbonBlack Structure and Specific Surface Area on the Mechanical and Di-electric Properties of Filled Rubber Composites,? International Journalof Polymer Science, vol. 2011, pp. 1?8, 2011.63Bibliography[42] T. Koga, T. Hashimoto, M. Takenaka, K. Aizawa, N. Amino, M. Naka-mura, D. Yamaguchi, and S. Koizumi, ?New Insight into HierarchicalStructures of Carbon Black Dispersed in Polymer Matrices: A Com-bined Small-Angle Scattering Study,? Macromolecules, vol. 41, pp. 453?464, Jan. 2008.[43] T. Koga, M. Takenaka, K. Aizawa, M. Nakamura, and T. Hashimoto,?Structure factors of dispersible units of carbon black filler in rubbers.,?Langmuir : the ACS journal of surfaces and colloids, vol. 21, pp. 11409?13, Nov. 2005.[44] S. Kohjiya, A. Katoh, T. Suda, J. Shimanuki, and Y. Ikeda, ?Visual-isation of carbon black networks in rubbery matrix by skeletonisationof 3D-TEM image,? Polymer, vol. 47, pp. 3298?3301, May 2006.[45] N. Kchit and G. Bossis, ?Electrical resistivity mechanism in magne-torheological elastomer,? Journal of Physics D: Applied Physics, vol. 42,p. 105505, May 2009.[46] W. F. Verhelst, K. G. Wolthuis, A. Voet, P. Ehrburger, and J. B.Donnet, ?The Role of Morphology and Structure of Carbon Blacksin the Electrical Conductance of Vulcanizates,? Rubber Chemistry andTechnology, vol. 50, pp. 735?746, Sept. 1977.[47] T. Brown, H. Lemay, and B. Bursten, Chemistry: The Central Science.Simon & Schuster, 1997.[48] Cabot Corporation, ?Vulcan XC72 Data Sheet,? 2004.[49] C. W. Lee, ?Effect of Binary Conductive Agents in Spinel LiMn2O4Cathodes on Electrochemical Performance of Li-ion Batteries,? Journalof industrial and engineering chemistry, vol. 12, no. 6, pp. 967?971,2006.[50] A. Clague, J. Donnet, T. Wang, and J. Peng, ?A comparison of dieselengine soot with carbon black,? Carbon, vol. 37, pp. 1553?1565, Jan.1999.[51] A. Kato, J. Shimanuki, S. Kohjiya, and Y. Ikeda, ?Three-DimensionalMorphology of Carbon Black in NR Vulcanizates as Revealed by 3D-Tem and Dielectric Measurements,? Rubber Chemistry and Technology,vol. 79, pp. 653?673, Sept. 2006.64[52] D. S. McLachlan, ?A grain consolidation model for the critical or per-colation volume fraction in conductor-insulator mixtures,? Journal ofApplied Physics, vol. 70, no. 7, p. 3681, 1991.[53] D. S. McLachlan, ?Analytical Functions for the dc and ac Conductivityof Conductor-Insulator Composites,? Journal of Electroceramics, vol. 5,no. 2, pp. 93?110, 2000.[54] A. Jean, D. Jeulin, S. Forest, S. Cantournet, and F. N?guyen, ?A mul-tiscale microstructure model of carbon black distribution in rubber,?Journal of microscopy, vol. 241, pp. 243?60, Mar. 2011.[55] S. Kumar and S. K. Kurtz, ?Simulation of material microstructure usinga 3D voronoi tesselation: Calculation of effective thermal expansion co-efficient of polycrystalline materials,? Acta Metallurgica et Materialia,vol. 42, pp. 3917?3927, Dec. 1994.[56] Z. Fan, Y. Wu, X. Zhao, and Y. Lu, ?Simulation of polycrystallinestructure with Voronoi diagram in Laguerre geometry based on randomclosed packing of spheres,? Computational Materials Science, vol. 29,pp. 301?308, Mar. 2004.[57] D. F. Watson, ?Computing the n-dimensional Delaunay tessellationwith application to Voronoi polytopes,? The Computer Journal, vol. 24,pp. 167?172, Feb. 1981.[58] P. Dong, ?Generating and updating multiplicatively weighted Voronoidiagrams for point, line and polygon features in GIS,? Computers &Geosciences, vol. 34, pp. 411?421, Apr. 2008.[59] G. Ambrosetti, N. Johner, C. Grimaldi, A. Danani, and P. Ryser, ?Per-colative properties of hard oblate ellipsoids of revolution with a softshell,? Physical Review E, vol. 78, p. 061126, Dec. 2008.[60] I. Balberg and N. Binenbaum, ?Invariant properties of the percola-tion thresholds in the soft-corehard-core transition,? Physical ReviewA, vol. 35, pp. 5174?5177, June 1987.65Appendix ASimulation parametersContained herein is a list of relevant modelling parameters for figures gen-erated using the model described throughout this document.Parameters for Figure 2.5? = 3%r = 65 nm (log normal deviation target = 15%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVvdens = 5? 106 mm?3vrMax = 1.0625 ?mri = 100 nmParameters for Figure 4.1? = 3%r = 62.5 nm (log normal deviation target = 10%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVvdens = 10? 106 mm?3vrMax = 187.5 nmri = 275 nm66Appendix A. Simulation parametersParameters for Figure 4.2? = 3%r = 62.5 nm (log normal deviation target = 10%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVParameters for Figure 4.3? = 3%r = 62.5 nm (log normal deviation target = 10%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVri = varies as per x-axisParameters for Figure 4.4? = 3%r = 62.5 nm (log normal deviation target = 10%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVParameters for Figure 4.5? = 3%r = 62.5 nm (log normal deviation target = 10%)G0 = 4.62 ?S??1 = varies as per figure description67Appendix A. Simulation parameters? = varies as per figure descriptionvdens = 108 mm?3vrMax = 187.5 nmri = varies as per x-axisParameters for Figure 4.8? = 3%r = 65 nm (log normal deviation target = 15%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVri = varies as per x-axisenforced proximity = 25 nmprobability of enforced proximity = varies as per legendParameters for Figure 4.9? = 3%r = 65 nm (log normal deviation target = 15%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVri = varies as per x-axisenforced proximity = 100 nmprobability of enforced proximity = varies as per legend68Appendix A. Simulation parametersParameters for Figure 4.10? = 3%r = 65 nm (log normal deviation target = 15%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVvdens = 5? 106 mm?3vrMax = 187.5 nmri = 100 nmParameters for Figure 4.11? = varies as per x-axisr = 65 nm (log normal deviation target = 15%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meVvdens = 5? 106 mm?3vrMax = varies as per legendri = 100 nmParameters for Figure 4.14? = varies as per legendr = 65 nm (log normal deviation target = 15%)G0 = 4.62 ?S??1 = 9.76 nm? = 0.1 meV69Appendix A. Simulation parametersvdens = 5? 106 mm?3vrMax = 1.0625 ?mri = varies as per x-axis70
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Modelling electrical transport in nano-particle composites...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Modelling electrical transport in nano-particle composites using 3D network analysis Williams, Kaan Owen 2014
pdf
Page Metadata
Item Metadata
Title | Modelling electrical transport in nano-particle composites using 3D network analysis |
Creator |
Williams, Kaan Owen |
Publisher | University of British Columbia |
Date Issued | 2014 |
Description | Conductive composites consist of a conductive filler dispersed within an insulating matrix. These composite materials have been known for many years and are regularly produced experimentally and commercially for a variety of applications. Novel techniques are now being found for creating composites that exhibit conductivity with less conductive filler material than classical physics suggests is sufficient if the particles are uniformly distributed. Several parties have offered physical explanations for the characteristics of their composites by incorporating a blend of classical and quantum physics but few attempts have been made to compare explanations or develop any mechanism to simulate the physics. The model presented in the present work incorporates first principles physics and semi-empirical theory to account for the distribution of particles within a composite and calculate resultant conductivity using three dimensional network analysis. Results from several model iterations are presented and they are compared with published experimental results. The model demonstrates that a random distribution of spherical particles smaller than 200 nm at 3% loading, given realistic wave function decay rates and reasonable tunnelling barrier heights, cannot explain experimentally observed conductivities in these composite materials. The final model, using a Voronoi tessellation approach, duplicates the behaviour trend of the composites being simulated and illustrates some gaps in the present material science knowledge of conductive composites. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2014-02-28 |
Provider | Vancouver : University of British Columbia Library |
Rights | CC0 1.0 Universal |
DOI | 10.14288/1.0165886 |
URI | http://hdl.handle.net/2429/46120 |
Degree |
Master of Applied Science - MASc |
Program |
Electrical and Computer Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2014-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/publicdomain/zero/1.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2014_spring_williams_kaan.pdf [ 1.09MB ]
- Metadata
- JSON: 24-1.0165886.json
- JSON-LD: 24-1.0165886-ld.json
- RDF/XML (Pretty): 24-1.0165886-rdf.xml
- RDF/JSON: 24-1.0165886-rdf.json
- Turtle: 24-1.0165886-turtle.txt
- N-Triples: 24-1.0165886-rdf-ntriples.txt
- Original Record: 24-1.0165886-source.json
- Full Text
- 24-1.0165886-fulltext.txt
- Citation
- 24-1.0165886.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0165886/manifest