Finite Differences Newton’s Method Example: Rectangular Duct With Incompressible Flow MATLAB for Sparse Matrix Adding Convection10.34, Numerical Methods Applied to Chemical Engineering Professor William H. Green Lecture #19: Boundary Value Problems (BVPs) Lecture 2. Finite Differences (2)(2)()(0xOxxxfxxfdxdfooxΔ+ΔΔ−−Δ+≡) Relatively good accuracy, better convergence ()xOxxfxxfdxdfooxΔ+Δ−Δ+≡)()(0 ONE SIDED upwind differencing: The error leads to numerical stability but is a mathematical trick. Adds in error Deff = Dtrue + VxΔx/2 and Pelocal, eff < 2. Still wrong, because artificially increased. V2φ + v·Vφ + S(φ) = 0 ()0)()(21211=+Δ−+Δ+−+−+xfxvxiiiiiφφφφφφ linear linear linear or nonlinear ()()00)(0)(21=⎟⎟⎟⎠⎞⎜⎜⎜⎝⎛+×%xfxfMφφφ ()()00021=⎟⎟⎟⎠⎞⎜⎜⎜⎝⎛+×%φφφffM approx. to differential operators Newton’s Method F(φ) = 0 JΔφ = -F Í Newton update J = M + ⎟⎟⎟⎟⎟⎟⎟⎟⎠⎞⎜⎜⎜⎜⎜⎜⎜⎜⎝⎛∂∂∂∂%0021φφφφff Cite as: William Green, Jr., course materials for 10.34 Numerical Methods Applied to Chemical Engineering, Fall 2006. MIT OpenCourseWare (http://ocw.mit.edu), Massachusetts Institute of Technology. Downloaded on [DD Month YYYY].Example: Rectangular Duct With Incompressible Flow Rectangular duct with incompressible fluid being pulled by gravity: V2vz(x,y) = ρg/μ No slip at walls: vz(boundaries) = 0 2222yx ∂∂+∂∂ jNinnyxyii+−=→(*)1()(),ϕφ ⎟⎟⎟⎟⎟⎟⎠⎞⎜⎜⎜⎜⎜⎜⎝⎛=),(),(),(),(1212111yxvyxvyxvyxvzNzzzy#ϕ •••••••• •↓•→•←•↑• Ny yÇ Δy x Æ Nx Δx 11234+••••••••Ny (1,1) Figure 1. Rectangular duct with incompressible flow. NxNy = # points Interior: 211211)(),(),(2),()(),(),(2),(yyxvyxvyxvxyxvyxvyxvjizjizjizjizjizjizΔ+−+Δ+−−+−+ Rows of M look like this: 10.34, Numerical Methods Applied to Chemical Engineering Lecture 19 Prof. William Green Page 2 of 4 Cite as: William Green, Jr., course materials for 10.34 Numerical Methods Applied to Chemical Engineering, Fall 2006. MIT OpenCourseWare (http://ocw.mit.edu), Massachusetts Institute of Technology. Downloaded on [DD Month YYYY].⎟⎟⎠⎞⎜⎜⎝⎛ΔΔ⎟⎟⎠⎞⎜⎜⎝⎛Δ+Δ−ΔΔ= 00)(100)(100)(1)(12)(1000)(1.....0222222yxyxyxM Mφ = ρg/μ {Mφ = b} Shown in MATLAB function makeAforLaplacian(Nx,Ny,Xmax,Ymax) Non-sparse: A(n, n+1) = 1/(Δy)2Ymax Ç Æ Nx Sparse format: mvec(k) = n ηvec(k) = n+1 A vec(k) = 1/(Δy)2 MATLAB for Sparse Matrix Nx = 20; Ny = 30; Xmax = 40; Ymax = 10; A_sparse = makeA_sparse(Nx,Ny,Xmax,Ymax); b = zeros(Nx*Ny,1); % unknown vector b = b+1; phi = A_sparse\b; check = A_sparse*phi; check(400) ans = 1.0 bbig = 1e7*b phi = A_sparse\bbig; check = A_sparse*phi; ans = 1e7 How do know if correct? V_shaped = reshape(phi, Ny, Nx); surf(V_shaped) Makes a beautiful plot of the solution Adding Convection Peclet number Pe = vL/D Pelocal = vxΔx/D 10.34, Numerical Methods Applied to Chemical Engineering Lecture 19 Prof. William Green Page 3 of 4 Cite as: William Green, Jr., course materials for 10.34 Numerical Methods Applied to Chemical Engineering, Fall 2006. MIT OpenCourseWare (http://ocw.mit.edu), Massachusetts Institute of Technology. Downloaded on [DD Month YYYY].V(vφ) + DV2φ + S(φ) = 0 ⎟⎟⎠⎞⎜⎜⎝⎛−∂∂−−=⎟⎟⎠⎞⎜⎜⎝⎛∂∂=+∂∂+∂∂+∂∂DSxvDDvxSxDxvxvxxxxϕϕϕϕϕϕϕϕϕ1'''0)(22 x∂∂ϕ Recall from ODE discussion: 10.34, Numerical Methods Applied to Chemical Engineering Lecture 19 Prof. William Green Page 4 of 4 Cite as: William Green, Jr., course materials for 10.34 Numerical Methods Applied to Chemical Engineering, Fall 2006. MIT OpenCourseWare (http://ocw.mit.edu), Massachusetts Institute of Technology. Downloaded on [DD Month YYYY]. λλϕϕ2<Δ−=∂∂tt numerical stability of explicit solvers 2Pe')'(local<2<Δ⇒2<Δ=+−=∂∂DvxxDvstuffDvxxxxλλϕϕ Achieve Pelocal < 2 by making Δx smaller. This leads to stiffness. Difficulties when implicit solving. Use adaptive meshing with Gear predictor-corrector. Deff = Dtrue + )()()(200xOxxfxxfxfxvxDxΔ+Δ−Δ+=∂∂Δ 2Pelocal< Method of lines for flow only in x-direction .....)(),,(),,(2),,(.....2222+⎟⎟⎠⎞⎜⎜⎝⎛ΔΔ−+−Δ++∂∂+∂∂+∇+∂∂yzyyxzyxzyyxDxDxvDxvxxϕϕϕϕϕϕϕ PDE Æ ODE {works only if: vy, vz ~ negligible; Pelocal< 2} Equations like this for all discrete z and y values in
View Full Document
Unlocking...