dss020=function(xl,xu,n,u,v) { # # An extensive set of documentation comments detailing the derivation # of the following five point pwind finite differences (FDs) is not given # here to conserve space. The derivation is also detailed in Schiesser, # W. E., The Numerical Method of Lines Integration of Partial Differential # Equations, Academic Press, San Diego, 1991. # # Declare arrays ux=rep(0,n); # # Grid spacing # # 1/(12*dx) for subsequent use dx=(xu-xl)/(n-1); r12dx=1/(12*dx); # # (1) Finite difference appoximation for positive v # if(v>0){ ux[ 1]=r12dx*(-25*u[ 1]+48*u[ 2]-36*u[ 3]+16*u[ 4] -3*u[ 5]); ux[ 2]=r12dx*( -3*u[ 1]-10*u[ 2]+18*u[ 3] -6*u[ 4] +u[ 5]); ux[ 3]=r12dx*( u[ 1] -8*u[ 2]+ +8*u[ 4] -u[ 5]); for(i in 4:(n-1)){ ux[ i]=r12dx*( -u[i-3] +6*u[i-2]-18*u[i-1]+10*u[i ] +3*u[i+1]);} ux[ n]=r12dx*( 3*u[n-4]-16*u[n-3]+36*u[n-2]-48*u[n-1]+25*u[n ]); } # # (2) Finite difference approximation for negative v # if(v<0){ ux[ 1]=r12dx*( -25*u[ 1]+48*u[ 2]-36*u[ 3]+16*u[ 4] -3*u[ 5]); for(i in 2:(n-3)){ ux[ i]=r12dx*( -3*u[i-1]-10*u[i ]+18*u[i+1] -6*u[i+2] +u[i+3]);} ux[n-2]=r12dx*( u[n-4] -8*u[n-3] +8*u[n-1] -u[n ]); ux[n-1]=r12dx*( -u[n-4] +6*u[n-3]-18*u[n-2]+10*u[n-1] +3*u[n ]); ux[ n]=r12dx*( 3*u[n-4]-16*u[n-3]+36*u[n-2]-48*u[n-1]+25*u[n ]); } # # All points concluded (x=xl,...,x=xu) return(c(ux)); }