Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:40:35 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
a456aa407c Andr*0001 #include "FIZHI_OPTIONS.h"
e4ce4355da Jean*0002        SUBROUTINE FIZHI_WRAPPER (myTime, myIter, myThid)
e337e4ca8c Andr*0003 c-----------------------------------------------------------------------
                0004 c  Subroutine fizhi_wrapper - 'Wrapper' routine to interface
6b923cbf47 Jean*0005 c        with physics driver.
                0006 c        1) Set up "bi, bj loop"  and some timers and clocks.
e337e4ca8c Andr*0007 c        2) Call do_fizhi - driver for physics which computes tendencies
                0008 c        3) Interpolate tendencies to dynamics grid in vertical
                0009 c        4) Convert u,v tendencies to C-Grid
                0010 c
                0011 c Calls: do_fizhi (get u,v,t,s tend, step tke, etc and tc, etc. forward)
                0012 c        phys2dyn (4 calls - all physics tendencies)
                0013 c        AtoC (u and v tendencies)
                0014 c-----------------------------------------------------------------------
e4ce4355da Jean*0015        IMPLICIT NONE
e337e4ca8c Andr*0016 #include "SIZE.h"
                0017 #include "GRID.h"
eb35f89333 Andr*0018 #include "EEPARAMS.h"
5cbcb662ee Andr*0019 #include "PARAMS.h"
eb35f89333 Andr*0020 #include "SURFACE.h"
                0021 #include "DYNVARS.h"
f4a0368053 Andr*0022 #include "fizhi_land_SIZE.h"
e337e4ca8c Andr*0023 #include "fizhi_SIZE.h"
                0024 #include "fizhi_coms.h"
                0025 #include "gridalt_mapping.h"
d9ce6b5984 Andr*0026 #include "fizhi_land_coms.h"
                0027 #include "fizhi_earth_coms.h"
                0028 #include "fizhi_ocean_coms.h"
                0029 #include "fizhi_chemistry_coms.h"
afe658afe3 Andr*0030 #ifdef ALLOW_DIAGNOSTICS
                0031 #include "fizhi_SHP.h"
                0032 #endif
e337e4ca8c Andr*0033 
e4ce4355da Jean*0034        INTEGER myIter, myThid
3768927683 Andr*0035        _RL myTime
e4ce4355da Jean*0036        LOGICAL  diagnostics_is_on
                0037        EXTERNAL diagnostics_is_on
e337e4ca8c Andr*0038 
                0039 c pe on dynamics and physics grid refers to bottom edge
e4ce4355da Jean*0040        _RL pephy4fiz(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys+1,nSx,nSy)
                0041        _RL pephy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys+1,nSx,nSy)
                0042        _RL pedyn(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1,nSx,nSy)
                0043        _RL tempphy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys,nSx,nSy)
                0044        _RL fracland(sNx,sNy,nSx,nSy)
6b923cbf47 Jean*0045        _RL tempLdiag(sNx,sNy,Nrphys+1)
                0046        _RL tempLdiag2(sNx,sNy,Nrphys)
                0047        _RL tempdiag(sNx,sNy)
                0048        _RL slp(sNx,sNy)
e337e4ca8c Andr*0049 
e4ce4355da Jean*0050        INTEGER i, j, L, Lbotij, bi, bj
                0051        INTEGER im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
8e402e4ff2 Andr*0052        _RL grav, getcon
e337e4ca8c Andr*0053 
8e402e4ff2 Andr*0054        grav = getcon('GRAVITY')
eb35f89333 Andr*0055        idim1 = 1-OLx
                0056        idim2 = sNx+OLx
                0057        jdim1 = 1-OLy
                0058        jdim2 = sNy+OLy
                0059        im1 = 1
                0060        im2 = sNx
                0061        jm1 = 1
                0062        jm2 = sNy
e337e4ca8c Andr*0063 
afe658afe3 Andr*0064 #ifdef ALLOW_DIAGNOSTICS
                0065       if ( useDiagnostics ) then
                0066        if(diagnostics_is_on('TENDUFIZ',myThid) .or.
e4ce4355da Jean*0067      &       diagnostics_is_on('CORRDU  ',myThid) ) then
afe658afe3 Andr*0068         do bj = myByLo(myThid), myByHi(myThid)
                0069         do bi = myBxLo(myThid), myBxHi(myThid)
                0070         do L = 1,Nrphys
                0071         do j = 1,sNy
                0072         do i = 1,sNx
                0073          ubef(i,j,L,bi,bj) = uphy(i,j,L,bi,bj)
                0074         enddo
                0075         enddo
                0076         enddo
                0077         do L = 1,Nr
                0078         do j = 1,sNy
                0079         do i = 1,sNx+1
                0080          udynbef(i,j,L,bi,bj) = uvel(i,j,L,bi,bj)
                0081         enddo
                0082         enddo
                0083         enddo
                0084         enddo
                0085         enddo
                0086        endif
                0087        if(diagnostics_is_on('TENDVFIZ',myThid) .or.
e4ce4355da Jean*0088      &       diagnostics_is_on('CORRDV  ',myThid) ) then
afe658afe3 Andr*0089         do bj = myByLo(myThid), myByHi(myThid)
                0090         do bi = myBxLo(myThid), myBxHi(myThid)
                0091         do L = 1,Nrphys
                0092         do j = 1,sNy
                0093         do i = 1,sNx
                0094          vbef(i,j,L,bi,bj) = vphy(i,j,L,bi,bj)
                0095         enddo
                0096         enddo
                0097         enddo
                0098         do L = 1,Nr
                0099         do j = 1,sNy+1
                0100         do i = 1,sNx
                0101          vdynbef(i,j,L,bi,bj) = vvel(i,j,L,bi,bj)
                0102         enddo
                0103         enddo
                0104         enddo
                0105         enddo
                0106         enddo
                0107        endif
                0108        if(diagnostics_is_on('TENDTFIZ',myThid) .or.
e4ce4355da Jean*0109      &       diagnostics_is_on('CORRDT  ',myThid) ) then
afe658afe3 Andr*0110         do bj = myByLo(myThid), myByHi(myThid)
                0111         do bi = myBxLo(myThid), myBxHi(myThid)
                0112         do L = 1,Nrphys
                0113         do j = 1,sNy
                0114         do i = 1,sNx
                0115          thbef(i,j,L,bi,bj) = thphy(i,j,L,bi,bj)
                0116         enddo
                0117         enddo
                0118         enddo
                0119         do L = 1,Nr
                0120         do j = 1,sNy
                0121         do i = 1,sNx
                0122          thdynbef(i,j,L,bi,bj) = theta(i,j,L,bi,bj)
                0123         enddo
                0124         enddo
                0125         enddo
                0126         enddo
                0127         enddo
                0128        endif
                0129        if(diagnostics_is_on('TENDQFIZ',myThid) .or.
e4ce4355da Jean*0130      &       diagnostics_is_on('CORRDQ  ',myThid) ) then
afe658afe3 Andr*0131         do bj = myByLo(myThid), myByHi(myThid)
                0132         do bi = myBxLo(myThid), myBxHi(myThid)
                0133         do L = 1,Nrphys
                0134         do j = 1,sNy
                0135         do i = 1,sNx
                0136          sbef(i,j,L,bi,bj) = sphy(i,j,L,bi,bj)
                0137         enddo
                0138         enddo
                0139         enddo
                0140         do L = 1,Nr
                0141         do j = 1,sNy
                0142         do i = 1,sNx
                0143          sdynbef(i,j,L,bi,bj) = salt(i,j,L,bi,bj)
                0144         enddo
                0145         enddo
                0146         enddo
                0147         enddo
                0148         enddo
                0149        endif
                0150       endif
                0151 #endif
                0152 
e337e4ca8c Andr*0153        do bj = myByLo(myThid), myByHi(myThid)
                0154        do bi = myBxLo(myThid), myBxHi(myThid)
                0155 
d0b11e35fb Andr*0156 c Construct the physics grid pressures
                0157 C  Note: Need one array to send to fizhi (top-down) and another
                0158 C        For the interpolations between physics and dynamics (bottom-up)
e337e4ca8c Andr*0159         do j = 1,sNy
                0160         do i = 1,sNx
00f44e1146 Andr*0161          pephy(i,j,1,bi,bj)=(Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj))
e337e4ca8c Andr*0162          do L = 2,Nrphys+1
                0163           pephy(i,j,L,bi,bj)=pephy(i,j,L-1,bi,bj)-dpphys(i,j,L-1,bi,bj)
                0164          enddo
                0165 c Do not use a zero field as the top edge pressure for interpolation
d0b11e35fb Andr*0166          do L = 1,Nrphys+1
                0167           pephy4fiz(i,j,Nrphys+2-L,bi,bj)=pephy(i,j,L,bi,bj)
                0168          enddo
e337e4ca8c Andr*0169          if(pephy(i,j,Nrphys+1,bi,bj).lt.1.e-5)
e4ce4355da Jean*0170      &                               pephy(i,j,Nrphys+1,bi,bj) = 1.e-5
e337e4ca8c Andr*0171         enddo
                0172         enddo
                0173 C Build pressures on dynamics grid
                0174         do j = 1,sNy
                0175         do i = 1,sNx
                0176          do L = 1,Nr
                0177           pedyn(i,j,L,bi,bj) = 0.
                0178          enddo
                0179         enddo
                0180         enddo
                0181         do j = 1,sNy
f98d2ec0f4 Andr*0182         do i = 1,sNx
e4ce4355da Jean*0183          Lbotij = kSurfC(i,j,bi,bj)
6b923cbf47 Jean*0184          if(Lbotij.ne.0.)
e4ce4355da Jean*0185      & pedyn(i,j,Lbotij,bi,bj) = (Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj))
e337e4ca8c Andr*0186         enddo
                0187         enddo
                0188         do j = 1,sNy
                0189         do i = 1,sNx
e4ce4355da Jean*0190          Lbotij = kSurfC(i,j,bi,bj)
e337e4ca8c Andr*0191          do L = Lbotij+1,Nr+1
6b923cbf47 Jean*0192           pedyn(i,j,L,bi,bj) = pedyn(i,j,L-1,bi,bj) -
e4ce4355da Jean*0193      &           drF(L-1)* rStarExpC(i,j,bi,bj)*hfacC(i,j,L-1,bi,bj)
e337e4ca8c Andr*0194          enddo
                0195 c Do not use a zero field as the top edge pressure for interpolation
                0196          if(pedyn(i,j,Nr+1,bi,bj).lt.1.e-5)
e4ce4355da Jean*0197      &                               pedyn(i,j,Nr+1,bi,bj) = 1.e-5
e337e4ca8c Andr*0198         enddo
                0199         enddo
180c245802 Andr*0200        enddo
                0201        enddo
                0202 
5cbcb662ee Andr*0203 #ifdef ALLOW_DIAGNOSTICS
                0204       if ( useDiagnostics ) then
                0205        if(diagnostics_is_on('FIZPRES ',myThid) ) then
                0206         do bj = myByLo(myThid), myByHi(myThid)
                0207         do bi = myBxLo(myThid), myBxHi(myThid)
                0208         do j = 1,sNy
                0209         do i = 1,sNx
                0210         do L = 1,Nrphys
                0211          tempphy(i,j,L,bi,bj) = pephy4fiz(i,j,L,bi,bj)
                0212         enddo
                0213         enddo
                0214         enddo
                0215         enddo
                0216         enddo
                0217         call diagnostics_fill(tempphy,'FIZPRES ',0,
e4ce4355da Jean*0218      &                                     Nrphys,0,1,1,myThid)
5cbcb662ee Andr*0219        endif
                0220       endif
                0221 #endif
                0222 
e4ce4355da Jean*0223        CALL TIMER_START ('DO_FIZHI          [FIZHI_WRAPPER]',myThid)
180c245802 Andr*0224        do bj = myByLo(myThid), myByHi(myThid)
                0225        do bi = myBxLo(myThid), myBxHi(myThid)
e4ce4355da Jean*0226         call get_landfrac(im2,jm2,nSx,nSy,bi,bj,maxtyp,
                0227      &        surftype,tilefrac,fracland(1,1,bi,bj))
76a78cf632 Andr*0228 
                0229 #ifdef ALLOW_DIAGNOSTICS
                0230       if ( useDiagnostics ) then
                0231        if(diagnostics_is_on('SLP     ',myThid) ) then
6b923cbf47 Jean*0232         L = Nrphys+1
76a78cf632 Andr*0233         do j = 1,sNy
                0234         do i = 1,sNx
6b923cbf47 Jean*0235          tempdiag(i,j) = topoZ(i,j,bi,bj)*grav
                0236          tempLdiag(i,j,L) = pephy4fiz(i,j,L,bi,bj)/100.
76a78cf632 Andr*0237         enddo
                0238         enddo
6b923cbf47 Jean*0239         do L = 1,Nrphys
76a78cf632 Andr*0240         do j = 1,sNy
                0241         do i = 1,sNx
6b923cbf47 Jean*0242          tempLdiag(i,j,L) = pephy4fiz(i,j,L,bi,bj)/100.
                0243          tempLdiag2(i,j,L) = thphy(i,j,L,bi,bj) *
e4ce4355da Jean*0244      &        (1.+0.609*sphy(i,j,L,bi,bj))
6b923cbf47 Jean*0245 
                0246         enddo
76a78cf632 Andr*0247         enddo
                0248         enddo
6b923cbf47 Jean*0249         call slprs( tempdiag, tempLdiag, tempLdiag2,
e4ce4355da Jean*0250      &              fracland(1,1,bi,bj),sNx,sNy,Nrphys,slp)
6b923cbf47 Jean*0251         call diagnostics_fill( slp,'SLP     ',1,1,3,bi,bj,myThid )
76a78cf632 Andr*0252        endif
                0253       endif
                0254 #endif
e337e4ca8c Andr*0255 c
                0256 c Compute physics increments
e4ce4355da Jean*0257         turbStart(bi,bj) = turbStart(bi,bj) .AND. ( myIter.EQ.nIter0 )
9524fe64b5 Andr*0258 
c3f10ffaf9 Andr*0259         call do_fizhi(myIter,myThid,
e4ce4355da Jean*0260      &  idim1,idim2,jdim1,jdim2,Nrphys,nSx,nSy,im1,im2,jm1,jm2,bi,bj,
                0261      &  turbStart(bi,bj), nchp,nchptot,nchpland,
                0262      &  uphy,vphy,thphy,sphy,pephy4fiz,xC,yC,topoZ,
                0263      &  ctmt,xxmt,yymt,zetamt,xlmt,khmt,tke,
                0264      &  tgz,sst,sice,phis_var,landtype,fracland,emiss,albnirdr,albnirdf,
                0265      &  albvisdr,albvisdf,ityp,chfr,alai,agrn,igrd,chlt,chlon,
                0266      &  tcanopy,tdeep,ecanopy,swetshal,swetroot,swetdeep,snodep,capac,
                0267      &  o3,qstr,co2,cfc11,cfc12,cfc22,n2o,methane,
                0268      &  iras,nlwcld,cldtot_lw,cldras_lw,cldlsp_lw,nlwlz,lwlz,
                0269      &  nswcld,cldtot_sw,cldras_sw,cldlsp_sw,nswlz,swlz,
                0270      &  imstturbsw,imstturblw,qliqavesw,qliqavelw,fccavesw,fccavelw,
                0271      &  raincon,rainlsp,snowfall,
                0272      &  duphy,dvphy,dthphy,dsphy)
180c245802 Andr*0273        enddo
                0274        enddo
046c3794c1 Andr*0275 
e4ce4355da Jean*0276        CALL TIMER_STOP ('DO_FIZHI          [FIZHI_WRAPPER]',myThid)
180c245802 Andr*0277 
e4ce4355da Jean*0278        CALL TIMER_START ('PHYS2DYN          [FIZHI_WRAPPER]',myThid)
180c245802 Andr*0279        do bj = myByLo(myThid), myByHi(myThid)
                0280        do bi = myBxLo(myThid), myBxHi(myThid)
e337e4ca8c Andr*0281 c Interpolate (A-Grid) physics increments to dynamics grid
6b923cbf47 Jean*0282 C   First flip the physics arrays (which are top-down)
d0b11e35fb Andr*0283 C   into bottom-up arrays for interpolation to dynamics grid
                0284         do j = 1,sNy
                0285         do i = 1,sNx
                0286          do L = 1,Nrphys
28e90e66e5 Andr*0287           tempphy(i,j,Nrphys+1-L,bi,bj)=duphy(i,j,L,bi,bj)
d0b11e35fb Andr*0288          enddo
                0289         enddo
                0290         enddo
28e90e66e5 Andr*0291         call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
e4ce4355da Jean*0292      & nSx,nSy,im1,im2,jm1,jm2,bi,bj,pedyn,kSurfC,Nr,nlperdyn,guphy)
d0b11e35fb Andr*0293         do j = 1,sNy
                0294         do i = 1,sNx
                0295          do L = 1,Nrphys
28e90e66e5 Andr*0296           tempphy(i,j,Nrphys+1-L,bi,bj)=dvphy(i,j,L,bi,bj)
d0b11e35fb Andr*0297          enddo
                0298         enddo
                0299         enddo
28e90e66e5 Andr*0300         call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
e4ce4355da Jean*0301      & nSx,nSy,im1,im2,jm1,jm2,bi,bj,pedyn,kSurfC,Nr,nlperdyn,gvphy)
d0b11e35fb Andr*0302         do j = 1,sNy
                0303         do i = 1,sNx
                0304          do L = 1,Nrphys
28e90e66e5 Andr*0305           tempphy(i,j,Nrphys+1-L,bi,bj)=dthphy(i,j,L,bi,bj)
d0b11e35fb Andr*0306          enddo
                0307         enddo
                0308         enddo
28e90e66e5 Andr*0309         call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
e4ce4355da Jean*0310      & nSx,nSy,im1,im2,jm1,jm2,bi,bj,pedyn,kSurfC,Nr,nlperdyn,gthphy)
d0b11e35fb Andr*0311         do j = 1,sNy
                0312         do i = 1,sNx
                0313          do L = 1,Nrphys
28e90e66e5 Andr*0314           tempphy(i,j,Nrphys+1-L,bi,bj)=dsphy(i,j,L,bi,bj)
d0b11e35fb Andr*0315          enddo
                0316         enddo
                0317         enddo
28e90e66e5 Andr*0318         call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
e4ce4355da Jean*0319      & nSx,nSy,im1,im2,jm1,jm2,bi,bj,pedyn,kSurfC,Nr,nlperdyn,gsphy)
e337e4ca8c Andr*0320 
                0321        enddo
                0322        enddo
e184c4ee7b Andr*0323 
e4ce4355da Jean*0324        CALL TIMER_STOP ('PHYS2DYN          [FIZHI_WRAPPER]',myThid)
e337e4ca8c Andr*0325 
                0326 c Convert guphy and gvphy from A-grid to C-grid for use by dynamics
e4ce4355da Jean*0327        CALL TIMER_START ('ATOC              [FIZHI_WRAPPER]',myThid)
eb35f89333 Andr*0328        call AtoC(myThid,guphy,gvphy,maskC,idim1,idim2,jdim1,jdim2,Nr,
e4ce4355da Jean*0329      &                      nSx,nSy,im1,im2,jm1,jm2,guphy,gvphy)
                0330        CALL TIMER_STOP ('ATOC              [FIZHI_WRAPPER]',myThid)
e337e4ca8c Andr*0331 
e4ce4355da Jean*0332        CALL TIMER_START ('EXCHANGES         [FIZHI_WRAPPER]',myThid)
f98d2ec0f4 Andr*0333 c Call the c-grid exchange routine to fill in the halo regions (du,dv)
                0334        call exch_uv_xyz_RL(guphy,gvphy,.TRUE.,myThid)
                0335 c Call the a-grid exchange routine to fill in the halo regions (dth,ds)
6637358eea Jean*0336        _EXCH_XYZ_RL(gthphy,myThid)
                0337        _EXCH_XYZ_RL(gsphy,myThid)
e4ce4355da Jean*0338        CALL TIMER_STOP ('EXCHANGES         [FIZHI_WRAPPER]',myThid)
f98d2ec0f4 Andr*0339 
e337e4ca8c Andr*0340       return
                0341       end