Consider a (1-dimensional) advection-diffusion process, with (constant) diffusivity D > 0, 
for some concentration u(x,t) whose initial profile is u(x,0) = uinit(x), driven by a given velocity V, 
in an interval a ≤ x ≤ b whose ends are impermeable, during some time  0 ≤ t ≤ tend.

1. State precisely and fully the mathematical problem (PDE, where it holds, IC, BCs) modeling this process.
   (the PDE itself, not conservation-law form)

2. Briefly describe a physical problem that could be modeled by this mathematical problem.

Implement the explicit Finite Volume scheme for such advection-diffusion process 
by modifying appropriately your (debugged! Lab7) diffusion code. No exact solution here, so comment out call to COMPARE.
Take as initial profile the "square bump" :   uinit(x) = 5.0  for 1 ≤ x ≤ 2,   uinit(x) = 0  otherwise.
Assume constant velocity V > 0 and diffusivity D > 0;  take tend=4.0, dtout=2, 
and dx = 1/MM with MM = 16. The endpoints a, b will be chosen for each run later. 

Starting at time=0, and every dtout, your code should output the time, number of time-steps, 
and the entire profile of U for plotting, including at time 0 and tend.

Also, as a check of how well it is conserving, your code should compute and output the total area 
under the concentration profile at each output time. What does this area represent? 
You may use the Trapezoidal Rule routine below (or write your own, but very carefully, for our mesh)  

First, consider the case of pure advection ( D = 0 , V > 0 ).

3. State precisely and fully the mathematical problem modeling pure advection (PDE, where it holds, IC, BC).
   Describe what you expect to happen qualitatively and sketch  the initial profile uinit(x) 
   and the expected profile at time t=4 if V=1 and if V=5.  How far will the pulse travel in each case ?

Next, consider the case of pure diffusion ( V = 0 , D > 0 ).

4. State precisely and fully the mathematical problem modeling pure diffusion (PDE, where it holds, IC, BCs).
   Describe what you expect to happen qualitatively for low and high diffusivity and sketch some profiles by hand.

Now we want to examine how the presence of diffusion affects the advection profiles and vice versa.

For each (sub)case listed bellow, do the following:

 Plot the profiles at times 0, dtout, tend on one plot 
     ( Set y-axis range to [0:5.2] to see the top clearly, e.g. in gnuplot: set yrange [0 : 5.2] ). 

 On the plot, mark the parameter values that generated it, mark which curve is at what time  
     
 Below the plot, write the computed area, and comment on what you observe (including the "also try"
    and other explorations), any surprises/unexpected behavior, what you think is happening and why, 
    how this case compares with other cases.

Here are the cases to examine.

5. Pure advection:
          (5a)  V=1., D=0., a=0., b= 7., factor=1.0
          (5b)  V=5., D=0., a=0., b=25., factor=1.0
          (5c)  V=5., D=0., a=0., b=25., factor=1.1
          (5d)  V=5., D=0., a=0., b=25., factor=0.97 (also try 0.99)

6. Pure diffusion:
          (6a)  V=0., D=0.1, a=0., b= 4., factor=0.99 (also try 1.0)
          (6b)  V=0., D=0.5, a=0., b= 4., factor=0.99 (also try 1.0)
          (6c) same as (6a) but with  dx = 1/32  (MM=32).

7. Advection-Diffusion: Low velocity:
          (7a)  V=1., D=0.1, a= 0., b=10., factor=1.0
          (7b)  V=1., D=0.1, a= 0., b=10., factor=0.98
          (7c)  V=1., D=0.5, a=−2., b=13., factor=1.0

8. Advection-Diffusion: High velocity:
          (8a)  V=5., D=0.1, a=0., b=30., factor=1.0
          (8b)  V=5., D=0.5, a=0., b=30., factor=1.0

	double precision FUNCTION TrapzRule( M, dx, U )
		dimension U(0:M+1)
		sum = ( U(0)+3*(U(1)+U(M))+U(M+1) )/4
	 DO  i= 2,M-1
		sum = sum + U(i)
	 ENDDO
	 	TrapzRule = dx * sum
	end