260 likes | 386 Views
MA2213 Lecture 11. PDE. Topics. Introduction p. 451-452. Poisson equation p. 453-466. Boundary conditions p. 453. Finite difference grid p. 454-456. MATLAB Program p. 456-458. Visualization of numerical results p. 459-466. One-dimensional heat equation p. 466-481. Introduction.
E N D
Topics Introduction p. 451-452 Poisson equation p. 453-466 Boundary conditions p. 453 Finite difference grid p. 454-456 MATLAB Program p. 456-458 Visualization of numerical results p. 459-466 One-dimensional heat equation p. 466-481
Introduction “Many phenomena in sciences and engineering depend on more than one variable. For example, an unknown function of a real-world problem usually depends on both time t and the location of the point (x,y,z).” p. 451 Physical laws, including the conservation of energy, momentum and mass, the laws of electricity and magnetism, thermodynamics, and chemical kinetics, require that the partial derivatives of these functions satisfy certain (partial differential) equations.
Introduction Increasingly, PDE’s are used to model biological and social phenomena. The models include the “law of supply and demand” in economics that determines equilibrium prices of goods and services, the Black-Sholes equation for options prices in arbitrage-free financial markets, and laws that describe the evolution of population densities that are used in epidemiology, ecology, and population genetics. http://www.imbs.uci.edu/index.html
Introduction Examples Poisson equation Heat equation Wave equation
Poisson Equation Boundary Conditions Let be a planar domain, and denote its boundary by The boundary value problem is called a Dirichlet problem because the value of u is specified on the boundary. For simplicity we will assume that and therefore
Finite Difference Grid For the square domain we choose a grid of points example boundary points interior points
Finite Difference Approximation We approximate the PDE by where at interior points For boundary points
Finite Difference Approximation We approximate by the solutions of Where, at every boundary point, we define
Finite Difference Equations We now choose to obtain equations for and equations for
Solution Using Gauss-Seidel for k = 1,…,maxit for i = 2,…,n for j = 2,…,n end % j loop end % i loop end % k loop Compare this with the Jacobi method that you used for problems in slides 37,41 of Lecture 7.
MATLAB Program n = 10; maxit = 60;h = 1/n; [X,Y] = meshgrid(0:h:1,0:h:1); u = zeros(n+1); f = -2*pi^2*sin(pi*X).*sin(pi*Y); figure(1); surf(X,Y,f); title([‘f = -source distribution’]) for k = 1:maxit for i = 2:n for j = 2:n u(i,j) = 0.25*(u(i-1,j)+u(i+1,j)+u(i,j+1)+u(i,j-1)-h^2*f(i,j)); end end end figure(2); surf(X,Y,u); title([‘computed solution’]) utrue = sin(pi*X).*sin(pi*Y); figure(3); surf(X,Y,utrue-u); title([‘error with ’ num2str(maxit) ‘ iterations’])
One Dimensional Heat Equation We consider the initial value problem for a function The heat equation arises in the modeling of a number of phenomena and is often used in financial mathematics in the modeling of options. The famous Black-Scholes option pricing model's differential equation can be transformed into the heat equation allowing relatively easy solutions.
One Dimensional Heat Equation The domain of the function boundary value boundary value initial value
Discretization We divide the interval [0,L] into equal parts to obtain a grid with size and (space) grid points and divide We also choose maximum time the interval [0,T] into equal parts to obtain a grid with size and (time) grid points This provides a two-dimensional grid with points
Discretization We define and will compute approximations We first use approximations
Discretization and combine them with the Heat PDE to obtain approximate equations We replace by
Discretization to obtain linear equations that can be solved directly by
Matrix Representation Define sequences Then where and the matrix in slide 38, L7 has eig.val. the eigenvalues satisfy
Stability The numerical roundoff error sequence satisfies therefore the algorithm is stable if and only if the iteration matrix has spectral radius < 1 Since this occurs for large when
MATLAB CODE L = 1; nx = 20; hx = L/nx; a = 1; ht = hx^2/3; % for stable % ht = hx^2; for unstable niter = 500; F = zeros(nx-1,1); F(nx-1) = ht/(hx^2); u(:,1) = zeros(nx-1,1); A = -diag(2*ones(nx-1,1)) + diag(ones(nx-2,1),1) + diag(ones(nx-2,1),-1); I = eye(nx-1); c = a*ht/(hx^2); M = I + c*A; for k = 1:niter u(:,k+1) = M*u(:,k) + F; end