pkg/cost/cost_atlantic_heat.F

5ed655852f Jean*0001 #include "COST_OPTIONS.h"
720be40b89 Patr*0002
                0003       subroutine cost_atlantic_heat( myThid )
                0004 C     /==========================================================\
                0005 C     | subroutine cost_atlantic_heat                            |
                0006 C     | o This routine computes the meridional heat transport.   |
                0007 C     |   The current indices are for North Atlantic 29N         |
                0008 C     |   2x2 global setup.                                      |
                0009 C     \==========================================================/
                0010        implicit none
                0011
                0012 C     == Global variables ===
                0013 #include "SIZE.h"
                0014 #include "EEPARAMS.h"
                0015 #include "PARAMS.h"
                0016 #include "GRID.h"
                0017 #include "DYNVARS.h"
                0018 #include "cost.h"
                0019
                0020 C     ======== Routine arguments ======================
                0021 C     myThid - Thread number for this instance of the routine.
                0022       integer myThid
                0023
                0024 #ifdef ALLOW_COST_ATLANTIC_HEAT
                0025 C     ========= Local variables =========================
                0026       integer    isecbeg      , isecend      , jsec
449f3c105b Patr*0027       integer    jsecbeg      , jsecend      , isec
720be40b89 Patr*0028       integer    kmaxdepth
                0029       integer i, j, k
                0030       integer ig, jg
                0031       integer bi, bj
                0032       _RL     locfc
781d8676b2 Patr*0033       _RL     uVel_bar(Nr), vVel_bar(Nr), theta_bar(Nr)
                0034       _RL     thetaUvel_bar(Nr), thetaVvel_bar(Nr)
449f3c105b Patr*0035       _RL     countU(Nr), countV(Nr), countT(Nr)
781d8676b2 Patr*0036       _RL     countTU(Nr), countTV(Nr)
720be40b89 Patr*0037       _RL     petawatt
                0038       _RL     sum
c1d89bddf8 Jean*0039       parameter( petawatt = 1. _d +15 )
720be40b89 Patr*0040
                0041 C     80W - 0W at 24N
449f3c105b Patr*0042       parameter( isecbeg = 69, isecend = 87, jsec = 28 )
                0043 cph      parameter( isecbeg = 70, isecend = 90, jsec = 30 )
                0044       parameter( jsecbeg = 10, jsecend = 27, isec = 59 )
a4efec1784 Patr*0045 #ifdef ALLOW_COST_ATLANTIC_HEAT_DOMASS
                0046       parameter ( kmaxdepth = 8 )
                0047 #else
                0048       parameter ( kmaxdepth = 14 )
                0049 #endif
720be40b89 Patr*0050
                0051 C------------------------------------------------------
                0052 C     Accumulate meridionally integrated transports
                0053 C     Note bar(V)*bar(T) not bar(VT)
e305438401 Mart*0054 C     Attention pYFaceA [m^2*gravity*rhoConst]
720be40b89 Patr*0055 C------------------------------------------------------
                0056
                0057       DO bj=myByLo(myThid),myByHi(myThid)
                0058        DO bi=myBxLo(myThid),myBxHi(myThid)
                0059
                0060         locfc = 0.0
8347a61619 Patr*0061         sum = 0.0
                0062
449f3c105b Patr*0063 #define MERID_TRANSPORT
720be40b89 Patr*0064
449f3c105b Patr*0065 #ifdef MERID_TRANSPORT
720be40b89 Patr*0066
449f3c105b Patr*0067 #undef ENERGYNORM
720be40b89 Patr*0068 #ifdef ENERGYNORM
                0069
                0070         do i=1,sNx
                0071          ig = myXGlobalLo-1+(bi-1)*sNx+i
                0072          if ((ig .ge. isecbeg) .and. (ig .le. isecend)) then
                0073           sum = 0.0
                0074           do k = 1, kmaxdepth
                0075            sum = sum
                0076      &            + (vVel(i,j,k,bi,bj) * maskS(i,j,k,bi,bj)
                0077      &            * drF(k))**2
                0078           end do
                0079           locfc = locfc + sum*dxG(i,j,bi,bj)
                0080          end if
                0081         end do
                0082
                0083 #else
                0084
449f3c105b Patr*0085         do j=1,sNy
                0086          jg = myYGlobalLo-1+(bj-1)*sNy+j
                0087          if (jg .eq. jsec) then
                0088
720be40b89 Patr*0089           do k = 1, Nr
                0090            vVel_bar(k) = 0.0
781d8676b2 Patr*0091            thetaVvel_bar(k) = 0.0
8347a61619 Patr*0092            countV(k) = 0.0
781d8676b2 Patr*0093            countTV(k) = 0.0
720be40b89 Patr*0094            do i=1,sNx
                0095             ig = myXGlobalLo-1+(bi-1)*sNx+i
8347a61619 Patr*0096 c
720be40b89 Patr*0097             if ((ig .ge. isecbeg) .and. (ig .le. isecend)) then
951926fb9b Jean*0098                vVel_bar(k) = vVel_bar(k)
a4efec1784 Patr*0099      &              + cMeanVVel(i,j,k,bi,bj)*dxG(i,j,bi,bj)
                0100      &                 *maskS(i,j,k,bi,bj)
781d8676b2 Patr*0101
                0102                thetaVvel_bar(k) = thetaVvel_bar(k)
                0103      &            + cMeanThetaVVel(i,j,k,bi,bj)*dxG(i,j,bi,bj)
                0104      &                 *maskS(i,j,k,bi,bj)*maskC(i,j,k,bi,bj)
                0105
951926fb9b Jean*0106              countTV(k) = countTV(k) +
781d8676b2 Patr*0107      &            maskS(i,j,k,bi,bj)*maskC(i,j,k,bi,bj)
951926fb9b Jean*0108              countV(k) = countV(k) +
781d8676b2 Patr*0109      &            maskS(i,j,k,bi,bj)
720be40b89 Patr*0110             end if
8347a61619 Patr*0111
720be40b89 Patr*0112            enddo
                0113           enddo
8347a61619 Patr*0114 c
                0115           do k = 1, Nr
781d8676b2 Patr*0116 #ifdef ALLOW_COST_ATLANTIC_HEAT_DOMASS
                0117            if ( k .LE. kmaxdepth .AND. countV(k) .NE. 0) then
8347a61619 Patr*0118             sum = sum
781d8676b2 Patr*0119      &            + vVel_bar(k)*drF(k)/countV(k)
8347a61619 Patr*0120            end if
781d8676b2 Patr*0121 #else
                0122            if ( k .LE. kmaxdepth .AND. countTV(k) .NE. 0) then
                0123             sum = sum
                0124      &            + thetaVVel_bar(k)*drF(k)/countTV(k)
                0125            end if
                0126 #endif
720be40b89 Patr*0127           end do
                0128
449f3c105b Patr*0129 #endif /* ENERGYNORM */
                0130
                0131 #else
                0132
951926fb9b Jean*0133 cph need to change this part to go from
781d8676b2 Patr*0134 cph \bar{u}*\bar{T} to \bar{u*T}
                0135 cph (required store dir. are now in place)
                0136
449f3c105b Patr*0137         do i=1,sNx
                0138          ig = myXGlobalLo-1+(bi-1)*sNx+i
                0139          if (ig .eq. isec) then
                0140
                0141           do k = 1, Nr
                0142            uVel_bar(k) = 0.0
                0143            theta_bar(k) = 0.0
                0144            countT(k) = 0.0
                0145            countU(k) = 0.0
                0146            do j=1,sNy
                0147             jg = myYGlobalLo-1+(bj-1)*sNy+j
                0148 c
                0149             if ((jg .ge. jsecbeg) .and. (jg .le. jsecend)) then
                0150              uVel_bar(k) = uVel_bar(k)
                0151      &                      + cMeanUVel(i,j,k,bi,bj)
                0152      &                      *maskW(i,j,k,bi,bj)
                0153      &                      *maskC(i,j,k,bi,bj)*maskC(i-1,j,k,bi,bj)
                0154              theta_bar(k) = theta_bar(k) +
                0155      &            0.5*( cMeanTheta(i,j,k,bi,bj)
                0156      &                 +cMeanTheta(i-,j,k,bi,bj) )
                0157      &                 *maskW(i,j,k,bi,bj)*dyG(i,j,bi,bj)
                0158      &                 *maskC(i,j,k,bi,bj)*maskC(i-1,j,k,bi,bj)
                0159              countT(k) = countT(k) + maskW(i,j,k,bi,bj)
                0160      &                   *maskC(i,j,k,bi,bj)*maskC(i-1,j,k,bi,bj)
                0161              countU(k) = countU(k) + maskW(i,j,k,bi,bj)
                0162      &                   *maskC(i,j,k,bi,bj)*maskC(i-1,j,k,bi,bj)
                0163             end if
                0164
                0165            enddo
                0166           enddo
                0167 c
                0168           do k = 1, Nr
                0169            if ( k .LE. kmaxdepth .AND.
                0170      &          countT(k) .NE. 0 .AND. countU(k) .NE. 0) then
                0171             sum = sum
                0172      &            + uVel_bar(k) * theta_bar(k) * drF(k)
                0173      &            / ( countT(k) * countU(k) )
                0174            end if
                0175           end do
                0176
720be40b89 Patr*0177 #endif
                0178
                0179          end if
                0180         end do
                0181
781d8676b2 Patr*0182 #ifdef ALLOW_COST_ATLANTIC_HEAT_DOMASS
951926fb9b Jean*0183         objf_atl(bi,bj) =
781d8676b2 Patr*0184      &     sum*1.E-6
                0185 #else
951926fb9b Jean*0186         objf_atl(bi,bj) =
449f3c105b Patr*0187      &     sum*HeatCapacity_Cp*rhoConst/petawatt
781d8676b2 Patr*0188 #endif
720be40b89 Patr*0189
                0190 c--   end of bi,bj loop
                0191        end do
                0192       end do
                0193
                0194 #endif
951926fb9b Jean*0195
720be40b89 Patr*0196       end