mod_comm.f90

CCSM Research Tools: Community Atmosphere Model (CAM)

      module mod_comm
      implicit none

#include "misc.h"
#if defined ( SPMD ) 

#include "params.h"

      integer maxpro       ! Max no. of MLP PE allowed
      integer nbuf
      integer nghost
      integer max_nq       ! Be carefiul: max_nq = max(nc, 2)
                           ! nc is the total # of advected tracers

      parameter ( maxpro = PLAT/4 ) ! This is the max 1D decomp
      parameter ( nbuf = 2 )
      parameter ( nghost = 3 )
!     parameter ( max_nq = 2 )
      parameter ( max_nq = PCNST + 1 )

#if !defined(USE_MLP)
#include "mpif.h"
#define mp_precision MPI_DOUBLE_PRECISION

      integer max_call
      integer igosouth, igonorth
      integer idimsize
      parameter (max_call = 2)
      parameter (igosouth = 0)
      parameter (igonorth = 1)
      parameter (idimsize = PLON*nghost*PLEV*PCNST)
#if defined(AIX) && defined(MPI2)
      integer(kind=MPI_ADDRESS_KIND) intptr
      pointer (buff_r_ptr, buff_r(idimsize*nbuf*max_call))
      pointer (buff_s_ptr, buff_s(idimsize*nbuf*max_call))
      pointer (buff4d_ptr, buff4d(PLON*PLAT*(PLEV+1)*PCNST))
      pointer (buff4d_r4_ptr, buff4d_r4(PLON*PLAT*(PLEV+1)*PCNST))
      real :: buff_r
      real :: buff_s
      real :: buff4d
      real*4 :: buff4d_r4
#else
      real, SAVE:: buff_r(idimsize*nbuf*max_call)
      real, SAVE:: buff_s(idimsize*nbuf*max_call)
      real, SAVE:: buff4d(PLON*PLAT*(PLEV+1)*PCNST)
      real*4, SAVE:: buff4d_r4(PLON*PLAT*(PLEV+1)*PCNST)
#endif
      integer, SAVE:: ncall_s
      integer, SAVE:: ncall_r

#if defined(MPI2)
      integer(kind=MPI_ADDRESS_KIND) bsize, tdisp
      integer, SAVE:: buffwin      ! Communication window
      integer, SAVE:: buff4dwin    ! Communication window
      integer, SAVE:: buff4d_r4win ! Communication window
#else
      integer, SAVE:: tdisp
      integer, SAVE:: nsend                   ! Number of messages out-going
      integer, SAVE:: nrecv                   ! Number of messages in-coming
      integer, SAVE:: nread                   ! Number of messages read
      integer, SAVE:: sqest(nbuf*max_call)
      integer, SAVE:: rqest(nbuf*max_call)
#endif

      integer, SAVE:: commglobal   ! Global Communicator
      integer, SAVE:: Status(MPI_STATUS_SIZE)
      integer, SAVE:: Stats(nbuf*max_call*MPI_STATUS_SIZE)
      integer ierror

#else

#if defined (LAHEY)

#define PTR_INT TRUE
#define NOT_ASSIGNED

#include "mlp_ptr.h"

#undef  PTR_INT
#undef NOT_ASSIGNED

#else
!
! Main vars:
!
      pointer (wing_4d, g_4d)
      real :: g_4d(PLON, PLAT, PLEV, max_nq)

! Other work arrays:
!
! Type 1: For variables defined at layer edge (wz & pk)
!
      pointer (wing_t1, g_t1)
      real :: g_t1(PLON, PLAT, PLEV+1, nbuf)

!
! Type 2: For edge pressure (pe)
!
      pointer (wing_t2, g_t2)
      real :: g_t2(PLON, PLEV+1, PLAT)
!
! Type 3: 
!
      pointer (wing_t3, g_t3)
      real :: g_t3(PLEV+PLAT, maxpro)

!
! General purpose 2D (x-y) array
!
      pointer (wing_2d, g_2d)
      real :: g_2d(PLON,PLAT)

!
! General purpose 1D array
!
      pointer (wing_1d, g_1d)
      real :: g_1d(PLAT)
#endif
#endif

      integer, SAVE:: nowpro,numpro,numcps(maxpro),numcpu  
      integer, SAVE:: gid, gsize
      integer, allocatable, SAVE:: yfirst(:)  ! First latitude
      integer, allocatable, SAVE:: ylast(:)   ! Last latitude
      integer, allocatable, SAVE:: zfirst(:)  ! First level
      integer, allocatable, SAVE:: zlast(:)   ! Last level

      public mp_init, mp_exit, y_decomp, set_decomp
 
#if defined (SEMA)
      integer semid
#endif

!.................

      contains

      subroutine mp_init

#if !defined(USE_MLP)
#if !defined (SET_CPUS)
#if defined (IRIX64)
        integer mp_suggested_numthreads
#else
        integer omp_get_num_threads
#endif
#endif
        integer idimBuff, idimBuff4d
        integer n, nowpro, nowcpu
        integer npthreads
        character*80 evalue
        integer info
        logical flag

        integer mp_size

#if defined(MPI2) && !defined(AIX)
        call MPI_INIT_THREAD(MPI_THREAD_MULTIPLE, npthreads, ierror)
        call MPI_QUERY_THREAD(npthreads, ierror)
        if (npthreads == MPI_THREAD_SINGLE) then
          write(*,*) 'Provided MPI_THREAD_SINGLE on', gid
          call MPI_FINALIZE(ierror)
          stop
        elseif (npthreads == MPI_THREAD_FUNNELED) then
          write(*,*) 'Provided MPI_THREAD_FUNNELED on', gid
          call MPI_FINALIZE(ierror)
          stop
        elseif (npthreads == MPI_THREAD_SERIALIZED) then
          write(*,*) 'Provided MPI_THREAD_SERIALIZED on', gid
          call MPI_FINALIZE(ierror)
          stop
        elseif (npthreads == MPI_THREAD_MULTIPLE) then
!          write(*,*) 'Provided MPI_THREAD_MULTIPLE on', gid
        else
          write(*,*) gid,': Error in MPI_INIT_THREAD', npthreads, ':', ierror
          call MPI_FINALIZE(ierror)
          stop
        endif
#else
        call MPI_INITIALIZED( flag, ierror )
        if ( .not. flag ) then
          call MPI_INIT( ierror )
        endif
#endif
        call MPI_COMM_RANK (MPI_COMM_WORLD, gid, ierror)
        call MPI_COMM_SIZE (MPI_COMM_WORLD, numpro, ierror)
        call MPI_COMM_DUP (MPI_COMM_WORLD, commglobal, ierror)

#if defined(MPI2)
        call MPI_INFO_CREATE(info, ierror)
        call MPI_INFO_SET(info, "no_locks", "true", ierror)
#if defined(AIX)
        call MPI_TYPE_EXTENT(mp_precision, mp_size, ierror)
        bsize=idimsize*nbuf*max_call*mp_size
        call MPI_ALLOC_MEM(bsize, MPI_INFO_NULL, intptr, ierror)
        buff_r_ptr = intptr
        call MPI_WIN_CREATE(buff_r, bsize, mp_size, info, commglobal, &
                            buffwin, ierror)

        call MPI_ALLOC_MEM(bsize, MPI_INFO_NULL, intptr, ierror)
        buff_s_ptr = intptr

        bsize=PLON*PLAT*(PLEV+1)*PCNST*mp_size
        call MPI_ALLOC_MEM(bsize, MPI_INFO_NULL, intptr, ierror)
        buff4d_ptr = intptr
        call MPI_WIN_CREATE(buff4d, bsize, mp_size, info, commglobal, &
                            buff4dwin, ierror)

        call MPI_TYPE_EXTENT(MPI_REAL, mp_size, ierror)
        bsize=PLON*PLAT*(PLEV+1)*PCNST*mp_size
        call MPI_ALLOC_MEM(bsize, MPI_INFO_NULL, intptr, ierror)
        buff4d_r4_ptr = intptr
        call MPI_WIN_CREATE(buff4d_r4, bsize, mp_size, info, commglobal, &
                            buff4d_r4win, ierror)
#else
        call MPI_TYPE_EXTENT(mp_precision, mp_size, ierror)
        bsize=idimsize*nbuf*max_call
        call MPI_WIN_CREATE(buff_r, bsize, mp_size, info, commglobal, &
                            buffwin, ierror)

        bsize=PLON*PLAT*(PLEV+1)*PCNST
        call MPI_WIN_CREATE(buff4d, bsize, mp_size, info, commglobal, &
                            buff4dwin, ierror)

        bsize=PLON*PLAT*(PLEV+1)*PCNST
        call MPI_TYPE_EXTENT(MPI_REAL, mp_size, ierror)
        call MPI_WIN_CREATE(buff4d_r4, bsize, mp_size, info, commglobal, &
                            buff4d_r4win, ierror)
#endif
        call MPI_INFO_FREE(info, ierror)
#else
        nsend = 0
        nrecv = 0
        nread = 0
#endif

        ncall_r = 0
        ncall_s = 0

#if defined(SET_CPUS)
        call getenv('NUMBER_CPUS_PER_MLP_PROCESS',evalue)
        read(evalue,*) numcpu
#if defined (IRIX64)
       call mp_set_numthreads(numcpu)  !keep it for a while, :)
#else
       call omp_set_num_threads(numcpu)
#endif
#if  defined( IRIX64 ) && defined(PIN_CPUS)
!$omp parallel do private(n,nowcpu)
        nowpro = gid
        do n=1,numcpu
          nowcpu = n + (nowpro) * numcpu-1
          call mp_assign_to_cpu(nowcpu)
        enddo
#endif
#else
#if defined (IRIX64)
        numcpu = mp_suggested_numthreads(0)
#else
#if defined (_OPENMP)
        numcpu = omp_get_num_threads()
#else
        numcpu = 1
#endif
#endif
#endif
#else
        if ( max_nq < PCNST ) then
           write(*,*) "Buffer size for MLP is NOT large enough!"
           stop
        endif

        call gotmem
        call forkit
#if defined (SEMA)
        call semcreate(semid)
#endif

#endif
        allocate( yfirst( numpro ) )
        allocate( ylast( numpro ) )
        allocate( zfirst( numpro ) )
        allocate( zlast( numpro ) )
      end subroutine mp_init

      subroutine mp_exit
        deallocate( yfirst )
        deallocate( ylast )
        deallocate( zfirst )
        deallocate( zlast )
#if !defined(USE_MLP)
#if defined(MPI2)
        call MPI_WIN_FREE( buffwin, ierror )
        call MPI_WIN_FREE( buff4dwin, ierror )
        call MPI_WIN_FREE( buff4d_r4win, ierror )
#endif
        call MPI_FINALIZE (ierror)
#endif
        return
      end subroutine mp_exit


#if defined(USE_MLP)
      subroutine gotmem

#define NOT_ASSIGNED
#include "mlp_ptr.h"
#undef  NOT_ASSIGNED


      integer n_svar
      integer*8 numvar       ! Total # of shared vars     
      parameter (n_svar=100)
      integer*8 isize(n_svar),ipnt(n_svar)
      integer n

      numvar    =  6
      isize(1)  =  PLON*PLAT*PLEV*max_nq
      isize(2)  =  PLON*PLAT*(PLEV+1)*nbuf
      isize(3)  =  PLON*PLAT*(PLEV+1)
      isize(4)  = (PLEV+PLAT)*maxpro
      isize(5)  =  PLON*PLAT
      isize(6)  =  PLAT

      do n=1,numvar
         isize(n) = isize(n) * 8
      enddo

      call mlp_getmem(numvar,isize,ipnt)

      wing_4d  = ipnt(1)
      wing_t1  = ipnt(2)
      wing_t2  = ipnt(3)
      wing_t3  = ipnt(4)
      wing_2d  = ipnt(5)
      wing_1d  = ipnt(6)

#if defined (LAHEY)
      ptrg_4d  = wing_4d
      ptrg_t1  = wing_t1
      ptrg_t2  = wing_t2
      ptrg_t3  = wing_t3
      ptrg_2d  = wing_2d
      ptrg_1d  = wing_1d
#endif

      return
      end subroutine gotmem


      subroutine forkit

#if defined(IRIX64)
#include <ulocks.h>
#endif
      integer fork,getpid
      integer master, n, nowpid, ierror, nowcpu
      character*80 evalue

!-----create mp environment
      call getenv('NUMBER_MLP_PROCESSES',evalue)
      read(evalue,*) numpro
      call getenv('NUMBER_CPUS_PER_MLP_PROCESS',evalue)
      read(evalue,*) numcpu

!-----get master pid
      master = getpid()
      nowpro = 1

!-----print fork message
!!!      write(*,510) numpro

#if defined(IRIX64)
!-----destroy mp environment
      call mp_destroy
#endif

!-----spawn the processes - manual forks
      do n=2,numpro
         nowpid = getpid()
         if(nowpid == master) then 
                              ierror=fork()
                              endif
         nowpid = getpid()
         if(nowpid /= master) then
                              nowpro=n
                              go to 200
                              endif
      enddo

!-----write note
  200 if(nowpro == 1) nowpid = master
!!!     write(*,500) nowpro,nowpid

      call omp_set_num_threads(numcpu)

#if  defined( IRIX64 ) && defined(PIN_CPUS)

!$omp parallel do private(n,nowcpu)
      do n=1,numcpu
         nowcpu = n+(nowpro-1)*numcpu-1
         call mp_assign_to_cpu(nowcpu)
      enddo
#endif

!******************************
!*    I/O formats             *
!******************************
  500 format('FORKIT: Current process:',i3,'    PID:',i10)
  510 format('FORKIT: Total active processes spawned:',i3)

      gid = nowpro-1
      gsize = numpro

      return
      end subroutine forkit
#endif

      subroutine y_decomp(jm, km, jfirst, jlast, kfirst, klast, myid)

      implicit none
      integer jm     ! Dimensions
      integer km     ! Levels
      integer myid

! OUTPUT PARAMETERS:
      integer jfirst, jlast, kfirst, klast

! Local
      integer p, p1, p2, lats, pleft
      integer, allocatable:: ydist(:)

      if (myid == 0) print *, "numpro", numpro, "numcpu", numcpu
      allocate( ydist( numpro ) )

      lats = jm / numpro
      pleft = jm - lats * numpro

      if( lats < 3 ) then
         write(*,*) 'Number of Proc is too large for jm=',jm
         stop
      endif

      do p=1,numpro
         ydist(p) = lats
      enddo

      if ( pleft .ne. 0 ) then
          p1 = (numpro+1) / 2 
          p2 = p1 + 1
        do while ( pleft .ne. 0 )
           if( p1 .eq. 1 ) p1 = numpro
               ydist(p1) = ydist(p1) + 1
               pleft = pleft - 1
               if ( pleft .ne. 0 ) then
                    ydist(p2) = ydist(p2) + 1
                    pleft = pleft - 1