anchor.diff:              HFlux(nsrflist(iBCB(i,2))) =
        &      HFlux(nsrflist(iBCB(i,2))) + heat(i)
             enddo
  -c... END DEBUGGING
  -c
  -c          if (countWrong .gt. 0) write(*,*) 'countWrong = ', countWrong
  -c... OLD lines
  -c          Force(1,nsrflist(iBCB(:,2))) = Force(1,nsrflist(iBCB(:,2))) + sum(tau1n)
  -c          Force(2,nsrflist(iBCB(:,2))) = Force(2,nsrflist(iBCB(:,2))) + sum(tau2n)
  -c          Force(3,nsrflist(iBCB(:,2))) = Force(3,nsrflist(iBCB(:,2))) + sum(tau3n)
  -c          PresFor(1,nsrflist(iBCB(:,2))) =
  -c     &    PresFor(1,nsrflist(iBCB(:,2))) + sum(pres1n)
  -c          PresFor(2,nsrflist(iBCB(:,2))) =
  -c     &    PresFor(2,nsrflist(iBCB(:,2))) + sum(pres2n)
  -c          PresFor(3,nsrflist(iBCB(:,2))) =
  -c     &    PresFor(3,nsrflist(iBCB(:,2))) + sum(pres3n)
  -c          HFlux(nsrflist(iBCB(:,2)))   = HFlux(nsrflist(iBCB(:,2)))   + sum(heat)
  -c... End OLD lines
           endif
   c
   c.... end of integration loop
  
positive.code: c
  c.... compute the forces on the body
  c
  c.... DEBUGGING
  c          where (nsrflist(iBCB(:,2)).ge.1)
  c            tau1n = ( pres * bnorm(:,1) - tau1n ) * WdetJb
  c            tau2n = ( pres * bnorm(:,2) - tau2n ) * WdetJb
  c            tau3n = ( pres * bnorm(:,3) - tau3n ) * WdetJb
  c            pres1n =  pres * bnorm(:,1) * WdetJb
  c            pres2n =  pres * bnorm(:,2) * WdetJb
  c            pres3n =  pres * bnorm(:,3) * WdetJb
  c            heat  = - heat * WdetJb
  c          elsewhere
  c            tau1n = zero
  c            tau2n = zero
  c            tau3n = zero
  c            pres1n = zero
  c           pres2n = zero
  c            pres3n = zero
  c            heat  = zero
  c          endwhere
  c
  c.... HARDCODED
  c        loc(1:300)   = 0.0
  c        loc(301:377) = 1.6924 * 1e-6
  c        loc(378:390) = 1.8338 * 1e-6
  c        loc(391:402) = 1.8738 * 1e-6
  c        loc(403:415) = 2.0562 * 1e-6
  c        loc(416:427) = 2.3090 * 1e-6
  c        loc(428:440) = 2.2736 * 1e-6
  c        loc(441:452) = 1.7806 * 1e-6
  
  c        dyn_org = 0.0
  c        do i = 300,lstep      
  c          dyn_org = dyn_org + loc(i-1) 
  c        enddo
  c        countWrong = 0; 
  c.... HARDCODED
  c
          do i = 1, npro
            if (nsrflist(iBCB(i,2)).ge.1) then
  c            if ((xlb(i,1,1) .ge. dyn_org) .and. (iBCB(i,2) .eq. 10)) then
  c               iBCB(i,2)  = 20
  c               countWrong = countWrong + 1 
  c            endif
  c            if ((xlb(i,1,1). lt. dyn_org) .and. (iBCB(i,2) .eq. 20)) then
  c               iBCB(i,2)  = 10
  c               countWrong = countWrong + 1
  c            endif
              tau1n(i) = ( pres(i) * bnorm(i,1) - tau1n(i) ) * WdetJb(i)
              tau2n(i) = ( pres(i) * bnorm(i,2) - tau2n(i) ) * WdetJb(i)
              tau3n(i) = ( pres(i) * bnorm(i,3) - tau3n(i) ) * WdetJb(i)
  
negative.code:      &            um1(npro),                   um2(npro),  
       &            um3(npro)
  c
          dimension lnode(27),               sgn(npro,nshl),
       &            shape(npro,nshl),        shdrv(npro,nsd,nshl)
  c
  