# Method of lines/example implementation/dss006

% File: dss006.m % function [ux]=dss006(xl,xu,n,u) % % Function dss006 computes the first derivative, u , of a % x % Variable u over the spatial domain xl le x le xu from classical % seven-point, sixth-order finite difference approximations % % Argument list % % xl Lower boundary value of x (input) % % xu Upper boundary value of x (input) % % n Number of grid points in the x domain including the % boundary points (input) % % u One-dimensional array containing the values of u at % the n grid point points for which the derivative is % to be computed (input) % % ux One-dimensional array containing the numerical % values of the derivatives of u at the n grid points % (output) % % The mathematical details of the following Taylor series (or % polynomials) are given in routines dss002 and dss004. % % Seven-point formulas % % (1) Left end, point i = 1 % % 2 3 4 % a(u2 = u1 + u1 ( dx) + u1 ( dx) + u1 ( dx) + u1 ( dx) % x 1f 2x 2f 3x 3f 4x 4f % % 5 6 7 % + u1 ( dx) + u1 ( dx) + u1 ( dx) + ...) % 5x 5f 6x 6f 7x 7f % % 2 3 4 % b(u3 = u1 + u1 (2dx) + u1 (2dx) + u1 (2dx) + u1 (2dx) % x 1f 2x 2f 3x 3f 4x 4f % % 5 6 7 % + u1 (2dx) + u1 (2dx) + u1 (2dx) + ...) % 5x 5f 6x 6f 7x 7f % % 2 3 4 % c(u4 = u1 + u1 (3dx) + u1 (3dx) + u1 (3dx) + u1 (3dx) % x 1f 2x 2f 3x 3f 4x 4f % % 5 6 7 % + u1 (3dx) + u1 (3dx) + u1 (3dx) + ...) % 5x 5f 6x 6f 7x 7f % % 2 3 4 % d(u5 = u1 + u1 (4dx) + u1 (4dx) + u1 (4dx) + u1 (4dx) % x 1f 2x 2f 3x 3f 4x 4f % % 5 6 7 % + u1 (4dx) + u1 (4dx) + u1 (4dx) + ...) % 5x 5f 6x 6f 7x 7f % % 2 3 4 % e(u6 = u1 + u1 (5dx) + u1 (5dx) + u1 (5dx) + u1 (5dx) % x 1f 2x 2f 3x 3f 4x 4f % % 5 6 7 % + u1 (5dx) + u1 (5dx) + u1 (5dx) + ...) % 5x 5f 6x 6f 7x 7f % % 2 3 4 % f(u7 = u1 + u1 (6dx) + u1 (6dx) + u1 (6dx) + u1 (6dx) % x 1f 2x 2f 3x 3f 4x 4f % % 5 6 7 % + u1 (6dx) + u1 (6dx) + u1 (6dx) + ...) % 5x 5f 6x 6f 7x 7f % % Constants a, b, c, d, e and f are selected so that the coeffi- % cients of the u1 terms sum to one and the coefficients of % x % the u1 , u1 , u1 , u1 and u1 terms sum to zero. % 2x 3x 4x 5x 6x % % 1 1 1 1 1 % a + 2 b + 3 c + 4 d + 5 e + 6 f = 1 % % 2 2 2 2 2 % a + 2 b + 3 c + 4 d + 5 e + 6 f = 0 % % 3 3 3 3 3 % a + 2 b + 3 c + 4 d + 5 e + 6 f = 0 % % 4 4 4 4 4 % a + 2 b + 3 c + 4 d + 5 e + 6 f = 0 % % 5 5 5 5 5 % a + 2 b + 3 c + 4 d + 5 e + 6 f = 0 % % 6 6 6 6 6 % a + 2 b + 3 c + 4 d + 5 e + 6 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for u1 gives the following seven-point approxi- % 6x x % tion % % u1 = (1/6f)(-1764*u1 + 4320u2 - 5400u3 + 4800u4 % x 6 (1) % + 2700u5 + 864u6 - 120u7) + O(dx ) % % (2) Interior point, i = 2 % % The preceding Taylor series can be summarized for point i = 2 % as % % a(u1 = u2 + ...) % % b(u3 = u2 + ...) % % c(u4 = u2 + ...) % % d(u5 = u2 + ...) % % e(u6 = u2 + ...) % % f(u7 = u2 + ...) % % The corresponding algebraic equations are % % 1 1 1 1 1 1 % -1 a + 1 b + 2 c + 3 d + 4 e + 5 f = 1 % % 2 2 2 2 2 2 % -1 a + 1 b + 2 c + 3 d + 4 e + 5 f = 0 % % 3 3 3 3 3 3 % -1 a + 1 b + 2 c + 3 d + 4 e + 5 f = 0 % % 4 4 4 4 4 4 % -1 a + 1 b + 2 c + 3 d + 4 e + 5 f = 0 % % 5 5 5 5 5 5 % -1 a + 1 b + 2 c + 3 d + 4 e + 5 f = 0 % % 6 6 6 6 6 6 % -1 a + 1 b + 2 c + 3 d + 4 e + 5 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for u2 gives the following seven-point approxi- % 6x x % tion % % u2 = (1/6f)(-120u1 - 924u2 + 1800u3 - 1200u4 % x 6 (2) % + 600u5 - 180u6 + 24u7) + O(dx ) % % (3) Interior point, i = 3 % % The preceding Taylor series can be summarized for point i = 3 % as % % a(u1 = u3 + ...) % % b(u2 = u3 + ...) % % c(u4 = u3 + ...) % % d(u5 = u3 + ...) % % e(u6 = u3 + ...) % % f(u7 = u3 + ...) % % The corresponding algebraic equations are % % 1 1 1 1 1 1 % -2 a + -1 b + 1 c + 2 d + 3 e + 4 f = 1 % % 2 2 2 2 2 2 % -2 a + -1 b + 1 c + 2 d + 3 e + 4 f = 0 % % 3 3 3 3 3 3 % -2 a + -1 b + 1 c + 2 d + 3 e + 4 f = 0 % % 4 4 4 4 4 4 % -2 a + -1 b + 1 c + 2 d + 3 e + 4 f = 0 % % 5 5 5 5 5 5 % -2 a + -1 b + 1 c + 2 d + 3 e + 4 f = 0 % % 6 6 6 6 6 6 % -2 a + -1 b + 1 c + 2 d + 3 e + 4 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for u3 gives the following seven-point approxi- % 6x x % tion % % u3 = (1/6f)(24u1 - 288u2 - 420u3 + 960u4 % x 6 (3) % - 360u5 + 96u6 - 12u7) + O(dx ) % % (4) Interior point, i ne 2, 3, n-2, n-1 % % The preceding Taylor series can be summarized for point i = i % as % % a(ui-3 = ui + ...) % % b(ui-2 = ui + ...) % % c(ui-1 = ui + ...) % % d(ui+1 = ui + ...) % % e(ui+2 = ui + ...) % % f(ui+3 = ui + ...) % % The corresponding algebraic equations are % % 1 1 1 1 1 1 % -3 a + -2 b + -1 c + 1 d + 2 e + 3 f = 1 % % 2 2 2 2 2 2 % -3 a + -2 b + -1 c + 1 d + 2 e + 3 f = 0 % % 3 3 3 3 3 3 % -3 a + -2 b + -1 c + 1 d + 2 e + 3 f = 0 % % 4 4 4 4 4 4 % -3 a + -2 b + -1 c + 1 d + 2 e + 3 f = 0 % % 5 5 5 5 5 5 % -3 a + -2 b + -1 c + 1 d + 2 e + 3 f = 0 % % 6 6 6 6 6 6 % -3 a + -2 b + -1 c + 1 d + 2 e + 3 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for ui gives the following seven-point approxi- % 6x x % tion % % ui = (1/6f)(-12ui-3 + 108ui-2 - 540ui-1 + 0ui % x 6 (4) % + 540ui+1 - 108ui+2 + 12ui+3) + O(dx ) % % (5) interior point, i = n-2 % % The preceding Taylor series can be summarized for point i = n-2 % as % % a(un-6 = un-2 + ...) % % b(un-5 = un-2 + ...) % % c(un-4 = un-2 + ...) % % d(un-3 = un-2 + ...) % % e(un-1 = un-2 + ...) % % f(un = un-2 + ...) % % The corresponding algebraic equations are % % 1 1 1 1 1 1 % -4 a + -3 b + -2 c + -1 d + 1 e + 2 f = 1 % % 2 2 2 2 2 2 % -4 a + -3 b + -2 c + -1 d + 1 e + 2 f = 0 % % 3 3 3 3 3 3 % -4 a + -3 b + -2 c + -1 d + 1 e + 2 f = 0 % % 4 4 4 4 4 4 % -4 a + -3 b + -2 c + -1 d + 1 e + 2 f = 0 % % 5 5 5 5 5 5 % -4 a + -3 b + -2 c + -1 d + 1 e + 2 f = 0 % % 6 6 6 6 6 6 % -4 a + -3 b + -2 c + -1 d + 1 e + 2 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for un-2 gives the following seven-point approxi- % 6x x % tion % % un-2 = (1/6f)(12un-6 - 96un-5 + 360un-4 - 960un-3 % x 6 (5) % + 420un-2 + 288un-1 - 24un) + O(dx ) % % (6) Interior point, i = n-1 % % The preceding Taylor series can be summarized for point i = n-1 % as % % a(un-6 = un-1 + ...) % % b(un-5 = un-1 + ...) % % c(un-4 = un-1 + ...) % % d(un-3 = un-1 + ...) % % e(un-2 = un-1 + ...) % % f(un = un-1 + ...) % % The corresponding algebraic equations are % % 1 1 1 1 1 1 % -5 a + -4 b + -3 c + -2 d + -1 e + 1 f = 1 % % 2 2 2 2 2 2 % -5 a + -4 b + -3 c + -2 d + -1 e + 1 f = 0 % % 3 3 3 3 3 3 % -5 a + -4 b + -3 c + -2 d + -1 e + 1 f = 0 % % 4 4 4 4 4 4 % -5 a + -4 b + -3 c + -2 d + -1 e + 1 f = 0 % % 5 5 5 5 5 5 % -5 a + -4 b + -3 c + -2 d + -1 e + 1 f = 0 % % 6 6 6 6 6 6 % -5 a + -4 b + -3 c + -2 d + -1 e + 1 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for un-1 gives the following seven-point approxi- % 6x x % tion % % un-1 = (1/6f)(-24un-6 + 180un-5 - 600un-4 + 1200un-3 % x 6 (6) % - 1800un-2 + 924un-1 + 120un) + O(dx ) % % (7) Right end, point i = n % % The preceding Taylor series can be summarized for point i = n % as % % a(un-6 = un + ...) % % b(un-5 = un + ...) % % c(un-4 = un + ...) % % d(un-3 = un + ...) % % e(un-2 = un + ...) % % f(un-1 = un + ...) % % The corresponding algebraic equations are % % 1 1 1 1 1 1 % -6 a + -5 b + -4 c + -3 d + -2 e + -1 f = 1 % % 2 2 2 2 2 2 % -6 a + -5 b + -4 c + -3 d + -2 e + -1 f = 0 % % 3 3 3 3 3 3 % -6 a + -5 b + -4 c + -3 d + -2 e + -1 f = 0 % % 4 4 4 4 4 4 % -6 a + -5 b + -4 c + -3 d + -2 e + -1 f = 0 % % 5 5 5 5 5 5 % -6 a + -5 b + -4 c + -3 d + -2 e + -1 f = 0 % % 6 6 6 6 6 6 % -6 a + -5 b + -4 c + -3 d + -2 e + -1 f = 0 % % Simultaneous solution for a, b, c, d, e and f followed by the % solution of the preceding Taylor series, truncated after the % u terms, for un gives the following seven-point approxi- % 6x x % tion % % un = (1/6f)(120un-6 - 864un-5 + 2700un-4 - 4800un-3 % x 6 (7) % + 5400un-2 - 4320un-1 + 1764un) + O(dx ) % % The weighting coefficients for equations (1) to (7) can be % summarized as % % -1764 4320 -5400 4800 -2700 864 -120 % % -120 -924 1800 -1200 600 -180 24 % % 24 -288 -420 960 -360 96 -12 % % 1/6f -12 108 -540 0 540 -108 12 % % 12 -96 360 -960 420 288 -24 % % -24 180 -600 1200 -1800 924 120 % % 120 -864 2700 -4800 5400 -4320 1764 % % which are the coefficients reported by Bickley for n = 6, m = % 1, p = 0 to 6 (Bickley, W. G., Formulae for Numerical Differ- % entiation, Math. Gaz., vol. 25, 1941). % % Equations (1) to (7) can now be programmed to generate the % derivative u (x) of function u(x) (arguments u and ux of % x % function dss006, respectively). % % Compute the spatial increment dx=(xu-xl)/(n-1); r6fdx=1./(720.*dx); nm3=n-3; % % Equation (1) ux( 1)=r6fdx*... ( -1764.*u( 1) +4320.*u( 2) -5400.*u( 3) +4800.*u( 4)... -2700.*u( 5) +864.*u( 6) -120.*u( 7)); % % Equation (2) ux( 2)=r6fdx*... ( -120.*u( 1) -924.*u( 2) +1800.*u( 3) -1200.*u( 4)... +600.*u( 5) -180.*u( 6) +24.*u( 7)); % % Equation (3) ux( 3)=r6fdx*... ( +24.*u( 1) -288.*u( 2) -420.*u( 3) +960.*u( 4)... -360.*u( 5) +96.*u( 6) -12.*u( 7)); % % Equation (4) for i=4:nm3 ux( i)=r6fdx*... ( -12.*u(i-3) +108.*u(i-2) -540.*u(i-1) +0.*u( i)... +540.*u(i+1) -108.*u(i+2) +12.*u(i+3)); end % % Equation (5) ux(n-2)=r6fdx*... ( +12.*u(n-6) -96.*u(n-5) +360.*u(n-4) -960.*u(n-3)... +420.*u(n-2) +288.*u(n-1) -24.*u( n)); % % Equation (6) ux(n-1)=r6fdx*... ( -24.*u(n-6) +180.*u(n-5) -600.*u(n-4) +1200.*u(n-3)... -1800.*u(n-2) +924.*u(n-1) +120.*u( n)); % % Equation (7) ux( n)=r6fdx*... ( +120.*u(n-6) -864.*u(n-5) +2700.*u(n-4) -4800.*u(n-3)... +5400.*u(n-2) -4320.*u(n-1) +1764.*u( n)); %