On the generalized method of cells and the prediction of effective elastic properties of polymer bonded explosives

0.61 MB·English

by  Biswajit Banerjee

On the generalized method of cells and the prediction of
effective elastic properties of polymer bonded explosives
2
1 BiswajitBanerjee 1 and Daniel O.Adams
0
2

n
a
Dept. ofMechanicalEngineering,UniversityofUtah,
J

1
50SCentralCampusDrive#2202,SaltLakeCity,UT84112,USA.
1

]
i
c
s
-
l
r
t
m
.
t
a
m
-
d
n
o
c
[

1
v
5
3
4
2
.
1
0
2
1
:
v
i
X
r
a
1Correspondingauthor.E-mail:[email protected] Phone:[phone] Fax:[phone]
Abstract
Thepredictionoftheeffectiveelasticpropertiesofpolymerbondedexplosivesusingdirectnumericalsimulationsis
computationallyexpensivebecause of the high volume fraction of particlesin these particulate composites(∼0.90)
and the strong modulus contrast between the particles and the binder (∼20,000). The generalized method of cells
(GMC)isanalternativetodirectnumericalsimulationsforthedeterminationofeffectiveelasticpropertiesofcompos-
ites. GMChasbeenshownto bemorecomputationallyefficientthanfiniteelementanalysisbasedapproachesfora
rangeofcomposites.Inthisinvestigation,theapplicabilityofGMCtothedeterminationofeffectiveelasticproperties
ofpolymerbondedexplosivesisexplored.GMCisshowntogenerateexcellentestimatesofeffectivemoduliforcom-
positescontainingsquarearraysofdisksatvolumefractionslessthan0.60andamoduluscontrastofapproximately
100. However,forhighvolumefractionandstrongmoduluscontrastpolymerbondedexplosivessuchasPBX9501,
the elastic propertiespredictedbyGMCare foundto beconsiderablylowerthanfinite elementbasedestimatesand
experimentaldata. Simulationsofmodelmicrostructuresareperformedtoshowthatnormalstiffnessesareunderes-
timatedbyGMCwhenstress-bridgingduetocontactbetweenparticlesisdominant. Additionally,thecomputational
efficiency of GMC decreases rapidly with an increase in the number of subcells used to discretize a representative
volumeelement.TheresultspresentedinthisworksuggestthatGMCmaynotbesuitableforcalculatingtheeffective
elastic propertiesof high volume fractionand strong moduluscontrast particulatecomposites. Finally, a real-space
renormalization group approach called the recursive cell method (RCM) is explored as an alternative to GMC and
showntoprovideimprovedestimatesoftheeffectivepropertiesofmodelsofpolymerbondedexplosives.
Keywords:EffectiveProperties,ParticulateComposites,HighVolumeFraction,StrongModulusContrast,Stress-
Bridging,MethodofCells,Real-SpaceRenormalizationGroup
1 Introduction
Thegeneralizedmethodofcells(GMC) (Aboudi1996,Paley&Aboudi1992)isasemi-analyticalmethodofdeter-
miningtheeffectivepropertiesofcomposites.Inthismethod,arepresentativevolumeelement(RVE)ofthecomposite
1
underconsiderationisdiscretizedintoaregulargridofsubcells. Equilibriumandcompatibilityaresatisfiedonanav-
eragebasisacrosssubcellsusingintegralsoversubcellboundaries. GMCgeneratesamatrixofalgebraicexpressions
containinginformationaboutsubcellmaterialproperties. Theeffectivestiffnessofthecompositecanbeobtainedby
invertingthismatrix.
OneadvantageofGMCoverothernumericaltechniquesisthatthefullsetofeffectiveelasticpropertiesofacom-
positecanbecalculatedinonestepinsteadofsolvinganumberofboundaryvalueproblemswithdifferentboundary
conditions. GMChasalso beenfoundto be morecomputationallyefficientthatfinite elementcalculationsforfiber
reinforcedcomposites(Aboudi1996,Wilt1995),sincefarfewerGMCsubcellsthanfiniteelementsarenecessaryto
obtainthesamedegreeofaccuracy. Theproblemofdiscretizationisalsominimizedsincearegularrectangulargrid
isusedinGMC.
Thegeneralizedmethodofcellsisdiscussedbrieflyinthiswork.Effectivestiffnessespredictedbythismethodare
comparedwithaccuratenumericalpredictionsforsquarearraysofdisks(Greengard&Helsing1998).Themethodis
thenappliedtotwo-dimensionalmodelsofageneralpolymerbondedexplosiveandtothemicrostructureofPBX9501
using a two-step proceduresimilar to thatof Low et al. (Low et al. 1994). GMC estimates of elastic propertiesare
comparedwithpredictionsfromdetailedfiniteelementcalculations.TheperformanceofGMCisexploredforseveral
microstructureswithcontactingparticlesandsomeshortcomingsofthemethodareidentified. Finally,therecursive
cellmethod(RCM)isexploredasanalternativetodirectGMCcalculationsinthepredictionofeffectivepropertiesof
polymerbondedexplosives.
2 The generalized method of cells
Figure 1 shows a schematic of the RVE, the subcells and the notation (Aboudi 1991) used in GMC. In the figure,
(X ,X ,X )istheglobalcoordinatesystemoftheRVEand(x(α),x(β),x(γ))isthecoordinatesystemlocaltoasubcell
1 2 3 1 2 3
denotedby(αβγ). Itisassumedthatthedisplacementfunctionu(αβγ)varieslinearlywithinasubcell(αβγ)andcan
i
bewrittenintheform
u(αβγ)(x(α),x(β),x(γ))=w(αβγ)(X ,X ,X )+Φ(αβγ)x(α)+Θ(αβγ)x(β)+Ψ(αβγ)x(γ) (1)
i 1 2 3 i 1 2 3 i 1 i 2 i 3
2
where i represents the coordinate direction and takes the values 1, 2 or 3; w(αβγ) is the mean displacement at the
i
centerofthesubcell(αβγ);andΦ(αβγ),Θ(αβγ),andΨ(αβγ)areconstantslocaltothesubcellthatrepresentgradients
i i i
ofdisplacementacrossthesubcell. Thestrain-displacementrelationsforthesubcellaregivenby
1
ǫ(αβγ) = (∂ u(αβγ)+∂ u(αβγ)) (2)
ij 2 i j j i
where∂ = ∂/∂x(α),∂ = ∂/∂x(β),and∂ = ∂/∂x(γ). Sincepolymerbondedexplosivesareisotropicparticulate
1 1 2 2 3 3
composites,thefollowingbriefdescriptionofGMCassumesthattheRVEiscubicandallsubcellsareofequalsize. If
eachsubcell(αβγ)hasthesamedimensions(2h,2h,2h)thentheaveragestraininthesubcellisdefinedasavolume
averageofthestrainfieldoverthesubcellas
1 h h h
ǫ(αβγ) = ǫ(αβγ)dx(α)dx(β)dx(γ) (3)
D ij ES 8h3 Z−hZ−hZ−h ij 1 2 3
Theaveragestraininthesubcellcanbeobtainedintermsofthedisplacementfieldvariables. Itisassumedthatthere
iscontinuityoftractionattheinterfaceoftwosubcells. Thedisplacementsandtractionsareassumedtobeperiodic
at the boundariesthe RVE. Applying the displacementcontinuity equationson an average basis over the interfaces
betweensubcells,theaveragestrainintheRVEcanbeexpressedintermsofthesubcellstrains. Theaveragesubcell
stressescanbeobtainedfromthesubcellstrainsusingthetractioncontinuityconditionandthestress-strainrelations
ofthematerialsinthesubcells. ArelationshipbetweenthesubcellstressesandtheaveragestrainsintheRVEisthus
obtained.
Fororthotropic,transverselyisotropicorisotropicmaterials,theapproachdiscussedaboveleadstothedecoupling
ofthenormalandshearresponseoftheRVE.Thisdecouplingleadstotwosystemsofequationsrelatingthesubcell
stresses and the average strains in the RVE. For the normal components of strain, the system of equations can be
3
writtenas
M1 M1 M1 T H 0 0
 1 2 3 1      
M2 M2 M2 T =2N0hǫ i +2NHhǫ i +2N0hǫ i (4)
 1 2 3 2   11 V   22 V   33 V
        
        
M3 M3 M3 T  0 0 H
 1 2 3 3      
        
whereN isthenumberofsubcellspersideoftheRVE.Thecorrespondingsystemofequationsfortheshearcompo-
nentsisoftheform
M4 0 0 T H 0 0
  12      
 0 M5 0 T =2N20hǫ i +2N2Hhǫ i +2N20hǫ i (5)
  23   12 V   23 V   13 V
        
        
 0 0 M6T  0 0 H
  13      
        
In equations (4) and (5) the M matrices contain material compliance terms. The T matrices contain the average
subcellstresses. ThevectorHcontainsthedimensionsofthesubcells.Thus,asparsesystemofequationsofsize3N2
isproducedthatrelatesthesubcellstressestotheaveragestrainsintheRVE.Afterinvertingtheseequationsandwith
some algebraicmanipulation,explicitalgebraicexpressionsfor the individualtermsofthe effectivestiffnessmatrix
canbeobtained. Thesestress-strainequationsthatrelatetheaverageRVE stressestotheaverageRVEstrainsareof
theform
hσ i Ceff Ceff Ceff 0 0 0 hǫ i
 11 V  11 12 13  11 V 
hσ22iV C1ef2f C2ef2f C2ef3f 0 0 0  hǫ22iV 
    
    
hσ33iV C1ef3f C2ef3f C3ef3f 0 0 0  hǫ33iV 
 =   (6)
    
hσ23iV  0 0 0 C4ef4f 0 0 2hǫ23iV
    
    
hσ31iV  0 0 0 0 C5ef5f 0 2hǫ31iV
    
    
hσ12iV  0 0 0 0 0 C6ef6f2hǫ12iV
    
where Ceff are the terms of the effective stiffness matrix. Details of the algebraic expressions for these terms have
ij
beenpublishedbyotherresearchers(Pindera&Bednarcyk1997).
InGMC,thenumberofequationstobesolvedequalsthenumberofsubcellsraisedtothedthpower,wheredis
4
thenumberofdimensionsintheproblem. Asaresult,thecomputationalefficiencyofGMCdecreasesasthenumber
ofsubcellsincreases. Thisissuehasbeenpartiallyresolved(Orozco1997)byidentifyingthesparsitycharacteristics
of the system of equations and by using the Harwell-Boeing suite of sparse solvers. The computational efficiency
of GMC has been further improvedafter a reformulation(Pindera & Bednarcyk1997, Bednarcyk& Pindera 1997)
that takes advantage of the continuity of tractions across subcells to obtain a system of O(N2) equations in three
dimensions.
Due to decoupling of the normal and shear response of the RVE, the shear components of the stiffness matrix
obtainedfromGMCaretheharmonicmeansofthesubcellshearstiffnessesandoftheform
N N N
1/Ceff =1/N3 1/C(αβγ). (7)
66 66
XXX
α=1β=1γ=1
BednarcykandArnold(2001)suggestthatthislackofcouplingmakesforan“ultra-efficient”micromechanicsmodel.
However,thislackofcouplingcanleadtogrossunderestimationofshearmoduliforhighvolumefractionandhigh
moduluscontrastmaterialssuchaspolymerbondedexplosives.Recently,researchers(Williams&Aboudi1999,Gan
et al. 2000)have attempted to solve the problemby using higher order expansionsfor the displacementand by ex-
plicitlysatisfyingbothsubcellequilibriumandcompatibility. However,theseapproachesdecreasethecomputational
efficiencyofGMCconsiderablyandarenotexploredinthiswork.
3 Validation - square arrays of disks
Inthissection, estimatesofeffectivepropertiesfromGMCarecomparedwithaccuratenumericalresultsforsquare
arraysofdisks. SquareRVEscontainingsquarearraysofdisksexhibitsquaresymmetry. Thetwo-dimensionallinear
elasticstress-strainrelationfortheseRVEscanbewrittenas
hσ i K +µ(1) K −µ(1) 0 hǫ i
 11 V  eff eff eff eff  11 V 
hσ22iV=Keff−µ(ef1f) Keff+µ(ef1f) 0  hǫ22iV  (8)
    
    
hσ12iV  0 0 µ(ef2f)2hǫ12iV
    
5
whereK isthetwo-dimensionaleffectivebulkmodulus,µ(1) istheeffectiveshearmoduluswhenashearstressis
eff eff
appliedalongthediagonalsoftheRVE,andµ(2)istheeffectiveshearmoduluswhenashearstressisappliedalongthe
eff
edgesoftheRVE.Thesethreeeffectivemodulihavebeendeterminedaccurately,usinganintegralequationapproach,
byGreengardandHelsing(1998)forsquarearraysofdiskscontainingarangeofdiskvolumefractions.
TocomparetheeffectivemodulipredictedbyGMCwiththosefromtheintegralequationcalculations(Greengard
&Helsing1998),RVEscontainingdiskvolumefractionsfrom0.10to0.70werecreated.TheseRVEswerediscretized
into64×64equalsizedsubcells. TheeffectivestiffnessmatrixofeachRVEwascalculatedusingGMC.Finally,the
two-dimensionaleffectivemoduliforeachRVEwerecalculatedfromtheeffectivestiffnessmatrixusingtherelations
K =0.5(Ceff+Ceff), µ(1) =0.5(Ceff−Ceff), µ(2) =Ceff. (9)
eff 11 12 eff 11 12 eff 66
Figure 2 shows the moduli predicted by GMC and those from the integral equation method of Greengard and
Helsing(G&H)fordiskvolumefractionsfrom0.10to0.70. Thematerialpropertiesofthedisksandthebinderused
inthecalculationsareshowninTable1. Theeffectivebulkmoduli(K )anddiagonalshearmoduli(µ(1))obtained
eff eff
fromtheGMCcalculationsarewithin4%ofthoseobtainedbytheintegralequationmethodforallvolumefractions
up to 0.60. At a volume fraction of 0.70, the GMC predictions for bulk modulus and diagonal shear modulus are
4% and 11% less, respectively. For the shear modulusµ(2), the GMC predictions are around 4% to 10% less than
eff
theestimatesofGreengardandHelsingforvolumefractionsfrom0.10to 0.60. Thedifferenceisaround24%fora
volumefractionof0.70.
TheseresultsshowthatGMCestimatesarequiteaccurateforcompositescontainingsquarearraysofdiskswith
volumefractionsupto0.60,confirmingresultsreportedelsewhere(Aboudi1996). Inthenextsection,GMCisused
todeterminetheeffectivepropertiesofmodelsofpolymerbondedexplosivesandtheresultsarecomparedtodetailed
finiteelementcalculationsandexperimentaldata.
6
4 Modeling polymer bonded explosives
Polymer bonded explosive (PBX) materials typically contain around 90% by volume of particles surrounded by a
binder. The particles consist of a mixture of coarse and fine grains with the finer grains forming a filler between
coarsergrains.ModelingthemicrostructureofthesematerialsisdifficultduetothecomplexshapesofHMXparticles
andthelargerangeofparticlesizes. Two-dimensionalapproximationsofthemicrostructureofPBXsbasedondigital
images(Benson & Conley1999)have beenused to studysome aspects of the micromechanicsof PBXs. However,
suchmicrostructuresaredifficulttogenerateandrequirecompleximageprocessingtechniquesandexcellentimage
qualitytoaccuratelycapturedetailsofthematerial.AcombinationofMonteCarloandmoleculardynamicstechniques
havealsobeenusedtogeneratethree-dimensionalmodelsofPBXs(Baer2001). Microstructurescontainingspheres
and oriented cubes have been generated using these techniques and appear to represent PBX microstructures well.
However,thegenerationofmicrostructuresusingdynamics-basedmethodsisextremelytimeconsumingwhentight
particlepackingisrequired,asisthecaseforvolumefractionsabove0.70.
Comparisonsoffiniteelementpredictionswithexactrelationsfortheeffectivepropertiesofcomposites(Banerjee
2002)haveshownthatdetailedfiniteelementestimatescanbeusedasabenchmarktochecktheaccuracyofpredic-
tionsfromGMC.Inthisinvestigation,manuallygeneratedPBXmicrostructurescontainingsymmetricallydistributed
circularparticlesareusedinitiallytocompareGMCandfiniteelementpredictions. Thetwo-dimensionalmicrostruc-
turescontain90%particlesbyvolumeand use two particle lengthscales. Two-dimensionalmodelscontainingran-
domlydistributedcircularparticlesthatreflecttheactualparticlesizedistributionofPBX9501arenextmodeledwith
GMCandtheresultscomparedtofiniteelementestimates.
Thematerialpropertiesusedfortheparticles,thebinder,andPBX9501inthesecalculationsareshowninTable2.
ThesepropertiescorrespondtothoseofHMX(theexplosiveparticles),thebinder,andPBX9501at25oCandastrain
rateof0.05/s(Wetzel1999).
7
4.1 Simplified models ofPBX materials
GMCandfiniteelementcalculationswereperformedforthesix,manuallygenerated,simplifiedmodelmicrostructures
ofpolymerbondedexplosivesshowninFigure3. Theserepresentativevolumeelements(RVEs)containoneorafew
relativelylargeparticlessurroundedbysmaller particlesto reflectcommonparticlesize distributionsof PBXs. The
volumefractionofparticlesineachofthesemodelsisaround0.90±0.005.Thebindermaterialsurroundsallparticles
inthesixmicrostructures.
FortheGMCcalculations,asquaregridwasoverlaidontheRVEstogeneratesubcells. Twodifferentapproaches
were used to assign materials to subcells before the determination of effective properties of the RVE. In the first
approach,referredto as the “binarysubcellapproach”,a subcellwas assignedthe materialpropertiesofparticlesif
morethan50%ofthesubcellwasoccupiedbyparticles.Binderpropertieswereassignedotherwise.Figure4(a)shows
aschematicofthebinarysubcellapproach.Inthesecondapproach,calledthe“effectivesubcellapproach”,amethod
ofcellscalculation(Aboudi1991)wasusedtodeterminetheeffectivepropertiesofasubcellbasedonthecumulative
volumefractionofparticlesinthesubcell(Banerjeeetal.2000).Figure4(b)showsaschematicoftheeffectivesubcell
approach.Afterthesubcellswereassignedmaterialproperties,theGMCtechniquewasusedtocomputetheeffective
propertiesoftheRVE.
Notethattheparticlesarenotresolvedwellwhenmaterialsareassignedtosubcellsinthismannerifthenumber
of subcells is small. However, the large size of the matrix to be inverted in GMC limits the number of subcells
thatcanbeusedto discretizetheRVE. Ifthebinarysubcellapproachisusedtoassignsubcellmaterials, contacting
particlesarecreatedwheretherearenoneintheactualmicrostructure,leadingtothepredictionofhigherthanactual
stiffnessvalues.Theeffectivesubcellapproachimprovesuponthebinarysubcellapproachby“smearing”thematerial
propertiesattheboundariesofparticlesandthusreducingtheparticlecontactartifactscausedbydiscretizationerrors.
ForvalidatingtheGMCresults,detailedfiniteelement(FEM)calculationswereperformedusingsix-nodedtrian-
gularelementstoaccuratelymodelthegeometryoftheparticles.Around65,000nodeswereusedtodiscretizeeachof
themodels.ThevolumeaveragestressandstrainineachRVEwasdeterminedforappliednormalandsheardisplace-
ments. Periodicitywas enforcedthroughdisplacementboundaryconditions. Since these finite elementcalculations
servetovalidatetheGMCcalculations,furthermeshrefinementwasexploredandtheresultswerefoundtoconverge
8
thosefromthe65,000nodefiniteelementmodels.
Table 3 lists the effective stiffnesses of the six RVEs shown in Figure 3 from GMC and FEM calculations. On
average,theGMCcalculationsusingthebinarysubcellapproachpredictvaluesofCeff thatarearound2.5timesthe
11
FEMbasedvalues. ThevaluesofCeff fromtheeffectivesubcellapproachbasedGMCcalculationsareclosertothe
11
FEMestimatesthanthebinarysubcellbasedGMCestimates. TheGMCandFEMestimatesofCeff arequiteclose.
12
ThevaluesofCeff fromGMCareonly10%oftheFEM values. Thelow valuesofCeff areobtainedbecauseGMC
66 66
predicteffectiveshearstiffnessesthatareharmonicmeansofsubcellshearstiffnesses. Models3and4(Figure3)pro-
duceagreementinCeff andCeff betweenthebinarysubcellapproachbasedGMCcalculationsandFEM(within5%)
11 12
whereastheotherfourmodelsproduceconsiderabledifferencesinpredictions. Theselargedifferencesareproduced
bydiscretizationerrorsintroducedintheGMCapproachthatleadtocontinuousstress-bridgingpathsacrosstheRVEs
andhencetoincreasedstiffness.
However, if the predicted effective stiffnesses shown in Table 3 are compared with the experimental effective
stiffnessofPBX 9501(showninTable2),itcanbeobservedthatthemodelspredictvaluesofCeff andCeff thatare
11 12
around10%oftheexperimentalvalues. Hence,thesesimplifiedmodelsarenotappropriateforthemodelingofPBX
9501.ThenextsectionexploresmodelsbasedontheactualparticlesizedistributionofPBX9501.
4.2 Models ofPBX 9501
CoarseandfineparticlesofHMXareblendedinaratioof1:3byweightandcompactedintheprocessofmanufacturing
PBX 9501. Figure 5(a) shows four RVEs of PBX 9501 based on the particle size distribution of the dry blend of
HMX (Wetzel 1999) prior to compaction. Figure 5(b) shows four RVEs based on the particle size distribution of
pressed PBX 9501 (Skidmore et al. 1998). The larger particles are broken up in the pressing process leading to a
larger proportion of smaller particles in pressed PBX 9501. The models of the dry blend have been labeled “DB”
whilethoseofpressedPBX9501havebeenlabeled“PP”.
GMCcalculationswereperformedonthePBX9501RVEsafterdiscretizingeachRVEinto100×100subcellsand
assigningmaterialstosubcellsusingtheeffectivesubcellapproach. InordertovalidatetheGMCpredictions,FEM
calculationswerealsoperformedontheRVEsafterdiscretizingeachRVEinto350×350four-nodedsquareelements.
9