c
c
c =========================================================
      subroutine src1(maxmx,meqn,mbc,mx,xlower,dx,q,maux,aux,t,dt)
c =========================================================
      implicit double precision (a-h,o-z)
      dimension q(1-mbc:maxmx+mbc, meqn)
      dimension aux(1-mbc:maxmx+mbc, maux)
c
c
      parameter (maxmode = 2002)
      dimension qstar(-1:maxmode, 3)
      common /param/  gamma,gamma1
c
c     # source terms for radial symmetry
c
c     # aux(i,1) = x(i) = radial distance

c     # 2-step Runge-Kutta method

      ndim = 3
c
      dt2 = dt/2.d0
      press = 0.d0
      do 10 i=1-mbc,mx+mbc
	 rho = q(i,1)
	 u = q(i,2)/q(i,1)
 	 press = gamma1*(q(i,3) - 0.5d0*rho*u**2)
	 qstar(i,1) = q(i,1) - dt2*(ndim-1)/aux(i,1) * q(i,2)
	 qstar(i,2) = q(i,2) - dt2*(ndim-1)/aux(i,1) * 
     &				(rho*u**2)
	 qstar(i,3) = q(i,3) - dt2*(ndim-1)/aux(i,1) * 
     &                          u*(q(i,3) + press)
   10    continue
c
      do 20 i=1-mbc,mx+mbc
	 rho = qstar(i,1)
	 u = qstar(i,2)/qstar(i,1)
 	 press = gamma1*(qstar(i,3) - 0.5d0*rho*u**2)
	 q(i,1) = q(i,1) - dt*(ndim-1)/aux(i,1) * qstar(i,2)
	 q(i,2) = q(i,2) - dt*(ndim-1)/aux(i,1) * 
     &				(rho*u**2)
	 q(i,3) = q(i,3) - dt*(ndim-1)/aux(i,1) * 
     &                          u*(qstar(i,3) + press)
   20    continue
      return
      end