four_point_upwind_uni_d1.m

有关the method of line即MOL算法的matlab代码

      function [D]=four_point_upwind_uni_D1(z0,zL,n,v)
%...
%...  A. Vande Wouwer, P. Saucez and W.E. Schiesser (2002)
%...
%...  function four_point_upwind_uni_D1 returns the differentiation matrix 
%...  for computing the first derivative, xz, of a variable x over the spatial 
%...  domain z0 < z < zL from upwind four-point, first-order finite difference 
%...  approximations (this function replaces dss016)
%...
%...  argument list
%...
%...     z0      lower boundary value of z (input)
%...
%...     zL      upper boundary value of z (input)
%...
%...     n       number of grid points in the z domain including the
%...             boundary points (input)
%...
%...     v       fluid velocity (positive from left to right - only the sign is used) (input)
%...
%...
%...  origin of the approximation
%...
%...  this function is an application of fourth-order directional
%...  differencing in the numerical method of lines.  It is intended
%...  specifically for the analysis of convective systems modelled by
%...  first-order hyperbolic partial differential equations as dis-
%...  cussed in subroutine dss012. The coefficients of the finite
%...  difference approximations used herein are taken from Bickley, W.
%...  G., Formulae for numerical differentiation, The Mathematical
%...  Gazette, pp. 19-27, 1941, n = 3, m = 1, p = 0, 1, 2, 3.
%...
%...
%...  compute the spatial increment
      dz=(zL-z0)/(n-1);
      r3fdz=1/(6*dz);
%...
%...     (1)  finite difference approximation for positive v     
              if v > 0
%...
%...             sparse discretization matrix      
%...
%...             interior points      
                 D=diag(-2*ones(n-3,1),-3)+diag(+9*ones(n-2,1),-2)+diag(-18*ones(n-1,1),-1)+diag(+11*ones(n,1),0);      
%...
%...             boundary points      
                 D(1,1:4) = [-11 +18 -9 +2];
                 D(2,1:4) = [-2  -3  +6 -1];
                 D(3,1:4) = [+1  -6  +3 +2];
              end;
%...
%...     (2)  finite difference approximation for negative v
              if v < 0
%...
%...             sparse discretization matrix      
%...
%...             interior points      
                 D=diag(-11*ones(n,1),0)+diag(+18*ones(n-1,1),+1)+diag(-9*ones(n-2,1),+2)+diag(+2*ones(n-2,1),+3);      
%...
%...             boundary points      
                 D(n-2,(n-3):n) = [-2 -3 +6 -1];
                 D(n-1,(n-3):n) = [+1 -6 +3 +2];
                 D(n,(n-3):n) = [-2 +9 -18 +11];
              end;                
%...      
      D=r3fdz*D;
      D=sparse(D);