anchor.diff:      implicit none
       real(kind=kreal), intent(inout)  :: ZP(:)
       integer(kind=kint) :: in, i, j, isL, ieL, isU, ieU, k, iold, rcm,NDOF,NDOF2,idof,jdof
  -    real(kind=kreal) :: SW1, SW2, SW3, X1, X2, X3
       real(kind=kreal) :: SW(NDOF),X(NDOF)
       NDOF2=NDOF*NDOF
  -      !C-- FORWARD
  -      do i= 1, N
  +    !C-- FORWARD
  +    do i= 1, N
  +      do idof = 1, NDOF
  +        SW(idof) = ZP(NDOF*(i-1)+idof)
  +      end do 
  +      isL= inumFI1L(i-1)+1
  +      ieL= inumFI1L(i)
  +      do j= isL, ieL
  +        in= FI1L(j)
           do idof = 1, NDOF
  -          SW(idof) = ZP(NDOF*(i-1)+idof)
  -        end do 
  -        isL= inumFI1L(i-1)+1
  -        ieL= inumFI1L(i)
  -        do j= isL, ieL
  -          in= FI1L(j)
  -          do idof = 1, NDOF
  -            X(idof) = ZP(NDOF*(in-1)+idof)
  +          X(idof) = ZP(NDOF*(in-1)+idof)
  +        end do 
  +        do idof = 1, NDOF
  +          do jdof = 1, NDOF
  +            SW(idof) = SW(idof) - RIFL(NDOF2*(j-1)+NDOF*(idof-1)+jdof)*X(jdof)
             end do 
  -          do idof = 1, NDOF
  -            do jdof = 1, NDOF
  -              SW(idof) = SW(idof) - RIFL(NDOF2*(j-1)+NDOF*(idof-1)+jdof)*X(jdof)
  -            end do 
  +        end do 
  +      enddo
  +      X = SW
  +      do idof = 2,NDOF
  +        do jdof = 1, idof-1
  +           X(idof) = X(idof) - RIFD(NDOF2*(i-1)+NDOF*(idof-1)+jdof )*X(jdof)
  +        end do 
  +      end do
  +      ZP(NDOF*(i-1)+1:NDOF*(i-1)+NDOF) = X(1:NDOF)
  +    enddo
  +
  +    do i=1, N
  +      do idof=1,NDOF
  +        ZP(NDOF*(i-1)+idof)= ZP(NDOF*(i-1)+idof)*RIFD(NDOF2*(i-1)+(idof-1)*NDOF+idof)
  +      end do
  +    enddo
  +  !C-- BACKWARD
  +    do i= N, 1, -1
  +      do idof = 1, NDOF
  +        SW(idof) = ZP(NDOF*(i-1)+idof)
  +      end do 
  +      isL= inumFI1U(i-1) + 1
  +      ieL= inumFI1U(i)
  +      do j= ieL, isL, -1
  +        in= FI1U(j)
  +        do idof = 1, NDOF
  +          X(idof) = ZP(NDOF*(in-1)+idof)
  +        end do 
  +        do idof = 1, NDOF
  +          do jdof = 1, NDOF
  +            SW(idof) = SW(idof) - RIFU(NDOF2*(j-1)+NDOF*(idof-1)+jdof)*X(jdof)
             end do 
  -        enddo
  -        X = SW
  -        do idof = 2,NDOF
  -          do jdof = 1, idof-1
  -             X(idof) = X(idof) - RIFD(NDOF2*(i-1)+NDOF*(idof-1)+jdof )*X(jdof)
  -          end do 
  -        end do
  -        ZP(NDOF*(i-1)+1:NDOF*(i-1)+NDOF) = X(1:NDOF)
  -      enddo
  -
  -      do i=1, N
  -        do idof=1,NDOF
  -          ZP(NDOF*(i-1)+idof)= ZP(NDOF*(i-1)+idof)*RIFD(NDOF2*(i-1)+(idof-1)*NDOF+idof)
  -        end do
  -      enddo
  -    !C-- BACKWARD
  -      do i= N, 1, -1
  -        do idof = 1, NDOF
  -          SW(idof) = ZP(NDOF*(i-1)+idof)
  -        end do 
  -        isL= inumFI1U(i-1) + 1
  -        ieL= inumFI1U(i)
  -        do j= ieL, isL, -1
  -          in= FI1U(j)
  -          do idof = 1, NDOF
  -            X(idof) = ZP(NDOF*(in-1)+idof)
  -          end do 
  -          do idof = 1, NDOF
  -            do jdof = 1, NDOF
  -              SW(idof) = SW(idof) - RIFU(NDOF2*(j-1)+NDOF*(idof-1)+jdof)*X(jdof)
  -            end do 
  -          end do 
  -        enddo
  -        X=SW
  -        do idof = NDOF, 1, -1
  -          do jdof = NDOF, idof+1, -1
  -            X(idof) = X(idof) - RIFD(NDOF*NDOF*(i-1)+NDOF*(jdof-1)+idof)*X(jdof)          
  -          end do 
  -        end do
  -        do idof = 1, NDOF
  -          ZP(NDOF*(i-1)+idof) = X(idof)
  -        end do 
  -      enddo
  +        end do 
  +      enddo
  +      X=SW
  +      do idof = NDOF, 1, -1
  +        do jdof = NDOF, idof+1, -1
  +          X(idof) = X(idof) - RIFD(NDOF*NDOF*(i-1)+NDOF*(jdof-1)+idof)*X(jdof)          
  +        end do 
  +      end do
  +      do idof = 1, NDOF
  +        ZP(NDOF*(i-1)+idof) = X(idof)
  +      end do 
  +    enddo
     end subroutine hecmw_precond_RIF_nn_apply
   
   
  
positive.code:     subroutine hecmw_rif_nn(hecMAT)
      implicit none
      type (hecmwST_matrix)     :: hecMAT
      integer(kind=kint) :: i, j, jS, jE, in, itr, iitr, ind(9), PRECOND, NP, NDOF, NDOF2,idof,jdof,k
      real(kind=krealp) :: YV1, YV2, YV3, X1, X2, X3, dd, dd1, dd2, dd3, dtemp(3)
      real(kind=krealp) :: YV(hecMAT%NDOF), X(hecMAT%NDOF)
      real(kind=krealp) :: FILTER, SIGMA_DIAG
      real(kind=krealp), allocatable :: zz(:), vv(:)
  
  
negative.code:     allocate(zz(NDOF*NP))
  
      do itr=1,NP
        do iitr=1,NDOF
          zz(:) = 0.0d0
          vv(:) = 0.0d0
          !{v}=[A]{zi}
          do idof = 1,NDOF
            zz(NDOF*(itr-1)+idof)= SAINVD(NDOF2*(itr-1)+NDOF*(idof-1)+iitr)
          end do 
          zz(NDOF*(itr-1)+iitr)= 1.0d0
  
          jS= inumFI1L(itr-1) + 1
          jE= inumFI1L(itr  )
          do j= jS, jE
            in  = FI1L(j)
            do idof = 1, NDOF
              zz(NDOF*(in-1)+idof)=SAINVL(NDOF2*(j-1)+NDOF*(iitr-1)+idof)
            end do 
          enddo
  
          do i= 1, itr
            do idof = 1,NDOF
              X(idof)=zz(NDOF*(i-1)+idof)
            end do 
            do idof = 1, NDOF
              do jdof = 1, NDOF
                vv(NDOF*(i-1)+idof) = vv(NDOF*(i-1)+idof) + D(NDOF2*(i-1)+NDOF*(idof-1)+jdof)*X(jdof)
              end do 
            end do 
  
            jS= indexL(i-1) + 1
            jE= indexL(i  )
            do j=jS,jE
              in = itemL(j)
              do idof = 1, NDOF
                do jdof = 1, NDOF
                  vv(NDOF*(in-1)+idof) = vv(NDOF*(in-1)+idof) + AL(NDOF2*(j-1)+NDOF*(jdof-1)+idof)*X(jdof)
                end do 
              end do 
            enddo
            jS= indexU(i-1) + 1
            jE= indexU(i  )
            do j= jS, jE
              in = itemU(j)
              do idof = 1, NDOF
                do jdof = 1, NDOF
                  vv(NDOF*(in-1)+idof) = vv(NDOF*(in-1)+idof) + AU(NDOF2*(j-1)+NDOF*(jdof-1)+idof)*X(jdof)
                end do 
              end do 
            enddo
          enddo
  
          !{d}={v^t}{z_j}
          do idof = 1, NDOF
            Dtmp(idof)= SAINVD(NDOF2*(itr-1)+NDOF*(idof-1)+idof)
          end do 
          
          do i= itr,NP
            do idof = 1,NDOF
              SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof) = vv(NDOF*(i-1)+idof)
              do jdof = 1, idof-1
                SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof) = SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof) + &
                 & SAINVD(NDOF2*(i-1)+NDOF*(jdof-1)+idof)*vv(NDOF*(i-1)+jdof)
              end do 
            end do
            jS= inumFI1L(i-1) + 1
            jE= inumFI1L(i  )
            do j= jS, jE
              in  = FI1L(j)
              do idof = 1,NDOF
                X(idof)=vv(NDOF*(in-1)+idof)
              end do 
              do idof = 1, NDOF
                do jdof = 1, NDOF
                  SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof) = SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof) + &
                         & SAINVL(NDOF2*(j-1)+NDOF*(idof-1)+jdof)*X(jdof)
                end do 
              end do 
            enddo
          enddo
  
          !Update D
          dd = 1.0d0/SAINVD(NDOF2*(itr-1)+NDOF*(iitr-1)+iitr)
          do idof=1,iitr-1
            SAINVD(NDOF2*(itr-1)+NDOF*(idof-1)+idof) = Dtmp(idof)
          end do 
          SAINVD(NDOF2*(itr-1)+NDOF*(iitr-1)+iitr)=dd
          do idof = iitr+1, NDOF
            SAINVD(NDOF2*(itr-1)+NDOF*(idof-1)+idof) = SAINVD(NDOF2*(itr-1)+NDOF*(idof-1)+idof)*dd
          end do 
  
          do i =itr+1,NP
            do idof = 1, NDOF
              SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof) = SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof)*dd
            end do 
          enddo
  
          do idof = iitr, NDOF
            RIFD(NDOF2*(itr-1)+NDOF*(idof-1)+iitr) = SAINVD(NDOF2*(itr-1)+NDOF*(idof-1)+idof) !RIF UPDATE
          end do 
  
          jS= inumFI1U(itr-1) + 1
          jE= inumFI1U(itr  )
          do j= jS, jE
            do idof = 1, NDOF
              RIFU(NDOF2*(j-1)+NDOF*(iitr-1)+idof) = SAINVD(NDOF2*(FI1U(j)-1)+NDOF*(idof-1)+idof)
            end do 
          enddo
  
          !Update Z
          do k=iitr+1,NDOF
            dd = SAINVD(NDOF2*(itr-1)+NDOF*(k-1)+k)
            if(abs(dd) > FILTER)then
              do jdof = 1, iitr 
                SAINVD(NDOF2*(itr-1)+NDOF*(jdof-1)+k)= SAINVD(NDOF2*(itr-1)+NDOF*(jdof-1)+k) - dd*zz(NDOF*(itr-1)+jdof)
              end do
              jS= inumFI1L(itr-1) + 1
              jE= inumFI1L(itr  )
              do j= jS, jE
                in  = FI1L(j)
                do idof = 1, NDOF
                  SAINVL(NDOF2*(j-1)+NDOF*(k-1)+idof) = SAINVL(NDOF2*(j-1)+NDOF*(k-1)+idof)-dd*zz(NDOF*(in-1)+idof)
                end do
              enddo
            endif
          end do 
          
          do i= itr +1,NP
            jS= inumFI1L(i-1) + 1
            jE= inumFI1L(i  )
            do idof = 1, NDOF
              dd = SAINVD(NDOF2*(i-1)+NDOF*(idof-1)+idof)
              if(abs(dd) > FILTER)then
                do j= jS, jE
                  in  = FI1L(j)
                  if (in > itr) exit
                  do jdof=1,NDOF
                  SAINVL(NDOF2*(j-1)+NDOF*(idof-1)+jdof)=SAINVL(NDOF2*(j-1)+NDOF*(idof-1)+jdof)-dd*zz(NDOF*(in-1)+jdof)
                  end do 
                enddo
              endif
            end do 
          enddo
        end do
      enddo
      deallocate(vv)
      deallocate(zz)
  