C $Header: /u/gcmpack/MITgcm/pkg/aim_v23/phy_vdifsc.F,v 1.1 2002/11/22 17:17:03 jmc Exp $ C $Name: $ #include "AIM_OPTIONS.h" SUBROUTINE VDIFSC (dpFac,SE,RH,QA,QSAT, O TTENVD,QTENVD, I kGrd,bi,bj,myThid) C-- C-- SUBROUTINE VDIFSC (UA,VA,SE,RH,QA,QSAT,PHI, C-- & UTENVD,VTENVD,TTENVD,QTENVD) C-- C-- Purpose: Compute tendencies of momentum, energy and moisture C-- due to vertical diffusion and shallow convection C-- Input: UA = u-wind (3-dim) C-- VA = v-wind (3-dim) C dpFac = cell delta_P fraction (3-dim) C-- SE = dry static energy (3-dim) C-- RH = relative humidity [0-1] (3-dim) C-- QA = specific humidity [g/kg] (3-dim) C-- QSAT = saturation sp. humidity [g/kg] (3-dim) C-- PHI = geopotential (3-dim) C-- Output: UTENVD = u-wind tendency (3-dim) C-- VTENVD = v-wind tendency (3-dim) C-- TTENVD = temperature tendency (3-dim) C-- QTENVD = sp. humidity tendency [g/(kg s)] (3-dim) C Input: kGrd = Ground level index (2-dim) C bi,bj = tile index C myThid = Thread number for this instance of the routine C------- C Note: a) dry static energy (SE,input) has been replaced by Pot.Temp. C b) UA,VA & U,V_TENVD not used => removed C------- C In contrast to other Physics S/R, VDIFSC return a real tendency (dQ/dt,dT/dt) C nevertheless /dpFac is not applied here but later in AIM_AIM2DYN C------- IMPLICIT NONE C Resolution parameters C-- size for MITgcm & Physics package : #include "AIM_SIZE.h" #include "EEPARAMS.h" C Physical constants + functions of sigma and latitude #include "com_physcon.h" C Vertical diffusion constants #include "com_vdicon.h" C-- Routine arguments: c _RL UA(NGP,NLEV), VA(NGP,NLEV) _RL dpFac(NGP,NLEV) _RL SE(NGP,NLEV), RH(NGP,NLEV), QA(NGP,NLEV), QSAT(NGP,NLEV) c _RL PHI(NGP,NLEV) c _RL UTENVD(NGP,NLEV), VTENVD(NGP,NLEV) _RL TTENVD(NGP,NLEV), QTENVD(NGP,NLEV) INTEGER kGrd(NGP) INTEGER bi,bj,myThid #ifdef ALLOW_AIM C-- Local variables: INTEGER J, K, Ktmp, NL1 _RL RSIG(NLEV) _RL dSEdp(NGP,NLEV-1), DeltaPI(NLEV-1), factP C- jmc: declare all local variables: _RL CVDI(NGP), FSHCQ _RL DRH0, DRH, DMSE, FLUXSE, FLUXQ C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| C-- 1. Initalization C N.B. In this routine, fluxes of dry static energy and humidity C are scaled in such a way that: C d_T/dt = d_F'(SE)/d_sigma, d_Q/dt = d_F'(Q)/d_sigma c_FM NL1 = NLEV-1 c_FM CSHC = DSIG(NLEV)/3600. c_FM CVDI = (SIGH(NL1)-SIGH(1))/((NL1-1)*3600.) c_FM FSHCQ = CSHC/TRSHC c_FM FSHCSE = CSHC/(TRSHC*CP) c_FM FVDIQ = CVDI/TRVDI c_FM FVDISE = CVDI/(TRVDS*CP) DO J=1,NGP NL1 = kGrd(J)-1 CVDI(J) = 0. IF (NL1.GE.2) THEN CVDI(J) = (SIGH(NL1)-SIGH(1))/((NL1-1)*3600. _d 0) ENDIF ENDDO DO K=1,NLEV RSIG(K)=1./DSIG(K) ENDDO DO K=1,NLEV DO J=1,NGP c UTENVD(J,K) = 0. c VTENVD(J,K) = 0. TTENVD(J,K) = 0. QTENVD(J,K) = 0. ENDDO ENDDO C --------------------------------------------- C Write Conditional stability based on Pot.Temp : C dSEdp(K) = Delta[Static-Energy] between 2 Plevels(k,k+1) C and corresponds to SE(K+1)-SE(K) in the original code C ------- DO K=1,NLEV-1 factP = CP*SIGH(K)**(RD/CP) DO J=1,NGP dSEdp(J,K)=(SE(J,K+1)-SE(J,K))*factP ENDDO DeltaPI(K) = SIG(K+1)**(RD/CP) - SIG(K)**(RD/CP) ENDDO C-- 2. Shallow convection DO J=1,NGP Ktmp = kGrd(J) NL1 = Ktmp - 1 IF (Ktmp.GE.2) THEN DRH0=RHGRAD*(SIG(Ktmp)-SIG(NL1)) FSHCQ = DSIG(Ktmp)*dpFac(J,Ktmp)/(TRSHC*3600. _d 0) c_FM DMSE = (SE(J,NLEV)-SE(J,NL1))+ALHC*(QA(J,NLEV)-QSAT(J,NL1)) DMSE = dSEdp(J,NL1) + ALHC*(QA(J,Ktmp)-QSAT(J,NL1)) DRH = RH(J,Ktmp)-RH(J,NL1) IF (DMSE.GE.0.0) THEN c_FM FLUXSE = FSHCSE*DMSE FLUXSE = FSHCQ *DMSE/CP TTENVD(J,NL1) = FLUXSE*RSIG(NL1) TTENVD(J,Ktmp) =-FLUXSE*RSIG(Ktmp) IF (DRH.GE.0.0) THEN FLUXQ = FSHCQ*QSAT(J,Ktmp)*DRH QTENVD(J,NL1) = FLUXQ*RSIG(NL1) QTENVD(J,Ktmp) =-FLUXQ*RSIG(Ktmp) ENDIF ELSE IF (DRH.GE.DRH0) THEN c_FM FLUXQ = FVDIQ*QSAT(J,NL1)*DRH FLUXQ = QSAT(J,NL1)*DRH*CVDI(J)/TRVDI QTENVD(J,NL1) = FLUXQ*RSIG(NL1) QTENVD(J,Ktmp) =-FLUXQ*RSIG(Ktmp) ENDIF ENDIF ENDDO C-- 3. Vertical diffusion of moisture above the PBL DO J=1,NGP DO K=3,kGrd(J)-2 DRH0=RHGRAD*(SIG(K+1)-SIG(K)) DRH=RH(J,K+1)-RH(J,K) IF (DRH.GE.DRH0) THEN c_FM FLUXQ = FVDIQ*QSAT(J,K)*DRH FLUXQ = QSAT(J,K)*DRH*CVDI(J)/TRVDI QTENVD(J,K) = QTENVD(J,K) +FLUXQ*RSIG(K) QTENVD(J,K+1)= QTENVD(J,K+1)-FLUXQ*RSIG(K+1) ENDIF ENDDO ENDDO C-- 4. Damping of super-adiabatic lapse rate DO J=1,NGP DO K=1,kGrd(J)-1 c_FM SE0 = SE(J,K+1)+SEGRAD*(PHI(J,K)-PHI(J,K+1)) DMSE = dSEdp(J,K) & +SEGRAD*CP*DeltaPI(K)*(SE(J,K+1)+SE(J,K))*0.5 _d 0 c_FM IF (SE(J,K).LT.SE0) THEN c_FM FLUXSE = FVDISE*(SE0-SE(J,K)) IF (DMSE.GT.0.) THEN FLUXSE = DMSE*CVDI(J)/(TRVDS*CP) TTENVD(J,K ) = TTENVD(J,K )+FLUXSE*RSIG(K) TTENVD(J,K+1) = TTENVD(J,K+1)-FLUXSE*RSIG(K+1) ENDIF ENDDO ENDDO C-- C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| #endif /* ALLOW_AIM */ RETURN END