1Finite Element Solutions of Finite Element Solutions of BoundaryBoundary--value Problems in value Problems in ODEsODEsLarry CarettoMechanical Engineering 501ABSeminar in Engineering AnalysisNovember 30, 20042Outline• Review last class – Finite difference method for boundary value problems– Eigenvalue problems• Comparison of shoot and try method and finite-differences• Introduction to finite elements• Finite element solutions to boundary value problems3Review Example• Differential equation: d2T/dx2+ a2T =0• Finite-diff: Ti-1+ (–2 + h2a2)Ti+ Ti+1= 0• Apply boundary conditions from ODE to boundary nodes–T0= T(x=0) = TAand TN= T(x=L) = TB– (–2 + h2a2)T1+ T2= –TA–TN-1+ (–2 + h2a2)TN= –TB• Resulting system of equations forms tridiagonal matrix4Review Example II⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡+−+−+−+−+−−−BANNTTTTTTT00021000012010000121000012100001212321MMMMLLMMMMMOMMMMMOLLLααααα• Finite-difference equations in matrix form with α = a2h25Review Thomas Algorithm• Forward computations– Initial: E0= – C0/ B0F0= D0/ B0– Apply equations below for i = 1,… N-1: – At final point111 −−−+−=+−=iiiiiiiiiiiiEABFADFEABCE11−−+−==NNNNNNNNEABFADFx• Back substitute: xi= Fi+ Eixi+16Review Compact Expressions• Based on operator δ2fi= fi+1–2 fi+ fi-1()4222121'' hOhffi+⎟⎟⎠⎞⎜⎜⎝⎛+=δδ()01210''22222=⎟⎟⎠⎞⎜⎜⎝⎛++⇒=+iiiyahyyayδδ012126512112222122=⎟⎟⎠⎞⎜⎜⎝⎛++⎟⎟⎠⎞⎜⎜⎝⎛−+⎟⎟⎠⎞⎜⎜⎝⎛+−− iiiyahyahyah• Apply to y’’ + a2y = 027Review Compact Expressions II012126512112222122=⎟⎟⎠⎞⎜⎜⎝⎛++⎟⎟⎠⎞⎜⎜⎝⎛−+⎟⎟⎠⎞⎜⎜⎝⎛+−− iiiyahyahyah• Still have tridiagonal array• Results are fourth-order accurate• Have compact expression for first derivative but combination not “compact”()4211612' hOhfffiii+⎟⎟⎠⎞⎜⎜⎝⎛+−=−+δ8Error and Error Order• Get overall measure of error (like norm of a vector)• Typically use maximum error (in absolute value) or root-mean-squared (RMS) error• N = 10 has emax= 0.00242 and eRMS= 0.00183. For N = 100, emax= 0.0000241 and eRMS= 0.0000173.• FDM error is second- or fourth-order error depending on expressions0001.0100.000896.5E-061.77E-060.990.0016011.6E-063.19E-060.880.0021115.0E-064.19E-060.770.0023816.9E-064.73E-060.660.0024216.9E-064.81E-060.550.0022415.5E-064.46E-060.440.0018712.7E-063.72E-060.330.001349.0E-062.66E-060.220.000704.6E-061.39E-060.110000.00FDM 2 ErrorShoot ErrFDM 4 Errxii10Review Boundary Gradients• Use gradients with second order errorhTTTkdxdTkqxx24321000−+−−=−==hTTTkdxdTkqNNNxxNN24321 −−=+−−=−=.00021-.9155.01-.91531.01786-.9332.1-.91531.000362.1999.012.19950.036182.2357.12.19950Error-q/kh-qexact/kx11Review Boundary Conditions• General condition a dT/dx + bT = c– a = 1, b = 0 for Neumann (gradient given)– a = 0, b = 1 for Dirichlet (value given)– Write gradient using second order forward (x = x0) or backward difference (x = xN)– Combine with equation for first node in from the boundary to eliminate term with second node from boundary– Result conforms to tridiagonal system12Review Nonlinear Problems• Shoot-and-try requires no special procedures for nonlinear problems• For finite difference or finite elements, solve a linearized equation– Expand nonlinear terms in Taylor series about current iteration values and keep only linear terms– Example: replace sin θ by the expression sin θ(m)+ cos θ(m)(θ(m+1)– θ(m))313Review Eigenvalue Problems • Eigenvalue problems in differential equations are transformed into matrix eigenvalue problems• The accuracy of the eigenvalues depends on the grid spacing• Accuracy decreases as eigenvalue increases– Can only find as many eigenvalues as there are grid nodes 14Finite Element Methods• Designed for 2D and 3D geometries• Can use for 1D case as example• Basic idea is to divide region into small elements (line, area, volume)• Use interpolating polynomial for each element– Represent both geometry (independent variables) and dependent variable– Polynomials called basis functions15Finite Element Methods II• Analysis for individual elements is assembled into a set of nodal equations for the entire region– Result is set of algebraic equations for the dependent variable at nodes that are points on elements– Just as in finite-difference approach we convert a continuous differential equation to a set of algebraic equations at distinct points in the region16Two-dimensional Element(x4, y4)(x2, y2)(x1, y1)(x3, y3)ξ = -1, η = -1ξ = -1, η = 1ξ = 1, η = 1ξ = 1, η = -1• Use dimensionless ξ-η coordinate system for basis functions• Each element has several shape or basis functions, φi17Shape (Basis) Functions• Model geometry and dependent variable over the element• Use of same basis functions for both is called isoparametric element• Shape functions associated with element nodes such that φi(x(j)) = δij∑∑∑======414141 iiiiiiiiiTTyyxxϕϕϕ18Shape (Basis) Functions II• Simplest shape functions are linear for 1D or bilinear for 2D• For a linear element between nodes i (at ξ = -1) and i + 1 (at ξ = 1) we have φi= (1 - ξ)/2 and φi+1= (1 + ξ)/2•x = xiφi+ xi+1φi+1is correct at nodes• Bilinear functions for 2D element have the form (1 ±ξ) (1 ±η)/2419Bilinear Shape Functions(x4, y4)(x2, y2)(x1, y1)(x3, y3)ξ = -1, η = -1ξ = -1, η = 1ξ = 1, η = 1ξ = 1, η = -14)1)(1(4)1)(1(4)1)(1(4)1)(1(4321ηξϕηξϕηξϕηξϕ+−=++=−+=−−=Note:φi(x(j)) = δij20Modeling Differential Equation• Look at same example used for finite differences: d2T/dx2+ a2T = 0∑==NiiiTT1ˆϕ• Equation for T in terms of basis functions gives approximate value• Seek solution in which differential equation is satisfied in an average way over the region; wiis weighting functionNidxTadxTdxwLiK,00ˆˆ)(0222==⎥⎦⎤⎢⎣⎡+∫21Modeling Differential Equation II• Last equation on previous chart called method of weighted residuals (MWR)• Various choices are used for weighting functions, wi• Galerkin method uses wi= φi• Known to match variational results for linear problemsNidxTadxTdxLiK,00ˆˆ)(0222==⎥⎦⎤⎢⎣⎡+∫ϕ22Modeling