c    # Several functions needed elsewhere.
c     ===================================================================

      double precision function qtrue(x,t,u)
c     # true solution
      implicit double precision (a-h,o-z)

      qtrue = q0(x-u*t,u) * dexp((sigma(x-u*t)-sigma(x))/u)
c
      return
      end

c     ===================================================================

      double precision function q0(x,u)
c     # initial data
      implicit double precision (a-h,o-z)

      q0 = qss(x,u) 
c     q0 = q0 + 0.3d0*dexp(-8.d0*(x-1.d0)**2) * dcos(10.d0*(x-1.d0))
      if (x.gt.-11.d0 .and. x.lt.-10.5d0) q0 = q0 + 0.2d0
      return
      end

c     ===================================================================

      double precision function qss(x,u)
c     # steady state solution
      implicit double precision (a-h,o-z)

      qss = dexp(-sigma(x)/u)
      return
      end

c     ===================================================================

      double precision function sigma(x)
      implicit double precision (a-h,o-z)
c     # source term function

      A = 8.d0
      beta = 1.d0
      x0 = 5.d0
      B = 0.5d0*dlog(1.0d0/1.4d0)
      gamma = 1.d0
      sigma = A*dexp(-beta*(x-x0)**2) + B*(dtanh(gamma*(x-x0)) + 1.d0)
c
      return
      end