call pkez(nxu, im, km, jfirst, jlast, 1, km, ifirst, ilast,  &                pe, pk, akap, ks, peln, pkz, .true.)! ((((((((((((((((( compute globally integrated TE >>>>>>>>>>>>>>>>      if( consv ) then!$omp  parallel do         &!$omp  default(shared)     &!$omp  private(i,j,k)        do k=1,km          do j=jfirst,jlast             do i=ifirst,ilast                dz(i,j,k) = te(i,j,k) * delp(i,j,k)             enddo          enddo        enddo!$omp  parallel do        &!$omp  default(shared)    &!$omp  private(i,j,k,bte)! Perform vertical integration        do 4000 j=jfirst,jlast          if ( j == 1 ) then! SP            tte(1) = 0.              do k=1,km              tte(1) = tte(1) + dz(ifirst,1,k)            enddo          elseif ( j .eq. jm) then! NP            tte(jm) = 0.            do k=1,km              tte(jm) = tte(jm) + dz(ifirst,jm,k)            enddo          else! Interior            do i=ifirst,ilast              bte(i) = 0.            enddo            do k=1,km              do i=ifirst,ilast                bte(i) = bte(i) + dz(i,j,k)              enddo            enddo            xysum(j,1) = 0.            do i=ifirst,ilast              xysum(j,1) = xysum(j,1) + bte(i)#if !defined (OLDWAY)              tmpij(i,j,1) = bte(i)#endif            enddo          endif4000    continue#if defined (SPMD)        if (nprxy_x .gt. 1) then# if defined (OLDWAY)          call parcollective(commxy_x, sumop, jlast-jfirst+1, 1, xysum)# else          call par_xsum(tmpij, ifirst, ilast, im, jlast-jfirst+1, xysum)# endif        endif#endif!$omp  parallel do        &!$omp  default(shared)    &!$omp  private(j)        do j = max(jfirst,2), min(jlast,jm-1)           tte(j) = xysum(j,1)*cosp(j)        enddo        if ( jfirst == 1 ) tte(1)  = acap * tte(1)        if ( jlast == jm ) tte(jm) = acap * tte(jm)        te1 = 0.        call par_vecsum(jm, jfirst, jlast, tte, te1, comm_use, npry_use)      endif   ! consv#if defined( SPMD )      incount  = 0      outcount = 0! Send u southward      if ( jfirst > 1 ) then        call bufferpack3d( u, ifirst,ilast,jfirst-ngus,jlast+ngun,1,km,       &                           ifirst,ilast,jfirst,jfirst,1,km,buff_s )        incount  = itot*km      endif      if ( jlast < jm ) then        outcount = itot*km      endif      call mp_barrier()      call mp_send(iam-nprxy_x, iam+nprxy_x, incount, outcount, buff_s, buff_r)      call mp_barrier()      call mp_recv(iam+nprxy_x, outcount, buff_r)      if ( jlast < jm ) then        call bufferunpack3d( u,ifirst,ilast,jfirst-ngus,jlast+ngun,1,km,    &                                ifirst,ilast,jlast+1,jlast+1,1,km,buff_r )      endif#endif      if( consv ) then!$omp  parallel do       &!$omp& default(shared)   &!$omp& private(i,j)        do j=js2g0, jn2g0          xysum(j,1) = 0.          xysum(j,2) = 0.        do i=ifirst,ilast          xysum(j,1) = xysum(j,1) + ps(i,j)          xysum(j,2) = xysum(j,2) + peln(i,km+1,j)#if !defined (OLDWAY)          tmpij(i,j,1) = ps(i,j)          tmpij(i,j,2) = peln(i,km+1,j) #endif        enddo       enddo#if defined( SPMD )       if (nprxy_x .gt. 1) then# if defined (OLDWAY)          call parcollective(commxy_x, sumop, jlast-jfirst+1, 2, xysum)# else          call par_xsum(tmpij, ifirst, ilast, im, 2*(jlast-jfirst+1), xysum)# endif       endif#endif!$omp  parallel do       &!$omp  default(shared)   &!$omp  private(j)        do j=js2g0, jn2g0        tte(j) = cp*cosp(j)*(xysum(j,1) - ptop*float(im) -           &                 akap*ptop*(xysum(j,2) - peln(ifirst,1,j)*float(im)) )! peln(i,1,j) should be independent of i (AAM)       enddo       if ( jfirst .eq. 1 ) tte(1) = acap*cp * (ps(ifirst,1) - 2.*ptop -    &                akap*ptop*(peln(ifirst,km+1,1) - peln(ifirst,1,1) ) )       if ( jlast .eq. jm ) tte(jm)= acap*cp * (ps(ifirst,jm) -             &                akap*ptop*(peln(ifirst,km+1,jm) - peln(ifirst,1,jm) ) )      endif ! consv      if (consv) then       sum=0.       call par_vecsum(jm, jfirst, jlast, tte, sum, comm_use, npry_use)       dtmp = (te0 - te1) / sum       if( diag ) then         CPP_PRT_PREFIX write(6,*) 'te=',te0, ' Energy deficit in T = ', dtmp       endif      endif              ! end consv check! Single x-subdomain case (periodic)      do k = 1, km      do j = jfirst, jlast        veast(j,k) = v(ifirst,j,k)      enddo      enddo! Nontrivial x-decomposition#if defined( SPMD )      if (itot .ne. im) then        call bufferpack3d(v,ifirst,ilast,jfirst-ngvs,jlast+ngvn,1,km,    &         ifirst,ifirst,jfirst,jlast,1,km,buff_s)        dest = myidxy_y*nprxy_x + MOD(iam+nprxy_x-1,nprxy_x)        src  = myidxy_y*nprxy_x + MOD(iam+1,nprxy_x)        call mp_barrier()        call mp_send(dest, src, km*(jlast-jfirst+1), km*(jlast-jfirst+1),   &                     buff_s, buff_r)        call mp_barrier()        call mp_recv(src, km*(jlast-jfirst+1), buff_r)        call bufferunpack2d(veast, jfirst, jlast, 1, km,      &                            jfirst, jlast, 1, km, buff_r)      endif#endif!$omp  parallel do       &!$omp  default(shared)   &!$omp  private(i,j,k, u2, v2)      do 8000 k=1,km! Compute KE        do j=js2g0,jn1g1          do i=ifirst,ilast            u2(i,j) = u(i,j,k)**2          enddo        enddo        do j=js2g0,jn2g0          do i=ifirst,ilast            v2(i,j) = v(i,j,k)**2          enddo          v2(ilast+1,j) = veast(j,k)**2        enddo        do j=js2g0,jn2g0          do i=ifirst,ilast            te(i,j,k) = te(i,j,k) - 0.25 * ( u2(i,j) + u2(i,j+1)     &                                            +v2(i,j) + v2(i+1,j) )          enddo        enddo        if ( jfirst .eq. 1 ) then! South pole          do i=ifirst,ilast            u2_sp(i,k) = u2(i,2)            v2_sp(i,k) = v2(i,2)          enddo        endif        if ( jlast .eq. jm ) then! North pole          do i=ifirst,ilast            u2_np(i,k) = u2(i,jm)            v2_np(i,k) = v2(i,jm-1)          enddo        endif8000  continue      if ( jfirst .eq. 1 ) then!$omp  parallel do       &!$omp  default(shared)   &!$omp  private(i, k)         do k = 1, km            te_sp(k) = 0.            do i=ifirst,ilast#if defined (OLDWAY)              te_sp(k) = te_sp(k) + u2_sp(i,k) + v2_sp(i,k)#else              tmpik(i,k) = u2_sp(i,k) + v2_sp(i,k)              te_sp(k) = te_sp(k) + tmpik(i,k)#endif            enddo         enddo#if defined( SPMD )         if (nprxy_x .gt. 1) then# if defined (OLDWAY)            call parcollective(commxy_x, sumop, km, te_sp)# else            call par_xsum(tmpik, ifirst, ilast, im, km, te_sp)# endif         endif#endif!$omp  parallel do       &!$omp  default(shared)   &!$omp  private(i, k)         do k = 1, km            te_sp(k) = te(ifirst,1,k) - 0.5*te_sp(k)/float(im)            do i=ifirst,ilast              te(i,  1,k) = te_sp(k)            enddo         enddo      endif      if ( jlast .eq. jm ) then!$omp  parallel do       &!$omp  default(shared)   &!$omp  private(i, k)         do k = 1, km            te_np(k) = 0.            do i=ifirst,ilast#if defined (OLDWAY)              te_np(k) = te_np(k) + u2_np(i,k) + v2_np(i,k)#else              tmpik(i,k) = u2_np(i,k) + v2_np(i,k)              te_np(k) = te_np(k) + tmpik(i,k)#endif            enddo         enddo#if defined( SPMD )         if (nprxy_x .gt. 1) then# if defined (OLDWAY)            call parcollective(commxy_x, sumop, km, te_np)# else            call par_xsum(tmpik, ifirst, ilast, im, km, te_np)# endif         endif#endif!$omp  parallel do       &!$omp  default(shared)   &!$omp  private(i, k)         do k = 1, km            te_np(k) = te(ifirst,jm,k) - 0.5*te_np(k)/float(im)            do i=ifirst,ilast              te(i,jm,k) = te_np(k)            enddo         enddo      endif! Recover (virtual) temperature!$omp  parallel do        &!$omp  default(shared)    &!$omp  private(ixj, i1, i2, i, j, k, rg, gz, dlnp)!     do 9000 j=jfirst,jlast      do 9000 ixj=1,jp         j  = jfirst + (ixj-1) / nxu         i1 = ifirst + it * mod(ixj-1, nxu)         i2 = i1 + it - 1         rg = akap * cp         do i=i1,i2            gz(i) = hs(i,j)               enddo        do k=km,1,-1          do i=i1,i2            dlnp  = rg*(peln(i,k+1,j) - peln(i,k,j))            tvm(i,k,j)  = delp(i,j,k)*(te(i,j,k) - gz(i)) /     &                        ( cp*delp(i,j,k) - pe(i,k,j)*dlnp )! Update phis            gz(i) = gz(i) + dlnp*tvm(i,k,j)          enddo          if( consv ) then              do i=i1,i2                 tvm(i,k,j) = tvm(i,k,j) + dtmp              enddo          endif          if( .not. convt ) then              do i=i1,i2                 pt(i,j,k) = tvm(i,k,j) / pkz(i,j,k)              enddo          endif        enddo           ! end k-loop9000  continue      return!EOC      end!-----------------------------------------------------------------------