diff --git a/fvmbook/chap22/corner/Makefile b/fvmbook/chap22/corner/Makefile
new file mode 100644
index 00000000..5249fb58
--- /dev/null
+++ b/fvmbook/chap22/corner/Makefile
@@ -0,0 +1,65 @@
+
+# Makefile for Clawpack code in this directory.
+# This version only sets the local files and frequently changed
+# options, and then includes the standard makefile pointed to by CLAWMAKE.
+CLAWMAKE = $(CLAW)/clawutil/src/Makefile.common
+
+# See the above file for details and a list of make options, or type
+#   make .help
+# at the unix prompt.
+
+
+# Adjust these variables if desired:
+# ----------------------------------
+
+CLAW_PKG = classic                  # Clawpack package to use
+EXE = xclaw                         # Executable to create
+SETRUN_FILE = setrun.py             # File containing function to make data
+OUTDIR = _output                    # Directory for output
+SETPLOT_FILE = setplot.py           # File containing function to set plots
+PLOTDIR = _plots                    # Directory for plots
+
+OVERWRITE ?= True                   # False ==> make a copy of OUTDIR first
+RESTART ?= False                    # Should = clawdata.restart in setrun
+
+# Environment variable FC should be set to fortran compiler, e.g. gfortran
+
+# Compiler flags can be specified here or set as an environment variable
+FFLAGS ?=  
+
+# ---------------------------------
+# List of sources for this program:
+# ---------------------------------
+
+LIB = $(CLAW)/classic/src/2d
+
+MODULES = \
+
+SOURCES = \
+  qinit.f \
+  setprob.f \
+  fdisc.f \
+  setaux.f \
+  $(CLAW)/classic/src/2d/cellave.f \
+  $(CLAW)/riemann/src/rpn2_vc_elasticity.f90 \
+  $(CLAW)/riemann/src/rpt2_vc_elasticity.f90 \
+  $(CLAW)/classic/src/2d/driver.f90 \
+  $(CLAW)/classic/src/2d/claw2ez.f \
+  $(CLAW)/classic/src/2d/claw2.f \
+  $(CLAW)/classic/src/2d/bc2.f \
+  $(CLAW)/classic/src/2d/b4step2.f90 \
+  $(CLAW)/classic/src/2d/step2.f90 \
+  $(CLAW)/classic/src/2d/step2ds.f90 \
+  $(CLAW)/classic/src/2d/dimsp2.f \
+  $(CLAW)/classic/src/2d/flux2.f90 \
+  $(CLAW)/classic/src/2d/copyq2.f \
+  $(CLAW)/classic/src/2d/inlinelimiter.f90 \
+  $(CLAW)/classic/src/2d/src2.f90 \
+  $(CLAW)/classic/src/2d/out2.f \
+  $(CLAW)/classic/src/2d/restart2.f \
+  $(CLAW)/classic/src/2d/opendatafile.f
+
+#-------------------------------------------------------------------
+# Include Makefile containing standard definitions and make options:
+include $(CLAWMAKE)
+
diff --git a/fvmbook/chap22/corner/README.rst b/fvmbook/chap22/corner/README.rst
new file mode 100644
index 00000000..4102f46a
--- /dev/null
+++ b/fvmbook/chap22/corner/README.rst
@@ -0,0 +1,24 @@
+
+.. _fvmbook_chap22/corner:
+
+Elasticity equations in a region consisting of two materials with a corner
+--------------------------------------------------------------------------
+
+    
+Example [book/chap22/corner] to accompany the book 
+`Finite Volume Methods for Hyperbolic Problems
+<http://www.clawpack.org/book.html>`_
+by R. J. LeVeque.
+
+Converted to Clawpack 5.0 form in 2014.
+
+
+**Note:**
+
+ - `mx=400, my=400` was used for the plots in the book.
+ - The elastic parameters are defined in `setprob.data`
+ - and the inclusion geometry is specified by `fdisc.f`
+ - Elastic parameters in each cell are set in `setaux.f`
+ - The wave is generated by specifying a right-going P-wave in `qinit.f`
+
+    
diff --git a/fvmbook/chap22/corner/fdisc.f b/fvmbook/chap22/corner/fdisc.f
new file mode 100644
index 00000000..2b2da5e3
--- /dev/null
+++ b/fvmbook/chap22/corner/fdisc.f
@@ -0,0 +1,24 @@
+c
+c
+c
+c     =================================================
+      function fdisc(x,y)
+c     =================================================
+
+      implicit double precision (a-h,o-z)
+
+c
+c     # for computing cell averages for initial data that has a
+c     # discontinuity along some curve.  fdisc should be negative to the 
+c     # left of the curve and positive to the right
+
+c     # half wedge
+      if (x .gt. 0.d0 .and. (y.lt.(0.55d0*x))) then
+          fdisc = 1.d0
+        else
+          fdisc = -1.d0
+        endif
+c
+      return
+      end
+
diff --git a/fvmbook/chap22/corner/qinit.f b/fvmbook/chap22/corner/qinit.f
new file mode 100644
index 00000000..f73c6d3e
--- /dev/null
+++ b/fvmbook/chap22/corner/qinit.f
@@ -0,0 +1,38 @@
+
+c
+c
+c
+c     =====================================================
+       subroutine qinit(meqn,mbc,mx,my,xlower,ylower,
+     &                   dx,dy,q,maux,aux)
+c     =====================================================
+c
+c     # Set initial conditions for q.
+c
+       implicit double precision (a-h,o-z)
+       dimension q(meqn,1-mbc:mx+mbc,1-mbc:my+mbc)
+       dimension aux(maux,1-mbc:mx+mbc,1-mbc:my+mbc)
+       common /comaux/ rho1,amu1,alam1,rho2,amu2,alam2
+c
+       cp2 = dsqrt((alam2+2.d0*amu2)/rho2)
+       do 20 i=1,mx
+          xi = xlower + (i-0.5d0)*dx
+          do 20 j=1,my
+             yj = ylower + (j-0.5d0)*dy
+             q(1,i,j) = 0.d0
+             q(2,i,j) = 0.d0
+             q(3,i,j) = 0.d0
+             q(4,i,j) = 0.d0
+             q(5,i,j) = 0.d0
+     	     if (xi .gt. -0.35d0 .and. xi.lt.-0.2d0) then
+c                # set to be an eigenvector corresponding to right-going
+c                # P-wave:
+    	         q(1,i,j) = aux(2,i,j) + 2.d0*aux(3,i,j)  
+    	         q(2,i,j) = aux(2,i,j)
+    	         q(4,i,j) = -aux(4,i,j)
+        		 endif
+  20         continue
+
+       return
+       end
+
diff --git a/fvmbook/chap22/corner/setaux.f b/fvmbook/chap22/corner/setaux.f
new file mode 100644
index 00000000..fe0fc88a
--- /dev/null
+++ b/fvmbook/chap22/corner/setaux.f
@@ -0,0 +1,42 @@
+c     =========================================================
+      subroutine setaux(mbc,mx,my,xlower,ylower,dx,dy,maux,aux)
+c     =========================================================
+c
+c     # set auxiliary arrays 
+c     # variable coefficient acoustics
+c     #  aux(1,i,j) = density rho in (i,j) cell
+c     #  aux(2,i,j) = lambda in (i,j) cell
+c     #  aux(3,i,j) = mu in (i,j) cell
+c     #  aux(4,i,j) = cp in (i,j) cell
+c     #  aux(5,i,j) = cs in (i,j) cell
+c
+c     # Piecewise constant medium
+c     # Material parameters are set in setprob.f
+
+c
+c     
+      implicit double precision (a-h,o-z)
+      dimension aux(5,1-mbc:mx+mbc,1-mbc:my+mbc)
+      common /comaux/ rho1,amu1,alam1,rho2,amu2,alam2
+c
+
+      do 30 j=1-mbc,my+mbc
+       do 20 i=1-mbc,mx+mbc
+          xl = xlower + (i-1.0d0)*dx
+          yl = ylower + (j-1.0d0)*dy
+          call cellave(xl,yl,dx,dy,w1)
+          w2 = 1.d0 - w1
+
+          aux(1,i,j) = w1*rho1 + w2*rho2
+          aux(2,i,j) = w1*alam1 + w2*alam2
+          aux(3,i,j) = w1*amu1 + w2*amu2
+          bulk       = aux(2,i,j) + 2.d0*aux(3,i,j)
+          aux(4,i,j) = dsqrt(bulk/aux(1,i,j))
+          aux(5,i,j) = dsqrt(aux(3,i,j)/aux(1,i,j))
+
+   20     continue
+   30    continue
+
+       return
+       end
+
diff --git a/fvmbook/chap22/corner/setplot.py b/fvmbook/chap22/corner/setplot.py
new file mode 100644
index 00000000..05870e1d
--- /dev/null
+++ b/fvmbook/chap22/corner/setplot.py
@@ -0,0 +1,125 @@
+
+""" 
+Set up the plot figures, axes, and items to be done for each frame.
+
+This module is imported by the plotting routines and then the
+function setplot is called to set the plot parameters.
+    
+""" 
+
+#--------------------------
+def setplot(plotdata):
+#--------------------------
+    
+    """ 
+    Specify what is to be plotted at each frame.
+    Input:  plotdata, an instance of clawpack.visclaw.data.ClawPlotData.
+    Output: a modified version of plotdata.
+    
+    """ 
+
+
+    from clawpack.visclaw import colormaps
+    from numpy import linspace
+
+    plotdata.clearfigures()  # clear any old figures,axes,items data
+    
+    def plot_corner(current_data):
+        from pylab import plot
+        plot([.0,.0],[-1,0],'r',linewidth=2)
+        plot([.0,1],[0,.55],'r',linewidth=2)
+
+    def sigmatr(current_data):
+        # trace of sigma
+        q = current_data.q
+        return q[0,:,:] + q[1,:,:]
+
+
+
+    # Figure for trace of sigma
+    plotfigure = plotdata.new_plotfigure(name='trace(sigma)', figno=0)
+
+    # Set up for axes in this figure:
+    plotaxes = plotfigure.new_plotaxes()
+    plotaxes.xlimits = 'auto'
+    plotaxes.ylimits = 'auto'
+    plotaxes.title = 'trace(sigma)'
+    plotaxes.scaled = True
+    plotaxes.afteraxes = plot_corner
+
+    # Set up for item on these axes:
+    plotitem = plotaxes.new_plotitem(plot_type='2d_pcolor')
+    plotitem.plot_var = sigmatr
+    plotitem.pcolor_cmap = colormaps.red_yellow_blue
+    plotitem.pcolor_cmin = -1.
+    plotitem.pcolor_cmax = 1.
+    plotitem.add_colorbar = True
+    plotitem.show = True       # show on plot?
+    
+
+
+    # Figure for shear stress
+    plotfigure = plotdata.new_plotfigure(name='shear', figno=1)
+
+    # Set up for axes in this figure:
+    plotaxes = plotfigure.new_plotaxes()
+    plotaxes.xlimits = 'auto'
+    plotaxes.ylimits = 'auto'
+    plotaxes.title = 'shear stress'
+    plotaxes.scaled = True
+    plotaxes.afteraxes = plot_corner
+
+    # Set up for item on these axes:
+    plotitem = plotaxes.new_plotitem(plot_type='2d_pcolor')
+    plotitem.plot_var = 2  # sigma_12
+    plotitem.pcolor_cmap = colormaps.red_yellow_blue
+    plotitem.pcolor_cmin = -0.2
+    plotitem.pcolor_cmax = 0.2
+    plotitem.add_colorbar = True
+    plotitem.show = True       # show on plot?
+    
+
+    # Figure for contours 
+    plotfigure = plotdata.new_plotfigure(name='contours', figno=2)
+
+    # Set up for axes in this figure:
+    plotaxes = plotfigure.new_plotaxes()
+    plotaxes.xlimits = 'auto'
+    plotaxes.ylimits = 'auto'
+    plotaxes.title = 'pressure(black) and shear(green)'
+    plotaxes.scaled = True
+    plotaxes.afteraxes = plot_corner
+
+    # Set up for item on these axes:
+    plotitem = plotaxes.new_plotitem(plot_type='2d_contour')
+    plotitem.plot_var = sigmatr
+    plotitem.contour_levels = linspace(-2,8,50)
+    plotitem.contour_colors = 'k'
+    plotitem.show = True       # show on plot?
+
+    # Set up for item on these axes:
+    plotitem = plotaxes.new_plotitem(plot_type='2d_contour')
+    plotitem.plot_var = 2  # sigma_12
+    plotitem.contour_levels = linspace(-0.4,0.4,30)
+    plotitem.contour_colors = 'g'
+    plotitem.show = True       # show on plot?
+    
+    
+
+    # Parameters used only when creating html and/or latex hardcopy
+    # e.g., via clawpack.visclaw.frametools.printframes:
+
+    plotdata.printfigs = True                # print figures
+    plotdata.print_format = 'png'            # file format
+    plotdata.print_framenos = 'all'          # list of frames to print
+    plotdata.print_fignos = 'all'            # list of figures to print
+    plotdata.html = True                     # create html files of plots?
+    plotdata.html_homelink = '../README.html'   # pointer for top of index
+    plotdata.latex = True                    # create latex file of plots?
+    plotdata.latex_figsperline = 2           # layout of plots
+    plotdata.latex_framesperline = 1         # layout of plots
+    plotdata.latex_makepdf = False           # also run pdflatex?
+
+    return plotdata
+
+    
diff --git a/fvmbook/chap22/corner/setprob.f b/fvmbook/chap22/corner/setprob.f
new file mode 100644
index 00000000..9f3d45fd
--- /dev/null
+++ b/fvmbook/chap22/corner/setprob.f
@@ -0,0 +1,31 @@
+      subroutine setprob
+      implicit double precision (a-h,o-z)
+      character*12 fname
+      common /comaux/ rho1,amu1,alam1,rho2,amu2,alam2
+
+c
+c     # Set the material parameters for the elasticity equations
+c
+c
+      iunit = 7
+      fname = 'setprob.data'
+c     # open the unit with new routine from Clawpack 4.4 to skip over
+c     # comment lines starting with #:
+      call opendatafile(iunit, fname)
+                
+
+c
+c     # Piecewise constant medium 
+c     # Material parameters
+
+      read(7,*) rho1
+      read(7,*) alam1
+      read(7,*) amu1
+
+      read(7,*) rho2
+      read(7,*) alam2
+      read(7,*) amu2
+
+      return
+      end
+
diff --git a/fvmbook/chap22/corner/setrun.py b/fvmbook/chap22/corner/setrun.py
new file mode 100644
index 00000000..a4998ef4
--- /dev/null
+++ b/fvmbook/chap22/corner/setrun.py
@@ -0,0 +1,243 @@
+""" 
+Module to set up run time parameters for Clawpack -- classic code.
+
+The values set in the function setrun are then written out to data files
+that will be read in by the Fortran code.
+    
+""" 
+
+import os
+import numpy as np
+
+#------------------------------
+def setrun(claw_pkg='classic'):
+#------------------------------
+    
+    """ 
+    Define the parameters used for running Clawpack.
+
+    INPUT:
+        claw_pkg expected to be "classic" for this setrun.
+
+    OUTPUT:
+        rundata - object of class ClawRunData 
+    
+    """ 
+    
+    from clawpack.clawutil import data 
+    
+    
+    assert claw_pkg.lower() == 'classic',  "Expected claw_pkg = 'classic'"
+
+    num_dim = 2
+    rundata = data.ClawRunData(claw_pkg, num_dim)
+
+    #------------------------------------------------------------------
+    # Problem-specific parameters to be written to setprob.data:
+    #------------------------------------------------------------------
+    # Sample setup to write one line to setprob.data ...
+    probdata = rundata.new_UserData(name='probdata',fname='setprob.data')
+    probdata.add_param('rho1',     1.,  'rho1')
+    probdata.add_param('lam1',     4.,  'lam1')
+    probdata.add_param('mu1',     0.5,  'mu1')
+    probdata.add_param('rho2',     1.,  'rho2')
+    probdata.add_param('lam2',     2.,  'lam2')
+    probdata.add_param('mu2',     1.0,  'mu2')
+    
+    #------------------------------------------------------------------
+    # Standard Clawpack parameters to be written to claw.data:
+    #------------------------------------------------------------------
+
+    clawdata = rundata.clawdata  # initialized when rundata instantiated
+
+
+    # ---------------
+    # Spatial domain:
+    # ---------------
+
+    # Number of space dimensions:
+    clawdata.num_dim = num_dim
+    
+    # Lower and upper edge of computational domain:
+    clawdata.lower[0] = -1.000000e+00          # xlower
+    clawdata.upper[0] = 1.000000e+00          # xupper
+    clawdata.lower[1] = -1.000000e+00          # ylower
+    clawdata.upper[1] = 1.000000e+00          # yupper
+    
+    # Number of grid cells:
+    clawdata.num_cells[0] = 200      # mx
+    clawdata.num_cells[1] = 200      # my
+    
+
+    # ---------------
+    # Size of system:
+    # ---------------
+
+    # Number of equations in the system:
+    clawdata.num_eqn = 5
+
+    # Number of auxiliary variables in the aux array (initialized in setaux)
+    clawdata.num_aux = 5
+    
+    # Index of aux array corresponding to capacity function, if there is one:
+    clawdata.capa_index = 0
+    
+    
+    # -------------
+    # Initial time:
+    # -------------
+
+    clawdata.t0 = 0.000000
+    
+
+    # Restart from checkpoint file of a previous run?
+    # Note: If restarting, you must also change the Makefile to set:
+    #    RESTART = True
+    # If restarting, t0 above should be from original run, and the
+    # restart_file 'fort.qNNNN' specified below should be in 
+    # the OUTDIR indicated in Makefile.
+
+    clawdata.restart = False               # True to restart from prior results
+    clawdata.restart_file = 'fort.q0006'   # File to use for restart data
+    
+    
+    # -------------
+    # Output times:
+    #--------------
+
+    # Specify at what times the results should be written to fort.q files.
+    # Note that the time integration stops after the final output time.
+ 
+    clawdata.output_style = 1
+ 
+    if clawdata.output_style==1:
+        # Output ntimes frames at equally spaced times up to tfinal:
+        # Can specify num_output_times = 0 for no output
+        clawdata.num_output_times = 10
+        clawdata.tfinal = 0.500000
+        clawdata.output_t0 = True  # output at initial (or restart) time?
+        
+    elif clawdata.output_style == 2:
+        # Specify a list or numpy array of output times:
+        # Include t0 if you want output at the initial time.
+        clawdata.output_times =  [0., 0.1]
+ 
+    elif clawdata.output_style == 3:
+        # Output every step_interval timesteps over total_steps timesteps:
+        clawdata.output_step_interval = 2
+        clawdata.total_steps = 4
+        clawdata.output_t0 = True  # output at initial (or restart) time?
+        
+
+    clawdata.output_format = 'ascii'      # 'ascii', 'binary', 'netcdf'
+
+    clawdata.output_q_components = 'all'   # could be list such as [True,True]
+    clawdata.output_aux_components = 'none'  # could be list
+    clawdata.output_aux_onlyonce = True    # output aux arrays only at t0
+    
+
+    # ---------------------------------------------------
+    # Verbosity of messages to screen during integration:  
+    # ---------------------------------------------------
+
+    # The current t, dt, and cfl will be printed every time step
+    # at AMR levels <= verbosity.  Set verbosity = 0 for no printing.
+    #   (E.g. verbosity == 2 means print only on levels 1 and 2.)
+    clawdata.verbosity = 1
+    
+    
+
+    # --------------
+    # Time stepping:
+    # --------------
+
+    # if dt_variable==True:  variable time steps used based on cfl_desired,
+    # if dt_variable==False: fixed time steps dt = dt_initial always used.
+    clawdata.dt_variable = True
+    
+    # Initial time step for variable dt.  
+    # (If dt_variable==0 then dt=dt_initial for all steps)
+    clawdata.dt_initial = 1.000000e-02
+    
+    # Max time step to be allowed if variable dt used:
+    clawdata.dt_max = 1.000000e+99
+    
+    # Desired Courant number if variable dt used 
+    clawdata.cfl_desired = 0.900000
+    # max Courant number to allow without retaking step with a smaller dt:
+    clawdata.cfl_max = 1.000000
+    
+    # Maximum number of time steps to allow between output times:
+    clawdata.steps_max = 500
+
+
+    # ------------------
+    # Method to be used:
+    # ------------------
+
+    # Order of accuracy:  1 => Godunov,  2 => Lax-Wendroff plus limiters
+    clawdata.order = 2
+    
+    # Use dimensional splitting? (not yet available for AMR)
+    clawdata.dimensional_split = 'unsplit'
+    
+    # For unsplit method, transverse_waves can be 
+    #  0 or 'none'      ==> donor cell (only normal solver used)
+    #  1 or 'increment' ==> corner transport of waves
+    #  2 or 'all'       ==> corner transport of 2nd order corrections too
+    clawdata.transverse_waves = 2
+    
+    
+    # Number of waves in the Riemann solution:
+    clawdata.num_waves = 4
+    
+    # List of limiters to use for each wave family:  
+    # Required:  len(limiter) == num_waves
+    # Some options:
+    #   0 or 'none'     ==> no limiter (Lax-Wendroff)
+    #   1 or 'minmod'   ==> minmod
+    #   2 or 'superbee' ==> superbee
+    #   3 or 'vanleer'  ==> van Leer
+    #   4 or 'mc'       ==> MC limiter
+    clawdata.limiter = ['mc', 'mc', 'mc', 'mc']
+    
+    clawdata.use_fwaves = False    # True ==> use f-wave version of algorithms
+    
+    # Source terms splitting:
+    #   src_split == 0 or 'none'    ==> no source term (src routine never called)
+    #   src_split == 1 or 'godunov' ==> Godunov (1st order) splitting used, 
+    #   src_split == 2 or 'strang'  ==> Strang (2nd order) splitting used,  not recommended.
+    clawdata.source_split = 0
+    
+    
+    # --------------------
+    # Boundary conditions:
+    # --------------------
+
+    # Number of ghost cells (usually 2)
+    clawdata.num_ghost = 2
+    
+    # Choice of BCs at xlower and xupper:
+    #   0 or 'user'     => user specified (must modify bcNamr.f to use this option)
+    #   1 or 'extrap'   => extrapolation (non-reflecting outflow)
+    #   2 or 'periodic' => periodic (must specify this at both boundaries)
+    #   3 or 'wall'     => solid wall for systems where q(2) is normal velocity
+    
+    clawdata.bc_lower[0] = 'extrap'   # at xlower
+    clawdata.bc_upper[0] = 'extrap'   # at xupper
+
+    clawdata.bc_lower[1] = 'extrap'   # at ylower
+    clawdata.bc_upper[1] = 'extrap'   # at yupper
+                  
+    return rundata
+
+    # end of function setrun
+    # ----------------------
+
+
+if __name__ == '__main__':
+    # Set up run-time parameters and write all data files.
+    import sys
+    rundata = setrun(*sys.argv[1:])
+    rundata.write()
+