Changes between Version 7 and Version 8 of u/adebrech/code

Timestamp:

04/14/17 11:35:23 (9 years ago)

Author:

adebrech

Comment:

—

u/adebrech/code

@@ -1 +1 @@
 = [attachment:GodunovSolver.tar Godunov Solver] =
 
-
+This program solves either the inviscid Burgers equation or the linear advection equation with the given boundary or initial conditions. First, input is taken from a file:
+
+{{{#!fortran
+NAMELIST/ProblemData/ starttime,finaltime,frames,xL,xR,cells,alpha,beta,BC,Eq,IC,advectionConst,CFL
+        
+        ! Make sure problem.data exists. If it does, read in. If not, output error and quit.
+        inquire(file='GodunovProblem.data', exist=ex)
+        if (ex) then
+         OPEN(UNIT=30, FILE='GodunovProblem.data', STATUS="OLD")
+         READ(30,NML=ProblemData)                               ! read into variables in namelist
+         CLOSE(30)
+        else
+         print *, 'Missing GodunovProblem.data. Exiting.'
+         stop           ! exit if missing
+        end if
+}}}
+
+The domain is then initialized according to the selected conditions:
+
+{{{#!fortran
+select case (IC)
+ case (Gaussian)
+  do i = 1, cells
+   ProblemRegion%u(i) = alpha*exp(-1.*beta*((i-0.5)*dx+xL)**2)
+  end do
+ case (SquareWave)
+  do i = 1, cells
+  x=(i-0.5)*dx + xL
+   if (x <= 0.3) then
+    ProblemRegion%u(i) = 0
+   else if (x >= 0.7) then
+    ProblemRegion%u(i) = 0
+   else
+    ProblemRegion%u(i) = 1
+   end if
+  end do
+ case (BurgersShock)
+  do i = 1, cells
+  x=(i-1)*dx + xL
+   if (x <= 0.5) then
+    ProblemRegion%u(i) = -0.5
+   else if (x >= 1.) then
+    ProblemRegion%u(i) = 0
+   else
+    ProblemRegion%u(i) = 1
+   end if
+  end do
+end select
+}}}
+
+The boundary conditions are then applied:
+
+{{{#!fortran
+select case (BC)
+ case(periodic) ! Wraps around
+  ProblemRegion%u(0) = ProblemRegion%u(cells)
+  ProblemRegion%u(cells+1) = ProblemRegion%u(1)
+ case(transparent)      ! Wave flows out
+  ProblemRegion%u(0) = ProblemRegion%u(1)
+  ProblemRegion%u(cells+1) = ProblemRegion%u(cells)
+end select
+}}}
+
+After this, the time step is determined using the maximum speed in the domain and the given CFL coefficient.
+
+{{{#!fortran
+! Determine maximum speed in domain
+maxSpeed = 0.
+do i = 0, cells+1
+ if (abs(ProblemRegion%u(i)) > maxSpeed) maxSpeed = abs(ProblemRegion%u(i))     ! (See Toro, Ch. 5.6)
+end do
+
+! Calculate delta t
+dt = CFL*dx/maxSpeed
+
+! If time step would pass the next output frame, recalculate
+if (time + dt > nexttime) dt = nexttime - time
+}}}
+
+The intercell flux is then calculated and the velocity array is updated to the next time.
+
+{{{#!fortran
+do i = 0, cells
+ uL = ProblemRegion%u(i)
+ uR = ProblemRegion%u(i+1)   
+ select case (Eq)
+  case (advection)
+   ProblemRegion%flux(i) = LinearAdvectionFlux(RiemannSoln(uL,uR))
+
+!∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨∨
+function RiemannSoln(uL,uR)
+
+if (uL > uR) then
+ S = 0.5*(uL + uR)
+ if (S <= 0) then
+  RiemannSoln = uR
+ else
+  RiemannSoln = uL
+ end if
+else
+ if (uL >= 0) then
+  RiemannSoln = uL
+ else if (uR <= 0) then
+  RiemannSoln = uR
+ else
+  RiemannSoln = 0
+ end if
+end if
+!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+  case (inviscidBurgers)
+   ProblemRegion%flux(i) = BurgersFlux(RiemannSoln(uL,uR))
+ end select
+end do
+
+! Loop through cells using previously calculated flux
+do i = 1, cells
+ ProblemRegion%u(i) = ProblemRegion%u(i) + (dt/dx)*(ProblemRegion%flux(i-1) - ProblemRegion%flux(i))
+end do
+}}}
+
+Finally, at each requested output frame, the current system state is written to a file:
+
+{{{#!fortran
+write (40,*) 'x u'       ! write column headers
+do i = 1, cells
+ x = xL + (i-0.5)*dx    ! Data is located at center of each cell
+ write (40,*) x, ProblemRegion%u(i)
+end do
+}}}
 
 = [attachment:ExactSolver.tar Exact Riemann Solver] =