anchor.diff:      if(c_associated(dat))then
       end if
     end function
  +  
  +  !> solve the pressure, such that mass is conserved
  +  subroutine solvePressure()
  +    use modAdvection
  +    double precision,allocatable,save::tmpP(:)
  +    integer::info
  +    
  +    if(.not.allocated(tmpP))then
  +      allocate(tmpP(grid%nC))
  +    else if(size(tmpP)/=grid%nC)then
  +      deallocate(tmpP)
  +      allocate(tmpP(grid%nC))
  +    end if
  +    tmpP(:)=p(1:grid%nC)
  +    call setBC()
  +    call findMassFlow(grid,rhou,flowRho1)
  +    call pressureEq%solve(tmpP,info=info)
  +    needRetry=info<0
  +    call pressureEq%getNIt(n)
  +    nItPressure=max(nItPressure,n)
  +    p(1:grid%nC)=tmpP(:)
  +  end subroutine
     
     !> residual function of the pressure equation
     function pressureRHS(oldVector,newVector,dat)
  
positive.code:     end if
    end function
    
    !> RHS of the density equation in the form of a fix point problem
    function densityRHS(oldVector,newVector,dat)
      use iso_c_binding
  
negative.code:     end forall
    end subroutine
    
    !> RHS of the momentum equation in the form of a fix point problem
    function momentumRHS(oldVector,newVector,dat)
      use iso_c_binding
  