anchor.diff:  c
   c        # second order corrections:
   
  -         do 120 m=1,meqn
  +         do m = 1,meqn
               cqxx(m,i) = 0.d0
  -            do 119 mw=1,mwaves
  +            do mw=1,mwaves
   c
                  if (use_fwaves) then
                      abs_sign = dsign(1.d0,s(mw,i))
                    else
                      abs_sign = dabs(s(mw,i))
  -                 endif
  +               endif
   
                  cqxx(m,i) = cqxx(m,i) + abs_sign
        &             * (1.d0 - dabs(s(mw,i))*dtdxave) * wave(m,mw,i)
  -c
  -  119          continue
  -            faddm(m,i) = faddm(m,i) + 0.5d0 * cqxx(m,i)
  -            faddp(m,i) = faddp(m,i) + 0.5d0 * cqxx(m,i)
  -  120       continue
  +             end do
  +         end do
  +         faddm(m,i) = faddm(m,i) + 0.5d0 * cqxx(m,i)
  +         faddp(m,i) = faddp(m,i) + 0.5d0 * cqxx(m,i)
  +      end do
   c
   c
     130  continue
  
positive.code:                gupdate = 0.5d0*dtdx1d(i-1) * bpasdq(m,i)
                 gaddm(m,i-1,2) = gaddm(m,i-1,2) - gupdate
                 gaddp(m,i-1,2) = gaddp(m,i-1,2) - gupdate
    160          continue
  c
  c     # split the right-going flux difference into down-going and up-going:
        call rpt2(ixy,2,maxm,meqn,mwaves,maux,mbc,mx,
  
negative.code:       if (limit) call clawpack5_inlinelimiter(maxm,meqn,mwaves,mbc,mx,
       &      wave,s,mthlim)
  c
        do  i = 1, mx+1
  c
  c        # For correction terms below, need average of dtdx in cell
  c        # i-1 and i.  Compute these and overwrite dtdx1d:
  