dss042=function(xl,xu,n,u,ux,nl,nu) { # # An extensive set of documentation comments detailing the derivation # of the following second order finite differences (FDs) is not given # here to conserve space. The comments are available in the Matlab # version of dss042 in http://www.pdecomp.net. The derivation is also # detailed in Schiesser, W. E., The Numerical Method of Lines Integration # of Partial Differential Equations, Academic Press, San Diego, 1991. # # Preallocate arrays ux =rep(0,n); uxx=rep(0,n); # # Grid spacing dx=(xu-xl)/(n-1); # # 1/(dx**2) for subsequent use rdxs=1/(dx^2); # # uxx vector # # Boundaries (x=xl,x=xu) if(nl==1)uxx[1]=rdxs*( 2*u[1]-5*u[ 2]+4*u[ 3]-u[ 4]); if(nu==1)uxx[n]=rdxs*( 2*u[n]-5*u[n-1]+4*u[n-2]-u[n-3]); if(nl==2)uxx[1]=rdxs*(-7*u[1]+8*u[ 2] -u[ 3]-6*dx*ux[1])/2; if(nu==2)uxx[n]=rdxs*(-7*u[n]+8*u[n-1] -u[n-2]+6*dx*ux[n])/2; # # Interior points (x=xl+dx,...,x=xu-dx) for(i in 2:(n-1))uxx[i]=rdxs*(u[i+1]-2*u[i]+u[i-1]); # # All points concluded (x=xl,...,x=xu) return(c(uxx)); }