Unformatted text preview:

Lecture 22Advection EquationPeriodic DataDiscrete Fourier TransformTransformCode for the DFTCode for Inverse DFTFast Fourier TransformSpectral DerivativeDetailReal Component of N/2 ModeModified Derivative FormulaSpectral Differentiation SchemecontFinal TwistSpectral Differentiation CodeTwo-Dimensional Fourier TransformSlide 18Runge-Kutta Time IntegratorSide note: Jameson-Schmidt-Turkel Runge-Kutta IntegratorJST Runge-KuttaFactorized SchemeJST + Advection EquationNow Use Spectral RepresentationMatlab ImplementationSlide 26Slide 27FFT ResourcesLab ExerciseLecture 22MA471 Fall 2003Advection Equation•Recall the 2D advection equation:•We will use a Runge-Kutta time integrator and spectral representation in space.0x yC C Ca at x y� � �+ + =� � �Periodic Data•Let’s assume we are given N values of a function f at N data points on the unit interval.•For instance N=8 on a unit interval:0 1/8 2/8 3/8 4/8 5/8 6/8 7/8 8/8Discrete Fourier Transform•We can seek a trigonometric interpolation of a function f as a linear combination of N (even) trigonometric functions:•Such that:( )( )/ 22/ 2 1ˆk Nik xkk NIf x f ep==- +=�( ) ( )( )/ 22/ 2 1ˆ for 1,...,where jk Nik xj j kk Njf x If x f e j NjxNp==- += = ==�Transform•The interpolation formula defines a linear system for the unknown fhat coefficients:( )/ 22/ 2 121ˆ for 1,...,1ˆ for k 1,...,2 2jjk Nik xj kk Nj Nik xk jjf f e j NN Nf f eNpp==- +=-== =� �= =- +� �� ���Or:( )( )ˆˆifftfft==f ff fCode for the DFTCode for Inverse DFTFast Fourier Transform•See handoutSpectral Derivative•We can differentiate the interpolant by:( ) ( )( )/ 22/ 2 1ˆ2k Nik xkk NdIfx ik f edxpp==- +=�( )( )/ 22/ 2 1ˆk Ni kxkk NIf x f ep==- +=�Detail•We note that the derivative of the k=(N/2) mode•is technically:•However, as we show on the next slide – this mode has turning points at all the data points.222ˆiNxNe fp� �� �� �222ˆ22iNxNiNe fpp� �� �� �Real Component of N/2 Modei.e. the slope of the k=(N/2) mode is zeroat all the 8 points..Modified Derivative Formula•So we can create an alternative symmetric derivative formula:( ) ( )( )( )/ 2 12/ 2 11ˆ22k Nik xkk NdIfx ik f edxppp= -=- +@�Spectral Differentiation Scheme1) Use an fft to compute:2) Differentiate in spectral space. i.e. compute:ˆ for 1, 2,..., 1,0,1,...,2 2 2kN N Nf k� � � � � �=- + - + -� � � � � �� � � � � �( )ˆ2 for 1,..., 12 2ˆ0 for 2kkN Ni k f kgNkp�� � � �=- + -� � � ���� � � �=��=��cont3) Then:4) Summary:a) fft transform data to compute coefficients b) scale Fourier coefficients c) inverse fft (ifft) scaled coefficients( )( )( )/ 2/ 2 1ˆk Nikxkk NdIfx g edx==- +@�Final Twist•Matlab stores the coefficients from the fast Fourier transform in a slightly odd order:•The derivative matrix will now be a matrix with diagonal entries:0 1 11 1 22 2 22ˆ ˆ ˆ ˆ ˆ ˆ, ,.., , , , ,...,ˆNN N Nf f ff f f f-- - + - +2 0,1,.., 1, , 1, 2,..., 12 2 2N N Ni p� �- - + - + -� �� �0Spectral Differentiation Code1) DFT data2) Scale Fourier coefficients3) IFT scaled coefficientsTwo-Dimensional Fourier Transform•We can now construct a Fourier expansion in two variables for a function of two variables..• The 2D inverse discrete transform and discrete transform are:( )( )/ 2/ 22 2/ 2 1 / 2 12 221 1ˆ for , 1,...,1ˆ for j,k 1,...,2 2n mn mj Nk Nij x ik ynm jkk N j Nn N m Nij x ik yjk nmn mf f e e n m NN Nf f e eNp pp p===- + =- += =- -= == =� �= =- +� �� �� ���( )( )( )/ 2/ 22 2/ 2 1 / 2 1ˆ,j Nk Nij x ik yjkk N j NIf x y f e ep p===- + =- +=� �Advection Equation•Recall the 2D advection equation:•We will use a Runge-Kutta time integrator and spectral representation in space.0x yC C Ca at x y� � �+ + =� � �Runge-Kutta Time Integrator•We will now discuss a particularly simple Runge-Kutta time integrator introduced by Jameson-Schmidt-Turkel•The idea is each time step is divided into s substeps, which taken together approximate the update to s’th order.Side note: Jameson-Schmidt-Turkel Runge-Kutta Integrator•Taylor’s theorem tell’s us that•We will compute an approximate update as:( ) ( )( )( )[ ]221 1* *11 ...1! 2! ! for some ,1 !sss ssdt d dt d dt dC t dt C tdt dt s dtdt d Ct t t t dts dt+ ++� �� � � � � �+ = + + + +� �� � � � � �� ����� �� �+ � +� �+�( ) ( )221 ...1! 2! !ssdt d dt d dt dC t dt C tdt dt s dt� �� � � � � �+ @ + + + +� �� � � � � �� ����� �JST Runge-Kutta•The numerical scheme will look like•We then factorize the polynomial derivative term:( ) ( )222211 ...1! 2! !1 ...1! 2! !ssssn ndt d dt d dt dC t dt C tdt dt s dtdt d dt d dt dC Cdt dt s dt+� �� � � � � �+ @ + + + +� �� � � � � �� ����� �� �� � � � � �� = + + + +� �� � � � � �� ����� �Factorized Schemen+1..1 2 1nn n n ndt d dt d dt d dt dCC C C C Cs dt s dt s dt dt������ ��= + + + +���� �����- -� ������2211 ...1! 2! !ssn ndt d dt d dt dC Cdt dt s dt+� �� � � � � �= + + + +� �� � � � � �� ����� �nnn+1Set C=Cfor : 1:1 C CendC =Ck sdt dCk dt= -� +JST + Advection Equation•We now use the PDE definition nnn+1Set C=Cfor : 1:1 C CendC =Ck sdt dCk dt= -� ++nnn+1Set C=Cfor : 1:1 C CendC =Cx yk sdt C Ca ak x y= -� �� �� - +� �� �� �0x yC C Ca at x y� � �+ + =� � �Now Use Spectral Representation( )( )( )( )j jSet for : 1:1ˆ = fftˆˆ ˆˆ ˆ ifftˆ ifft for j=1,..,Nend=x xy yx xy yx x y yj j j jk sdtk== -====� - +c cc cd D cd D cd dd dc c a d a dc c%%%%A time step now consists of s substages.Each stage involves:a) Fourier transforming the ctildeusing a fast Fourier transform (fft)b) Scaling the coefficients to differentiate in Fourier spacec) Transforming the derivatives back tophysical values at the nodes by inverse fast Fourier transform (ifft).d) Finally updating ctilde according tothe advection equation.At the end we update the concentration.Matlab Implementation•First set up the nodes and the differentiation scalings.24) Compute dt26) Initiate concentration28-29) compute advection vector32-45)


View Full Document

Rice CAAM 471 - Lecture Notes

Download Lecture Notes
Our administrator received your request to download this document. We will send you the file to your email shortly.
Loading Unlocking...
Login

Join to view Lecture Notes and access 3M+ class-specific study document.

or
We will never post anything without your permission.
Don't have an account?
Sign Up

Join to view Lecture Notes 2 2 and access 3M+ class-specific study document.

or

By creating an account you agree to our Privacy Policy and Terms Of Use

Already a member?