anchor.diff:      double precision,intent(in)::v(:,:) !< transported variables
       double precision,intent(in)::mFlow(:) !< mass flow through pairs
       double precision,allocatable,intent(inout)::flow(:,:) !< v flow rate output
  +    double precision,allocatable::gradV(:,:,:)
  +    double precision::vf(size(v,1)),dV(size(v,1)),r(size(v,1))
       integer::up,dn
       
       call grid%up()
  
positive.code:             up=n
              dn=m
            end if
          end if
          flow(:,i)=v(:,up)*mFlow(i)
        end if
      end do
      !$omp end parallel do
    end subroutine
    
    !> find flow rate of scalar v on polyFvGrid using mass flow rate
  
negative.code:       allocate(flow(size(v,1),grid%nP))
      end if
      flow(:,:)=0d0
      allocate(gradV(DIMS,size(v,1),grid%nC))
      call findGrad(grid,v(:,1:grid%nC),gradV)
      !$omp parallel do default(shared)&
      !$omp& private(m,n,up,dn,dV,r,vf)
      do i=1,grid%nP
        m=grid%iEP(1,i)
        n=grid%iEP(2,i)
  