Practical Finite Element Modelling in Earth Science using MatLab

247 Pages·2017·4.88 MB·english

by  Guy Simpson

Practical Finite Element Modeling
in Earth Science
Using Matlab
GuySimpson
DepartmentofEarthScience
UniversityofGeneva
Switzerland
Thiseditionfirstpublished2017©2017JohnWiley&SonsLtd
LibraryofCongressCataloging-in-Publicationdataappliedfor
ISBN:9781119248620
BriefContents
Preface xiii
Symbols xv
AbouttheCompanionWebsite xvii
PartI TheFiniteElementMethodwithMatlab 1
1 Preliminaries 3
2 BeginningwiththeFiniteElementMethod 13
3 ProgrammingtheFiniteElementMethodinMatlab 25
4 NumericalIntegrationandLocalCoordinates 35
5 TheFiniteElementMethodinTwoDimensions 49
6 TheFiniteElementMethodinThreeDimensions 67
7 GeneralizationofFiniteElementConcepts 81
PartII ApplicationsoftheFiniteElementMethodinEarthScience 119
8 HeatTransfer 121
9 LandscapeEvolution 137
10 FluidFlowinPorousMedia 151
11 LithosphericFlexure 167
12 DeformationofEarth’sCrust 183
13 GoingFurther 207
AppendixA DerivationoftheDiffusionEquation 217
AppendixB BasicsofLinearAlgebrawithMatlab 221
AppendixC ComparisonbetweenDifferentNumericalMethods 227
AppendixD IntegrationbyParts 237
AppendixE TimeDiscretization 239
References 241
Index 245
Contents
Preface xiii
Symbols xv
AbouttheCompanionWebsite xvii
PartI TheFiniteElementMethodwithMatlab 1
1 Preliminaries 3
1.1 MathematicalModels 3
1.2 BoundaryandInitialConditions 4
1.3 AnalyticalSolutions 5
1.4 NumericalSolutions 5
1.5 NumericalSolutionMethods 7
1.6 MatlabScript 8
1.7 Exercises 10
SuggestedReading 12
2 BeginningwiththeFiniteElementMethod 13
2.1 TheGoverningPDE 13
2.2 ApproximatingtheContinuousVariable 14
2.3 MinimizingtheResidual 15
2.4 EvaluatingtheElementMatrices 17
2.5 TimeDiscretization 18
2.6 Assembly 19
2.7 BoundaryandInitialConditions 21
2.8 SolutionoftheAlgebraicEquations 21
2.9 Exercises 22
SuggestedReading 23
3 ProgrammingtheFiniteElementMethodinMatlab 25
3.1 ProgramStructureandPhilosophy 25
3.2 SummaryoftheProblem 25
3.3 DiscretizedEquations 26
3.4 TheProgram 27
3.4.1 PreprocessorStage 27
3.4.2 SolutionStage 29
3.4.3 PostprocessorStage 30
3.5 MatlabScript 30
3.6 Exercises 33
SuggestedReading 34
4 NumericalIntegrationandLocalCoordinates 35
4.1 Gauss–LegendreQuadrature 36
4.2 LocalCoordinates 37
4.3 EvaluatingtheIntegrals 39
4.4 VariableMaterialProperties 40
4.5 ProgrammingConsiderations 41
4.6 MatlabScript 43
4.7 Exercises 45
SuggestedReading 47
5 TheFiniteElementMethodinTwoDimensions 49
5.1 Discretization 50
5.2 GeometryandNodalConnectivity 52
5.3 IntegrationofElementMatrices 54
5.4 MultielementAssembly 57
5.5 BoundaryConditionsandSolution 60
5.6 MatlabScript 61
5.7 Exercises 65
SuggestedReading 66
6 TheFiniteElementMethodinThreeDimensions 67
6.1 Discretization 67
6.2 ElementIntegration 69
6.3 AssemblyforMultielementMesh 72
6.4 BoundaryConditionsandSolution 73
6.5 MatlabProgram 74
6.6 Exercises 79
SuggestedReading 80
7 GeneralizationofFiniteElementConcepts 81
7.1 TheFEMforanEllipticProblem 84
7.2 TheFEMforaHyperbolicProblem 96
7.3 TheFEMforSystemsofEquations 102
7.4 Exercises 116
SuggestedReading 116
PartII ApplicationsoftheFiniteElementMethodinEarthScience 119
8 HeatTransfer 121
8.1 ConductiveCoolinginanErodingCrust 122
8.2 ConductiveCoolingofanIntrusion 126
SuggestedReading 135
9 LandscapeEvolution 137
9.1 Evolutionofa1DRiverProfile 138
9.2 EvolutionofaFluviallyDissectedLandscape 143
SuggestedReading 150
10 FluidFlowinPorousMedia 151
10.1 FluidFlowAroundaFault 152
10.2 ViscousFingering 157
SuggestedReading 166
11 LithosphericFlexure 167
11.1 GoverningEquations 167
11.2 FEMDiscretization 168
11.3 MatlabImplementation 171
SuggestedReading 181
12 DeformationofEarth’sCrust 183
12.1 GoverningEquations 183
12.2 RateFormulation 185
12.3 FEMDiscretization 186
12.4 Viscoelastoplasticity 188
12.5 MatlabImplementation 190
SuggestedReading 205
13 GoingFurther 207
13.1 Optimization 207
13.2 UsingOtherFEMs 213
13.3 UseofExistingFiniteElementSoftware 215
AppendixA DerivationoftheDiffusionEquation 217
AppendixB BasicsofLinearAlgebrawithMatlab 221
AppendixC ComparisonbetweenDifferentNumericalMethods 227
AppendixD IntegrationbyParts 237
AppendixE TimeDiscretization 239
References 241
Index 245
Preface
Over the past few decades, mathematical models have become an increasingly important tool for
Earth scientists to understand and make predictions about how our planet functions and evolves
throughtimeandspace.Thesemodelsoftenconsistofpartialdifferentialequations(PDEs)thatare
discretizedwithanumericalmethodandsolvedonacomputer.Themostcommonlyuseddiscretiza-
tionmethodsarethefinitedifferencemethod(FDM),thefinitevolumemethod,thefiniteelement
method (FEM), the discrete element method, the boundary element method, and various spectral
methods.Intheory,eachmethodprovidesthesamesolutiontotheoriginalPDEs.However,inprac-
tice,certainmethodsarebettersuitedtocertainproblemsthanothers.Often,onemethoddominates
withinanygivendisciplineandintheEarthsciences,theFDMisthemostprevalent,duetoitssim-
plicity.AlthoughtheFEMisarguablybettersuitedtomanyEarthscienceproblems—especiallythose
withcomplicatedgeometryand/ormaterialbehavior—Earthscientistshavebeenhesitanttowhole-
heartedlyembracethistechniquebecauseitisregardedasbeingcomplicatedtoimplementcompared
tootherschemes.However,thisperceiveddifficultylargelyreflectsthefactthatmosttextbookson
thismethodarewrittenbyengineersormathematiciansforengineerswhohaveadifferenteduca-
tionalbackgroundascomparedtoEarthscientistsandwhoareinterestedindifferentapplications.
ThisisunfortunatebecausetheFEMisaremarkablyflexibleandpowerfultoolwithenormouspoten-
tialintheEarthsciencesthatisnomoredifficult(oreveneasier)toimplementthanothernumerical
schemes.
ThetextisintendedforstudentsandresearchersinEarthscience,attemptingtheirfirststepswith
theFEM.ItprovidesapracticalguideonhowtheFEMcaneasilybeusedtosolvevariousEarthsci-
enceproblemsusingMatlab.Forthemostpart,Iassumethattheequationsgoverningtheprocesses
ofinterestareknown.EmphasisisonhowoneactuallycomputesthesolutionusingtheFEM.Thetext
doesnotdealindetailwithbenchmarkingandinterpretationofmodelresultsorwithapplicationof
themodelresultstospecificcasestudies.Toguidereaders,manysamplefiniteelementMatlabscripts
arepresented.Thesescriptsarewrittenwithanemphasisonsimplicityandclarity,notonmodularity
andefficiency.However,oncetheunderlyingconceptsareclear,thesestandalonecodescouldeas-
ilybemodularized,optimized,andtransportedtoothermoreefficientlanguagessuchasFortranor
C/C++.ItisassumedthatthereaderisfamiliarwithlinearalgebraandPDEsandhasbasicprogram-
mingexperience.SomeoftheseaspectsarecoveredbrieflyinChapter1andAppendixB.Thetextis
directedtowardgraduatestudents,advancedundergraduates,andEarthscienceresearchers.While
thetextisintendedtoshowhowfiniteelementprogramscanbewrittenfromscratch,itshouldalso
beofinteresttoresearcherswhouseexistingFEMsoftware(e.g.,ABAQUSandCOMSOL)butwho
wanttoknowmoreaboutwhatgoesonwithinthe“blackbox”.Becausethelevelofthematerialpre-
sentedisquitebasiccomparedtootherfiniteelementtexts,readersarestronglyadvisedtoconsult
othermoreadvancedbooksoncetheyunderstandthebasicsandseehowprogramsareconstructed
inpractice.Thefollowingtitlesaregoodstartingpoints:TheFiniteElementMethod(threevolumes,
byZienkiewiczandTaylor,2000a,b,c)andTheFiniteElementMethod(byHughes,2000).
This book is structured in two parts. Part I begins with a general introduction to numerical
modelingbeforepassingtoaseriesofchaptersthatshowhowanarchetypicalmathematicalmodel
(i.e., the diffusion equation) is discretized with the FEM, programmed in Matlab, and solved on a
computer.Eachchapterbuildsonthepreviousoneandintroducesone(ormore)keyaspectofthe
FEM.Chapter7generalizestheconceptsintroducedinChapters1–6byshowinghowtheFEMcan
extendedfromsingleparabolicequationstosystemsofequationsandalsotoellipticandhyperbolic
equations.BytheendofPartI,thereadershouldunderstandtheessentialsoftheFEMandbeable
to write their own Matlab scripts from scratch to solve the most commonly encountered PDEs in
onedimension(1D),twodimension(2Ds),andthreedimension(3Ds).Thismaterialcanbetaught
asaonesemestercourseonnumericalmodelinginEarthscienceformasterandPhDstudents.Part
IIcomprisesaseriesofindependentchapters,eachofwhichfocusesonhowtheFEMcanbeapplied
in different contexts in Earth science. The problems investigated are heat transfer in the crust,
landscapeevolutionmodeling,fluidflowinporousmedia,flexure,anddeformationofEarth’scrust.
Although readers can choose to read only the chapter(s) that fall closest to their topic of interest,
each chapter introduces a different aspect of the FEM, and so every chapter should be studied by
readersinterestedineventuallymasteringthetechnique.
I am very much indebted to Yuri Podladchikov who initially inspired my interest in the subject
presentedandwhocontributedinamajorwaytomyunderstandingoftheFEM.Duringthesame
period, I also benefited enormously from discussions and interaction with many other colleagues
fromtheETHinZurich,including,inparticular,AlanThompson,JamieConnolly,NeilMancktelow,
Jean-PierreBurg,SteveMiller,LuigiBurlini,KatjaPetrini,TarasGerya,StefanSchmalholz,Daniel
Schmid,BorisKaus,andDaveMay.LineProbstfromtheUniversityofGenevaisthankedformany
comments and corrections on the text. I gratefully acknowledge financial support for my research
fromtheSwissNationalScienceFoundation,theDepartmentofEarthScienceattheUniversityof
Geneva,andtheSchmidheinyFoundation.Finally,Ithanktheencouragementandsupportfrommy
parents,mywifeKatja,andmychildrenLuca,Lara,andFabio.
GuySimpson
Geneva
Symbols
Matlabvariable Dimensions Description
b [sdof, 1] Globalright-hand-sidevector
bee [nst, ntot] Kinematicstrain—displacementmatrix
bcdof [1 ndn] ArraycontaininglistofequationswhereDirichlet
boundaryconditionsareimposed
ndn scalar NumberofequationstowhichDirichletconditions
areapplied
bcval [1 ndn] Arraycontainingfixedboundaryvalues
bx0 [1 ny]in2D Arrayscontaininglistofnodesonx=0boundary
bxn [1 ny]in2D Arrayscontaininglistofnodesonx=lxboundary
coord [nod, ndim] Nodecoordinatesforoneelement
dee [nst, nst] Materialdeformationmatrix
der [ndim nod] Derivativesofshapefunctionsinlocalcoordinates,
evaluatedatanintegrationpoint
der_s [ndim nod nip] Derivativesofshapefunctionssavedforall
integrationpoints
deriv [ndim nod] Derivativesofshapefunctionsinglobal
coordinates,evaluatedatanintegrationpoint
detjac scalar DeterminantoftheJacobianmatrix
displ [sdof, 1] Globalsolutionvector
displ0 [sdof, 1] Globalsolutionvectorfromprevioustimestep
dt scalar Timeincrement
dx scalar Elementdimensioninthex-direction
F [ntot, 1] Elementloadvector
ff [sdof, 1] Globalright-hand-sideloadvector
fun [1 nod] Shapefunctions,evaluatedatanintegrationpoint
fun_s [nip nod] Shapefunctionssavedforallintegrationpoints
g [ntot 1] Equationlistforoneelement
g_coord [ndim nn] Nodecoordinatesforentiremesh
g_g [ntot nels] Equationnumbersofeachelementforentiremesh
g_num [nod nels] Nodenumbersofeachelementforentiremesh
invjac [ndim, ndim] InverseoftheJacobianmatrix
jac [ndim, ndim] Jacobianmatrix
KM [ntot, ntot] Elementstiffnessmatrix