cd_core.f90

CCSM Research Tools: Community Atmosphere Model (CAM)

              yzt(i,j,k) = pt(i,j,k)
            enddo
          enddo
        enddo
        call redistributestart(inter_ijk, .true., yzt)
        call redistributefinish(inter_ijk, .true., ptxy)

        call geopk(ptop, pexy, delpxy, pkxy, wzxy, hsxy, ptxy, 0, im, jm, km,    &
                   jfirstxy, jlastxy, ifirstxy, ilastxy,                    &
                   cp, akap, nx, ipe)

!
! Transpose back to yz decomposition
!

        call redistributestart(inter_ijkp, .false., pkxy)
        call redistributefinish(inter_ijkp, .false., pkkp)
        call redistributestart(inter_ijkp, .false., wzxy)
        call redistributefinish(inter_ijkp, .false., wzkp)

!$omp parallel do private(i, j, k)
        do k = kfirst, klastp
        do j = jfirst, jlast
        do i = 1, im
          pkc(i,j,k) = pkkp(i,j,k)
          wz(i,j,k) = wzkp(i,j,k)
        enddo
        enddo
        enddo

!
! pexy is an output quantity from geopk if ipe=1
! Save pexy and pkxy to avoid main transpose; pe not needed.
!

        if (npr_z .gt. 1) then
!
! Fill in klast+1
!
         incount = 0
         outcount = 0
         if (kfirst .gt. 1) then
          call bufferpack3d(pkc, 1, im, jfirst-1, jlast+1, kfirst, klast+1,    &
                            1, im, jfirst, jlast, kfirst, kfirst,              &
                            buff_s) 
          incount = im * (jlast-jfirst+1)
          call bufferpack3d(wz, 1, im, jfirst-1, jlast+1, kfirst, klast+1,     &
                            1, im, jfirst, jlast, kfirst, kfirst,              &
                            buff_s(incount+1)) 
          incount = incount + im * (jlast-jfirst+1)
         endif
         if (klast .lt. km) then
          outcount = 2 * im * (jlast - jfirst + 1)
         endif
         call mp_barrier()
         call mp_send(iam-npr_y,iam+npr_y,incount,outcount,buff_s,buff_r)
         call mp_barrier()
         call mp_recv(iam+npr_y,outcount,buff_r)
         if (klast .lt. km) then
          call bufferunpack3d(pkc,1,im,jfirst-1,jlast+1,kfirst,klast+1,     &
                              1, im, jfirst, jlast, klast+1, klast+1,       &
                              buff_r )
          outcount = im * (jlast - jfirst + 1)
          call bufferunpack3d(wz,1,im,jfirst-1,jlast+1,kfirst,klast+1,      &
                              1, im, jfirst, jlast, klast+1, klast+1,       &
                              buff_r(outcount+1) )
         endif
        endif       ! npr_z .gt. 1

       endif       ! geopkd16

      else
#endif

        if (geopkd16) then

!
! Use 16-byte reals (for compatibility with 2-D decomposition)
!

          call geopk16(ptop, pe, delp, pkcc, wzc, hs, pt, ng_d, im, jm, km,          &
                       jfirst, jlast, 1, im, cp, akap, 1, km)

        else

!
! Use 8-byte reals (standard)
!

          call geopk(ptop, pe, delp, pkcc, wzc, hs, pt, ng_d, im, jm, km,            &
                     jfirst, jlast, 1, im,                                     &
                     cp, akap, nx, ipe)

        endif

!
! Geopk does not need j ghost zones of pkc and wz
!

!$omp parallel do private(i, j, k)
        do k = kfirst, klast+1
          do j = jfirst, jlast
            do i = 1, im
              pkc(i,j,k) = pkcc(i,j,k)
              wz(i,j,k) = wzc(i,j,k)
            enddo
          enddo
        enddo

#if defined( SPMD )
      endif
#endif
!
! Upon exit from geopk, the quantities pe, pkc and wz will have been
!      updated at klast+1

      call t_stopf('d_geop')


#if defined( SPMD )
! Exchange boundary regions on north and south for pkc and wz
      call t_startf('d_mpi_1')
      call mp_send2_ns(im, jm, jfirst, jlast, kfirst, klast+1, &
                       1, pkc, wz)
      call t_stopf('d_mpi_1')
#endif

   if ( ipe .ne. 1 ) then          !  not the last call
!
! Perform some work while sending data on the way
!

!$omp parallel do private(i, j, k)

    do k=kfirst,klast
         do j=jfirst,jlast
            do i=1,im
               delpf(i,j,k) = delp(i,j,k)
             enddo
          enddo
          call pft2d( delpf(1,js2g0,k), sc(js2g0), dc(1,js2g0),  &
                      im, jn2g0-js2g0+1, ifax, trigs )
    enddo

   else

!$omp parallel do private(i, j, k)

    do k=kfirst,klast+1
        do j=jfirst,jlast
           do i=1,im
             pk(i,j,k) = pkc(i,j,k)
           enddo
         enddo
    enddo

   endif

#if defined( SPMD )
      call t_startf('d_mpi_2')
      call mp_recv2_ns(im, jm, jfirst, jlast, kfirst, klast+1, &
                       1, pkc, wz)
      call t_stopf('d_mpi_2')
#endif


#if defined ( HIGH_P )

      call t_startf('highp')

      call highp2(pkc,  wz,  wz3,    wzz,   dpt,          &
                  im,   jm,   km,  jfirst, jlast,         &
                  kfirst, klast, klastp, nx )

      call t_stopf('highp')

#else

!
! Beware k+1 references directly below (AAM)
!
!$omp parallel do private(i, j, k)

   do k=kfirst,klast
      do j=js1g1,jn1g1                  ! dpt needed NS
         do i=1,im                       ! wz, pkc ghosted NS
            dpt(i,j,k)=(wz(i,j,k+1)+wz(i,j,k))*(pkc(i,j,k+1)-pkc(i,j,k))
         enddo
      enddo
   enddo

#endif


#if defined( SPMD )

      call t_startf('d_mpi_3')

      if ( ipe .ne. 1 ) then          !  not the last call
         call mp_send3d_ns(im, jm, jfirst, jlast, kfirst, klast, &
                            ng_d, ng_d, delpf, 1 )
      endif

      call t_stopf('d_mpi_3')

#endif

 
!  GHOSTING:   wz (input) NS ; pkc (input) NS

      call t_startf('d-4500')

!$omp parallel do private(i, j, k, wk3, wk1)

#if defined ( HIGH_P )
   do 4500 k=kfirst,klast+1
#else
   do 4500 k=max(kfirst,2),klast+1
#endif

      do j=js2g1,jn2g0                             ! wk3 needed S

#if defined ( HIGH_P )
          do i=1,im
             wk3(i,j) = wz3(i,j,k) 
          enddo
#else
            wk3(1,j) = (wz(1,j,k)+wz(im,j,k)) *       &
                     (pkc(1,j,k)-pkc(im,j,k))
          do i=2,im
            wk3(i,j) = (wz(i,j,k)+wz(i-1,j,k)) *      & 
                       (pkc(i,j,k)-pkc(i-1,j,k))
          enddo
#endif
      enddo

      do j=js2g1,jn2g0                               
         do i=1,im-1
            wk1(i,j) = wk3(i,j) + wk3(i+1,j)        
         enddo
            wk1(im,j) = wk3(im,j) + wk3(1,j)      ! wk3 ghosted S
      enddo

      if ( jfirst == 1 ) then
          do i=1,im
            wk1(i, 1) = 0.
          enddo
      endif

      if ( jlast == jm ) then
          do i=1,im
            wk1(i,jm) = 0.
          enddo
      endif

      do j=js2g0,jlast                          ! wk1 ghosted S
         do i=1,im
            wz3(i,j,k) = wk1(i,j) + wk1(i,j-1)
         enddo
      enddo

! N-S walls

      do j=js2g0,jn1g1                     ! wk1 needed N

#if defined ( HIGH_P )
          do i=1,im                          ! wz, pkc ghosted NS
             wk1(i,j) = wzz(i,j,k)
          enddo
#else
          do i=1,im                         ! wz, pkc ghosted NS
             wk1(i,j) = (wz(i,j,k)+wz(i,j-1,k))*(pkc(i,j,k)-pkc(i,j-1,k))
          enddo
#endif
      enddo

      do j=js2g0,jn1g1                    ! wk3 needed N
            wk3(1,j) = wk1(1,j) + wk1(im,j)      ! wk1 ghosted N
          do i=2,im
            wk3(i,j) = wk1(i,j) + wk1(i-1,j)   ! wk1 ghosted N
          enddo
      enddo

      do j=js2g0,jn2g0
         do i=1,im
            wz(i,j,k) = wk3(i,j) + wk3(i,j+1)  ! wk3 ghosted N
         enddo
      enddo
!!!      call avgc( pkc(1,jfirst-1,k), pkc(1,jfirst,k), im, jm,  &
!!!                 jfirst, jlast, wk1 )
      do j=js1g1,jn1g1
        wk1(1,j) = pkc(1,j,k) + pkc(im,j,k)
        do i=2,im
          wk1(i,j) = pkc(i,j,k) + pkc(i-1,j,k)
        enddo
      enddo

      do j=js2g0,jn1g1
        do i=1,im
          pkc(i,j,k) = wk1(i,j) + wk1(i,j-1)
        enddo
      enddo

4500  continue

      call t_stopf('d-4500')

#if ( !defined HIGH_P )
      if (kfirst == 1) then
        do j=js2g0,jlast
          do i=1,im
            wz3(i,j,1) = 0.
            wz(i,j,1) = 0.
          enddo
        enddo

        pk4 = 4.*ptop**akap
        do j=js2g0,jn1g1
          do i=1,im
            pkc(i,j,1) = pk4
          enddo
        enddo
      endif
#endif

!     GHOSTING:   dpt (loop 4000) NS ; pkc (loop 4500) N
      call t_startf('d-6000')

!
! Beware k+1 references directly below (AAM)
!
!$omp parallel do private(i, j, k, ub, wk1, wk2, wk3)

   do 6000 k=kfirst,klast

!!!      call avgc(dpt(1,jfirst-1,k), wk2(1,jfirst), im, jm,          &
!!!                jfirst, jlast, wk1)

      do j=js1g1,jn1g1
        wk1(1,j) = dpt(1,j,k) + dpt(im,j,k)
        do i=2,im
          wk1(i,j) = dpt(i,j,k) + dpt(i-1,j,k)
        enddo
      enddo
 
      do j=js2g0,jn1g1
        do i=1,im
          wk2(i,j) = wk1(i,j) + wk1(i,j-1)
        enddo
      enddo

      do j=js2g0,jn1g1
         do i=1,im
            ub(i,j) = pkc(i,j,k+1) - pkc(i,j,k)
         enddo
      enddo

      do j=js2g0,jlast
         do i=1,im-1
            wk3(i,j) = uc(i,j,k)  +  dtdxe(j)/(ub(i,j) + ub(i+1,j))      &
                       * (wk2(i,j)-wk2(i+1,j)+wz3(i,j,k+1)-wz3(i,j,k))
         enddo
          wk3(im,j) = uc(im,j,k)  +  dtdxe(j)/(ub(im,j) + ub(1,j))       &
                      * (wk2(im,j)-wk2(1,j)+wz3(im,j,k+1)-wz3(im,j,k))
      enddo

      do j=js2g0,jn2g0                  ! Assumes wk2 ghosted on N
         do i=1,im
            wk1(i,j) = vc(i,j,k)  +  dtdy/(ub(i,j)+ub(i,j+1)) *          &
                       (wk2(i,j)-wk2(i,j+1)+wz(i,j,k+1)-wz(i,j,k))
         enddo
      enddo

#if ( !defined ALT_PFT )
      call pft2d( wk3(1,js2g0), se(js2g0), de(1,js2g0), im,           &
                  jlast-js2g0+1,  ifax, trigs )
      call pft2d( wk1(1,js2g0), sc(js2g0), dc(1,js2g0), im,           &
                  jn2g0-js2g0+1,  ifax, trigs )
#endif

      do j=js2g0,jlast
         do i=1,im
            u(i,j,k) = u(i,j,k) + wk3(i,j)
         enddo
      enddo

      do j=js2g0,jn2g0
         do i=1,im
            v(i,j,k) = v(i,j,k) + wk1(i,j)
         enddo
      enddo

#if defined ( ALT_PFT )
      call pft2d( u(1,js2g0,k), se(js2g0), de(1,js2g0), im,    &
                  jlast-js2g0+1,  ifax, trigs )
      call pft2d( v(1,js2g0,k), sc(js2g0), dc(1,js2g0), im,    &
                  jn2g0-js2g0+1,  ifax, trigs )
#endif

6000  continue

      call t_stopf('d-6000')

#if defined( SPMD )

      call t_startf('d_mpi_4')

      if ( ipe .ne. 1 ) then
         call mp_recv3d_ns(im, jm, jfirst, jlast, kfirst, klast, &
                           ng_d, ng_d, delpf, 1)
      endif

      call t_stopf('d_mpi_4')

#endif

      return
!EOC
 end
!-----------------------------------------------------------------------