Download Presentation
## MA557/MA578/CS557 Lecture 19

- - - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - - -

**MA557/MA578/CS557Lecture 19**Spring 2003 Prof. Tim Warburton timwar@math.unm.edu**DG Scheme for the Scalar ADE(Lax-Friedrichs flux –**equivalent to upwind)**DG Derivative Operator**• We are going to introduce a DG derivative operator to simplify the scheme definition. • The linear operator Dtilde is such that the following holds for all intervals Ij • With the choice of penalty terms tauL and tauR to be determined.**In Operator Notation**Advectionterm(Lax-Friedrichs flux ~ upwind fluxfor scalar case) Diffusionterm**General Dtilde Operator**• We previously defined Dtilde for the special case where the test function is compactly supported on the j’th cell. • We generalize by:**Skew-symmetry of Dtilde1,1**[integrate by parts]**Skew-symmetry**• Note:the negative sign when the operator moves. • Trivially:**Recap**• We already proved that the DG scheme, with Lax-Friedrichs fluxes, for advection (D=0) is stable in the sense that : • Adding the diffusion term (dropping boundary terms) gives: • In terms of norms:**Stability!!**• Assuming D>=0 (reasonable since particles do not jump randomly with negative rates)**D operator**• If we move to the Legendre basis we can write down a discrete description of Dtilde:**Implementation of Dtilde**function [dfdx] = LEGdgderiv(f, nodex, tauL, tauR, D, F, G, H, J) [p,N] = size(f); % p = maximum polynomial order p=p-1; % N = number of nodes N=N+1; % space for derivative dfdx = zeros(p+1,N-1); % cells 2 to N-2 ids = (2:N-2); dfdx(:,ids) = (D*f(:,ids) + tauL*F*f(:,ids) + tauL*G*f(:,ids-1) + tauR*H*f(:,ids) + tauR*J*f(:,ids+1)); % cell 1 (periodicity assumed) dfdx(:,1) = (D*f(:,1) + tauL*F*f(:,1) + tauL*G*f(:,N-1) + tauR*H*f(:,1) + tauR*J*f(:,2)); % cell N-1 (periodicity assumed) dfdx(:,N-1) = (D*f(:,N-1) + tauL*F*f(:,N-1) + tauL*G*f(:,N-2) + tauR*H*f(:,N-1) + tauR*J*f(:,1)); % apply chain rule for physical cell width dx = nodex(2:N)-nodex(1:N-1); coeff = ones(p+1,1)*(2./dx); dfdx = dfdx.*coeff;**In Action (advection eqn)**% TIME STEPPING tauL = (1/u)*(u+abs(u))/2; tauR = (1/u)*(u-abs(u))/2; for tstep = 1:Ntsteps sigma = rho; for rkstage = 3:-1:1 dsigmadx = LEGdgderiv(sigma, nodex, tauL, tauR, D, F, G, H, J); sigma = rho + (dt/rkstage)*(-u*dsigmadx); end rho = sigma; if ( ~mod(tstep, 100) ) q = V*rho; % plot(x(:),abs(q(:)-exp(-(x(:)-dt*tstep).^2))); plot(x,q); pause(0.1); tstep*dt end end**Discrete Scheme**• Step 1: q = x DG derivative of C • Step 2: time rate of change of C Advectionterm(Lax-Friedrichs flux ~ upwind fluxfor scalar case) Diffusionterm**Dropping q**• We can do away with q: • This should make the following observation clear. The dt dependence for stability is now: • The first term is due to the spectral radius of the advection operator. • The second term is due to the spectral radius of the diffusion operator.**In Action (ADE)**% TIME STEPPING for tstep = 1:Ntsteps sigma = rho; for rkstage = p:-1:1 % advection: upwind bias on flux terms dsigmadx = LEGdgderiv(sigma, nodex, tauL, tauR, D, F, G, H, J); % diffusion: centered fluxes dsdx = LEGdgderiv(sigma, nodex, 1, 1, D, F, G, H, J); d2sdx2 = LEGdgderiv(dsdx, nodex, 1, 1, D, F, G, H, J); sigma = rho + (dt/rkstage)*(-u*dsigmadx+Dcoeff*d2sdx2); end rho = sigma; if ( ~mod(tstep, 400) ) q = V*rho; plot(x,q); axis([xmin xmax 0 1]) pause(0.1); end end**Downloads**• You can download the Dtilde based scripts from: • http://www.math.unm.edu/~timwar/MA578S03/MatlabScripts/LEGadvectionLaxF.m • http://www.math.unm.edu/~timwar/MA578S03/MatlabScripts/LEGade.m • You will need the Dtilde routine (named LEGdgderiv.m): • http://www.math.unm.edu/~timwar/MA578S03/MatlabScripts/LEGdgderiv.m**Legendre Scheme To Lagrange Scheme(modal to nodal)**• So far we have used Legendre expansions to represent the solution in each cell. • This representation causes some problems if we wish to evaluate non-linear function of the solution. • For instance solving the inviscid Burger’s equation:requires the evaluation of u2**Transforming the Modal Advection DG Scheme To a Nodal**Advection DG Scheme • Recall that given a Legendre expansion we may evaluate it at a set of points by: • We substitute: into the scheme:**Nodal Scheme**• Multiplying both sides on the left with V • This clearly shows that the nodal scheme is a change of basis away from the modal scheme.**Summary of Nodal Scheme**Note: we do not need to transform the initial condition or solution for visualization…**Next Lecture**• The lecture on Friday 03/07 is cancelled in favour of the Analysis Conference • The lecture on Monday 03/10 will be the first time we generalize the DG to two spatial dimensions. • Do not miss the Monday lecture • Run the LEGade.m script (you will need the LEGdgderiv.m script as well) • How to use nodal methods in general.