#! /usr/bin/env
#
def lesson01 ( ):

#*****************************************************************************80
#
## lesson01 presents 1D linear convection.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    13 November 2015
#
#  Author:
#
#    Initial version by Lorena Barba
#    This version by John Burkardt
#
#-------------------------------------------------------------------------------
#  We'll start by importing a few libraries to help us out.
#
#  * numpy provides useful matrix operations akin to MATLAB.
#  * matplotlib is a 2D plotting library that we will use to plot our results.
#  * time and sys provide basic timing functions that we'll use to 
#    slow down animations for viewing
#-------------------------------------------------------------------------------
  import numpy
  from matplotlib import pyplot
  import time, sys
#-------------------------------------------------------------------------------
#  Now let's define a few variables; we want to define an evenly spaced grid of 
#  points within a spatial domain that is 2 units of length wide, x in (0,2). 
#
#  nx will be the number of grid points we want.
#  What happens if you rerun the problem with nx = 81?
#
#  dx will be the distance between any pair of adjacent grid points.
#
#  nt will be the number of time steps
#
#  dt will be the size of a single time step.
#
#  v will be the constant velocity of the flow.
#-------------------------------------------------------------------------------
  nx = 41
  dx = 2.0 / ( nx - 1 )
  nt = 25
  dt = 0.025
  v = 1.0
#-------------------------------------------------------------------------------
#  We set up our initial conditions. The initial velocity u0 is 
#  given as u=2 in the interval 0.5 <= x <= 1 and u=1 everywhere else in (0,2) 
#  (i.e., a hat function).  Here, we use the function ones() defining a numpy 
#  array which is nx elements long with every value equal to 1.
#-------------------------------------------------------------------------------
  u = numpy.ones(nx)
  u [ .5/dx : 1/dx+1] = 2.0
  print(u)
#-------------------------------------------------------------------------------
#  Now let's take a look at those initial conditions using a Matplotlib plot. 
#  We've imported the matplotlib plotting library pyplot and the plotting 
#  function is called plot, so we'll call pyplot.plot. To learn about the 
#  myriad possibilities of Matplotlib, explore the Gallery of example plots.
#  Here, we use the syntax for a simple 2D plot: plot(x,y), where the x 
#  values are evenly distributed grid points:
#-------------------------------------------------------------------------------
  x = numpy.linspace ( 0.0, 2.0, nx )
  pyplot.plot ( x, u )
#-------------------------------------------------------------------------------
#  Why doesn't the hat function have perfectly straight sides? Think for a bit.
#
#  Now it's time to implement the discretization of the convection equation 
#  using a finite-difference scheme.  For every element of our array u, we need 
#  to perform the operation 
#
#    u^(n+1)(i) = u^n(i) - v dt dx ( u^n(i) - u^n(i-1) )
#
#  We'll store the result 
#  in a new (temporary) array un, which will be the solution u for the next 
#  time-step. We will repeat this operation for as many time-steps as we specify 
#  and then we can see how far the wave has convected.  We first initialize 
#  our placeholder array un to hold the values we calculate for the n+1 
#  timestep, using once again the NumPy function ones().
#  Then, we may think we have two iterative operations: one in space and one 
#  in time (we'll learn differently later), so we'll start by nesting one 
#  loop inside the other. Note the use of the nifty range() function. 
#
#  When we write: for i in range(1,nx) we will iterate through the u array, 
#  but we'll be skipping the first element (the zero-th element). Why?
#-------------------------------------------------------------------------------
#  initialize a temporary array
#  loop for values of n from 0 to nt, so it will run nt times
#  copy the existing values of u into un
#-------------------------------------------------------------------------------
  un = numpy.ones ( nx ) 
  for n in range ( nt ):
    un = u.copy()
    for i in range ( 1, nx ):
      u[i] = un[i] - v * dt / dx * ( un[i] - un[i-1] )
#-------------------------------------------------------------------------------
#  We will learn later that the code as written above is quite inefficient, 
#  and there are better ways to write this, Python-style. But let's carry on.
#  Now let's try plotting our u array after advancing in time.
#-------------------------------------------------------------------------------
  
  pyplot.plot ( x, u );

#-------------------------------------------------------------------------------
#  OK! So our hat function has definitely moved to the right, but it's no 
#  longer a hat. What's going on?
#-------------------------------------------------------------------------------
#  LEARN MORE:
#  For a more thorough explanation of the finite difference method, including
#  topics like truncation error and order of convergence, watch Video Lessons
#  2 and 3 by Professor Barba on YouTube.
#-------------------------------------------------------------------------------
  return