anchor.diff:         END DO
         END DO
   END SUBROUTINE
  +
  +
  +!Calculate the matrices B, C and D
  +!  N        Figure of elements            (IN)
  +!  XNODE    Nodes vector                  (IN)
  +!  NW       Number of wake elements       (IN)
  +!  XWNODE   Wake's nodes vector           (IN)
  +!  B        Matrix B                      (OUT)
  +!  C        Matrix C                      (OUT)
  +!  D        Matrix D (wake doublet int.)  (OUT)
  +subroutine SrfMatBCD(N, Xnode, NW, XWnode, B, C, D)
  +      IMPLICIT NONE
  +      integer, intent(IN) :: N, NW
  +      real(kind=8), dimension(N,2), intent(IN) :: Xnode
  +      real(kind=8), dimension(NW,2), intent(IN) :: XWnode
  +      real(kind=8), dimension(N,N), intent(OUT) :: B, C
  +      real(kind=8), dimension(N,NW), intent(OUT) :: D
  +      real(kind=8), dimension(2) :: XI
  +      real(kind=8) :: CIJ, BIJ, DIJ
  +      integer :: i, j, NOE
  +      do i = 1, N
  +       !Midpoint of the i-th element
  +       XI = (Xnode(NOE(N,i,0),:) + Xnode(NOE(N,i,1),:))/2.
  +       C(i,i) = 0.
  +       do j = 1, N
  +        call BCIJ(XI, Xnode(NOE(N,j,0),:), Xnode(NOE(N,j,1),:),BIJ,CIJ)
  +        B(i,j) = BIJ
  +        if (i .ne. j) then
  +          C(i,j) = CIJ
  +        end if
  +       end do
  +       do j = 1, NW-1
  +        D(i,j) = DIJ(XI, XWnode(j,:), XWnode(j+1,:))
  +       end do
  +      end do
  +end subroutine
   
   
   !Calculate B and C matrices on the field points X
  
positive.code: 
  !Calculate the matrix B and C
  !  N        Figure of elements  (IN)
  !  B        Matrix B            (OUT)
  !  C        Matrix C            (OUT)
  !  XNODE    Nodes vector        (IN)
  SUBROUTINE SRFMATBC(N, XNODE, B, C)
        IMPLICIT NONE
        INTEGER, INTENT(IN) :: N
  
negative.code:       INTEGER :: I, J, NOE
        DO I = 1, N
         !Midpoint of I element
         XI = (XNODE(NOE(N,I,0),:) + XNODE(NOE(N,I,1),:))/2.
         C(I,I) = 0.
         DO J = 1, N
          CALL BCIJ(XI, XNODE(NOE(N,J,0),:), XNODE(NOE(N,J,1),:),BIJ,CIJ)
          B(I,J) = BIJ
          IF (I .NE. J) THEN
            C(I,J) = CIJ
  