
c---------------------------------------------------------------------
c---------------------------------------------------------------------

      subroutine pintgr

c---------------------------------------------------------------------
c---------------------------------------------------------------------

      implicit none

      include 'mpinpb.h'
      include 'applu.incl'

c---------------------------------------------------------------------
c  local variables
c---------------------------------------------------------------------
      integer i, j, k
      integer ibeg, ifin, ifin1
      integer jbeg, jfin, jfin1
      integer iglob, iglob1, iglob2
      integer jglob, jglob1, jglob2
      integer ind1, ind2
      double precision  phi1(0:isiz2+1,0:isiz3+1),
     >                  phi2(0:isiz2+1,0:isiz3+1)
      double precision  frc1, frc2, frc3
      double precision  dummy

      integer IERROR


c---------------------------------------------------------------------
c   set up the sub-domains for integeration in each processor
c---------------------------------------------------------------------
      ibeg = nx + 1
      ifin = 0
      iglob1 = ipt + 1
      iglob2 = ipt + nx
      if (iglob1.ge.ii1.and.iglob2.lt.ii2+nx) ibeg = 1
      if (iglob1.gt.ii1-nx.and.iglob2.le.ii2) ifin = nx
      if (ii1.ge.iglob1.and.ii1.le.iglob2) ibeg = ii1 - ipt
      if (ii2.ge.iglob1.and.ii2.le.iglob2) ifin = ii2 - ipt
      jbeg = ny + 1
      jfin = 0
      jglob1 = jpt + 1
      jglob2 = jpt + ny
      if (jglob1.ge.ji1.and.jglob2.lt.ji2+ny) jbeg = 1
      if (jglob1.gt.ji1-ny.and.jglob2.le.ji2) jfin = ny
      if (ji1.ge.jglob1.and.ji1.le.jglob2) jbeg = ji1 - jpt
      if (ji2.ge.jglob1.and.ji2.le.jglob2) jfin = ji2 - jpt
      ifin1 = ifin
      jfin1 = jfin
      if (ipt + ifin1.eq.ii2) ifin1 = ifin -1
      if (jpt + jfin1.eq.ji2) jfin1 = jfin -1

c---------------------------------------------------------------------
c   initialize
c---------------------------------------------------------------------
      do i = 0,isiz2+1
        do k = 0,isiz3+1
          phi1(i,k) = 0.
          phi2(i,k) = 0.
        end do
      end do

      do j = jbeg,jfin
         jglob = jpt + j
         do i = ibeg,ifin
            iglob = ipt + i

            k = ki1

            phi1(i,j) = c2*(  u(5,i,j,k)
     >           - 0.50d+00 * (  u(2,i,j,k) ** 2
     >                         + u(3,i,j,k) ** 2
     >                         + u(4,i,j,k) ** 2 )
     >                        / u(1,i,j,k) )

            k = ki2

            phi2(i,j) = c2*(  u(5,i,j,k)
     >           - 0.50d+00 * (  u(2,i,j,k) ** 2
     >                         + u(3,i,j,k) ** 2
     >                         + u(4,i,j,k) ** 2 )
     >                        / u(1,i,j,k) )
         end do
      end do

c---------------------------------------------------------------------
c  communicate in i and j directions
c---------------------------------------------------------------------
      call exchange_4(phi1,phi2,ibeg,ifin1,jbeg,jfin1)

      frc1 = 0.0d+00

      do j = jbeg,jfin1
         do i = ibeg, ifin1
            frc1 = frc1 + (  phi1(i,j)
     >                     + phi1(i+1,j)
     >                     + phi1(i,j+1)
     >                     + phi1(i+1,j+1)
     >                     + phi2(i,j)
     >                     + phi2(i+1,j)
     >                     + phi2(i,j+1)
     >                     + phi2(i+1,j+1) )
         end do
      end do

c---------------------------------------------------------------------
c  compute the global sum of individual contributions to frc1
c---------------------------------------------------------------------
      dummy = frc1
      call MPI_ALLREDUCE( dummy,
     >                    frc1,
     >                    1,
     >                    dp_type,
     >                    MPI_SUM,
     >                    MPI_COMM_WORLD,
     >                    IERROR )

      frc1 = dxi * deta * frc1

c---------------------------------------------------------------------
c   initialize
c---------------------------------------------------------------------
      do i = 0,isiz2+1
        do k = 0,isiz3+1
          phi1(i,k) = 0.
          phi2(i,k) = 0.
        end do
      end do
      jglob = jpt + jbeg
      ind1 = 0
      if (jglob.eq.ji1) then
        ind1 = 1
        do k = ki1, ki2
           do i = ibeg, ifin
              iglob = ipt + i
              phi1(i,k) = c2*(  u(5,i,jbeg,k)
     >             - 0.50d+00 * (  u(2,i,jbeg,k) ** 2
     >                           + u(3,i,jbeg,k) ** 2
     >                           + u(4,i,jbeg,k) ** 2 )
     >                          / u(1,i,jbeg,k) )
           end do
        end do
      end if

      jglob = jpt + jfin
      ind2 = 0
      if (jglob.eq.ji2) then
        ind2 = 1
        do k = ki1, ki2
           do i = ibeg, ifin
              iglob = ipt + i
              phi2(i,k) = c2*(  u(5,i,jfin,k)
     >             - 0.50d+00 * (  u(2,i,jfin,k) ** 2
     >                           + u(3,i,jfin,k) ** 2
     >                           + u(4,i,jfin,k) ** 2 )
     >                          / u(1,i,jfin,k) )
           end do
        end do
      end if

c---------------------------------------------------------------------
c  communicate in i direction
c---------------------------------------------------------------------
      if (ind1.eq.1) then
        call exchange_5(phi1,ibeg,ifin1)
      end if
      if (ind2.eq.1) then
        call exchange_5 (phi2,ibeg,ifin1)
      end if

      frc2 = 0.0d+00
      do k = ki1, ki2-1
         do i = ibeg, ifin1
            frc2 = frc2 + (  phi1(i,k)
     >                     + phi1(i+1,k)
     >                     + phi1(i,k+1)
     >                     + phi1(i+1,k+1)
     >                     + phi2(i,k)
     >                     + phi2(i+1,k)
     >                     + phi2(i,k+1)
     >                     + phi2(i+1,k+1) )
         end do
      end do

c---------------------------------------------------------------------
c  compute the global sum of individual contributions to frc2
c---------------------------------------------------------------------
      dummy = frc2
      call MPI_ALLREDUCE( dummy,
     >                    frc2,
     >                    1,
     >                    dp_type,
     >                    MPI_SUM,
     >                    MPI_COMM_WORLD,
     >                    IERROR )

      frc2 = dxi * dzeta * frc2

c---------------------------------------------------------------------
c   initialize
c---------------------------------------------------------------------
      do i = 0,isiz2+1
        do k = 0,isiz3+1
          phi1(i,k) = 0.
          phi2(i,k) = 0.
        end do
      end do
      iglob = ipt + ibeg
      ind1 = 0
      if (iglob.eq.ii1) then
        ind1 = 1
        do k = ki1, ki2
           do j = jbeg, jfin
              jglob = jpt + j
              phi1(j,k) = c2*(  u(5,ibeg,j,k)
     >             - 0.50d+00 * (  u(2,ibeg,j,k) ** 2
     >                           + u(3,ibeg,j,k) ** 2
     >                           + u(4,ibeg,j,k) ** 2 )
     >                          / u(1,ibeg,j,k) )
           end do
        end do
      end if

      iglob = ipt + ifin
      ind2 = 0
      if (iglob.eq.ii2) then
        ind2 = 1
        do k = ki1, ki2
           do j = jbeg, jfin
              jglob = jpt + j
              phi2(j,k) = c2*(  u(5,ifin,j,k)
     >             - 0.50d+00 * (  u(2,ifin,j,k) ** 2
     >                           + u(3,ifin,j,k) ** 2
     >                           + u(4,ifin,j,k) ** 2 )
     >                          / u(1,ifin,j,k) )
           end do
        end do
      end if

c---------------------------------------------------------------------
c  communicate in j direction
c---------------------------------------------------------------------
      if (ind1.eq.1) then
        call exchange_6(phi1,jbeg,jfin1)
      end if
      if (ind2.eq.1) then
        call exchange_6(phi2,jbeg,jfin1)
      end if

      frc3 = 0.0d+00

      do k = ki1, ki2-1
         do j = jbeg, jfin1
            frc3 = frc3 + (  phi1(j,k)
     >                     + phi1(j+1,k)
     >                     + phi1(j,k+1)
     >                     + phi1(j+1,k+1)
     >                     + phi2(j,k)
     >                     + phi2(j+1,k)
     >                     + phi2(j,k+1)
     >                     + phi2(j+1,k+1) )
         end do
      end do

c---------------------------------------------------------------------
c  compute the global sum of individual contributions to frc3
c---------------------------------------------------------------------
      dummy = frc3
      call MPI_ALLREDUCE( dummy,
     >                    frc3,
     >                    1,
     >                    dp_type,
     >                    MPI_SUM,
     >                    MPI_COMM_WORLD,
     >                    IERROR )

      frc3 = deta * dzeta * frc3
      frc = 0.25d+00 * ( frc1 + frc2 + frc3 )
c      if (id.eq.0) write (*,1001) frc

      return

 1001 format (//5x,'surface integral = ',1pe12.5//)

      end