{VERSION 6 0 "IBM INTEL NT" "6.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 1 }{CSTYLE "2D Comment" 2 18 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 } {PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 483 " This Maple program is used in MATH 36 5 - Computational\n Fluid Dynamics at James Madison University. This program was\n written by James Sochacki of the Department of Mathem atics\n at James Madison University as part of the NSF Grant - A Col laborative \n Computational Science Program between JMU and North Ca rolina Central\n University. His info is\n \n Jim Sochacki\n D epartment of Mathematics\n James Madison University\n Harrisonburg , VA 22807\n jim@math.jmu.edu\n" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "with(plots):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "plotsetup(gdi):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 82 "This spread sheet gives the nume rical solution to the advection diffusion equation" }}{PARA 0 "" 0 "" {TEXT -1 17 " " }{XPPEDIT 18 0 "diff(w,t)+c*diff(w,x)+ d*diff(w,`$`(x,2)) = 0;" "6#/,(-%%diffG6$%\"wG%\"tG\"\"\"*&%\"cGF*-F&6 $F(%\"xGF*F**&%\"dGF*-F&6$F(-%\"$G6$F/\"\"#F*F*\"\"!" }{TEXT -1 0 "" } }{PARA 0 "" 0 "" {TEXT -1 21 " " }{XPPEDIT 18 0 "w (x,0) = q*arctan(s*x)+r;" "6#/-%\"wG6$%\"xG\"\"!,&*&%\"qG\"\"\"-%'arct anG6#*&%\"sGF,F'F,F,F,%\"rGF," }{TEXT -1 12 " for -L " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 19 "Value of w(x,0) at " }{XPPEDIT 18 0 "-infinity;" "6#,$%)infinityG!\"\"" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "a_left := -1.0;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 19 "Value of w(x,0) a t " }{XPPEDIT 18 0 "infinity;" "6#%)infinityG" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "a_right := 1.0;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 23 "Slope factor for w(0,0)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "s := \+ 1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "L := 10.0;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "q := (a_right-a_left)/Pi;" } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "r := (a_left+a_right)/2;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "p := x -> q*arctan(s*x)+r ;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 27 "Plot the initial condition." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 37 "pl ot(p(x),x=-L..L,y=a_left..a_right);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 24 "The Leap Frog scheme for" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }{XPPEDIT 18 0 "Diff(w,t)+c*Diff(w,x)+d*Diff(w,`$`(x,2)) = 0; " "6#/,(-%%DiffG6$%\"wG%\"tG\"\"\"*&%\"cGF*-F&6$F(%\"xGF*F**&%\"dGF*-F &6$F(-%\"$G6$F/\"\"#F*F*\"\"!" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 11 "is given by" }}{PARA 0 "" 0 "" {XPPEDIT 18 0 "w[j]^(n+1) = w[j]^(n-1)-c*Delta*t*(w[j+1]^n-w[j-1]^n)/(Delta*x)+ 2*d*Delta*t*(w[j+1]^n-2*w[j]^n+w[j-1]^n)/(Delta*x^2);" "6#/)&%\"wG6#% \"jG,&%\"nG\"\"\"F+F+,()&F&6#F(,&F*F+F+!\"\"F+*,%\"cGF+%&DeltaGF+%\"tG F+,&)&F&6#,&F(F+F+F+F*F+)&F&6#,&F(F+F+F1F*F1F+*&F4F+%\"xGF+F1F1*.\"\"# F+%\"dGF+F4F+F5F+,()&F&6#,&F(F+F+F+F*F+*&FBF+)&F&6#F(F*F+F1)&F&6#,&F(F +F+F1F*F+F+*&F4F+*$F@FBF+F1F+" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 88 " We initialize jmax - maximum number of x grids, nmax - maximum number of t steps, dt = " }{XPPEDIT 18 0 "Delta*t;" "6#*&%&DeltaG\"\"\"%\"tG F%" }{TEXT -1 7 ", dx = " }{XPPEDIT 18 0 "Delta*x;" "6#*&%&DeltaG\"\" \"%\"xGF%" }{TEXT -1 0 "" }}{PARA 0 "" 0 "" {XPPEDIT 18 0 "x[0];" "6#& %\"xG6#\"\"!" }{TEXT -1 7 " = -L, " }{XPPEDIT 18 0 "x[jmax+1];" "6#&% \"xG6#,&%%jmaxG\"\"\"F(F(" }{TEXT -1 3 "= L" }}}{EXCHG {PARA 0 "" 0 " " {TEXT -1 9 "Time step" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "dt := 0. 25;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 4 "Time" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "T := 10;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "num_time_steps := ceil(T/dt);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 9 " Grid size" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "dx := 0.25;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "x[0] := -L;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "jmax := ceil(2*L/dx);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "x[jmax+1] := L;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 27 "Put in the initial values." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "for j from 0 by 1 to jmax+1 do" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 " " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 20 " x[j] := x[0]+j*dx:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " w0[j] := evalf(p(x[j])):" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 2 " " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 3 "od:" }} }{EXCHG {PARA 0 "" 0 "" {TEXT -1 48 "Plot the initial value for the le ap frog scheme." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "listplot([seq([x[j],w0[j]],j =1..jmax)]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 33 "Generate the numerical solution. " }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 11 "The speed c" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "c := 1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 " c1 := -c*dt/dx;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 15 "The diffusion \+ d" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "d := 0;" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 23 "c2 := -2.*d*dt/(dx*dx);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 71 "Start the calculations. The first time step is calcu lated using upwind." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 " if (c<0) then" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 32 " for j from 1 by \+ 1 to jmax do" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 76 " w1[j] := w0[j]+(c*dt/dx)*(w0[j+1]-w0[j])+c2*(w 0[j+1]-2.*w0[j]+w0[j-1]):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 7 " od:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 6 " else" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 32 " for j from 1 by 1 to jmax do" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 76 " w1[j] := w0[j]+(c*dt/dx)*(w0[j ]-w0[j-1])+c2*(w0[j+1]-2.*w0[j]+w0[j-1]):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 " od:" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 " fi:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 46 "Plot the first time step approximate solution." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 40 "listplot([seq([x[j],w1[j]],j=1..jmax)]);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 39 "Put the values for the end grid points." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 41 "w1[0] := w0[0]: w1[jmax+1] := w0[jmax+1]: " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 50 "Use leap frog to generate the numerical solutions." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "for n from 1 by 1 to num_ time_steps do" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 " for j from 1 by \+ 1 to jmax do" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 69 " w2[j] := w0[j]+c1*(w1[j+1]-w1[j-1])+c2*(w1[j+1]- 2.*w1[j]+w1[j-1])\n" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 5 " od:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 " \+ for j from 1 by 1 to jmax do" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 19 " \+ w0[j] := w1[j]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 19 " w1[j] := w 2[j]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 5 " od:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 51 " gr[n] := listplot([seq([x[j],w0[j]],j=1..jmax)]):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 3 "od:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 41 "Click on the graph to view the animation." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 58 "display([seq(gr[n],n=1..num_time_st eps)],insequence=true);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 " " }}}}{MARK "0 0 0" 8 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }