!> \file GFS_suite_interstitial_2.f90 !! Contains code related used to calculate radiation-based and PBL-based diagnostics that are executed after radiation time interpolation and before the surface layer. module GFS_suite_interstitial_2 use machine, only: kind_phys real(kind=kind_phys), parameter :: one = 1.0_kind_phys logical :: linit_mod = .false. contains !> \section arg_table_GFS_suite_interstitial_2_run Argument Table !! \htmlinclude GFS_suite_interstitial_2_run.html !! subroutine GFS_suite_interstitial_2_run (im, levs, lssav, ldiag3d, lsidea, flag_cice, shal_cnv, old_monin, mstrat, & do_shoc, frac_grid, imfshalcnv, dtf, xcosz, adjsfcdsw, adjsfcdlw, cice, pgr, ulwsfc_cice, lwhd, htrsw, htrlw, xmu, ctei_rm, & work1, work2, prsi, tgrs, prsl, qgrs_water_vapor, qgrs_cloud_water, cp, hvap, prslk, suntim, adjsfculw, adjsfculw_lnd, & adjsfculw_ice, adjsfculw_wat, dlwsfc, ulwsfc, psmean, dtend, dtidx, index_of_process_longwave, index_of_process_shortwave, & index_of_process_pbl, index_of_process_dcnv, index_of_process_scnv, index_of_process_mp, index_of_temperature, & ctei_rml, ctei_r, kinver, dry, icy, wet, frland, huge, use_LW_jacobian, htrlwu, errmsg, errflg) implicit none ! interface variables integer, intent(in ) :: im, levs, imfshalcnv logical, intent(in ) :: lssav, ldiag3d, lsidea, shal_cnv logical, intent(in ) :: old_monin, mstrat, do_shoc, frac_grid, use_LW_jacobian real(kind=kind_phys), intent(in ) :: dtf, cp, hvap logical, intent(in ), dimension(:) :: flag_cice real(kind=kind_phys), intent(in ), dimension(:) :: ctei_rm real(kind=kind_phys), intent(in ), dimension(:) :: xcosz, adjsfcdsw, adjsfcdlw, pgr, xmu, work1, work2 real(kind=kind_phys), intent(in ), dimension(:) :: ulwsfc_cice real(kind=kind_phys), intent(in ), dimension(:) :: cice real(kind=kind_phys), intent(in ), dimension(:,:) :: htrsw, htrlw, htrlwu, tgrs, prsl, qgrs_water_vapor, qgrs_cloud_water, prslk real(kind=kind_phys), intent(in ), dimension(:,:) :: prsi real(kind=kind_phys), intent(in ), dimension(:,:,:) :: lwhd integer, intent(inout), dimension(:) :: kinver real(kind=kind_phys), intent(inout), dimension(:) :: suntim, dlwsfc, ulwsfc, psmean, ctei_rml, ctei_r real(kind=kind_phys), intent(in ), dimension(:) :: adjsfculw_lnd, adjsfculw_ice, adjsfculw_wat real(kind=kind_phys), intent(inout), dimension(:) :: adjsfculw ! dtend is only allocated if ldiag3d is .true. real(kind=kind_phys), optional, intent(inout), dimension(:,:,:) :: dtend integer, intent(in), dimension(:,:) :: dtidx integer, intent(in) :: index_of_process_longwave, index_of_process_shortwave, & index_of_process_pbl, index_of_process_dcnv, index_of_process_scnv, & index_of_process_mp, index_of_temperature logical, intent(in ), dimension(:) :: dry, icy, wet real(kind=kind_phys), intent(in ), dimension(:) :: frland real(kind=kind_phys), intent(in ) :: huge character(len=*), intent( out) :: errmsg integer, intent( out) :: errflg ! local variables real(kind=kind_phys), parameter :: czmin = 0.0001_kind_phys ! cos(89.994) integer :: i, k, idtend real(kind=kind_phys) :: tem1, tem2, tem, hocp logical, dimension(im) :: invrsn real(kind=kind_phys), dimension(im) :: tx1, tx2 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys, one = 1.0_kind_phys real(kind=kind_phys), parameter :: qmin = 1.0e-10_kind_phys, epsln=1.0e-10_kind_phys ! Initialize CCPP error handling variables errmsg = '' errflg = 0 hocp = hvap/cp if (lssav) then ! --- ... accumulate/save output variables ! --- ... sunshine duration time is defined as the length of time (in mdl output ! interval) that solar radiation falling on a plane perpendicular to the ! direction of the sun >= 120 w/m2 do i = 1, im if ( xcosz(i) >= czmin ) then ! zenth angle > 89.994 deg tem1 = adjsfcdsw(i) / xcosz(i) if ( tem1 >= 120.0_kind_phys ) then suntim(i) = suntim(i) + dtf endif endif enddo ! --- ... sfc lw fluxes used by atmospheric model are saved for output if (frac_grid) then do i=1,im tem = (one - frland(i)) * cice(i) ! tem = ice fraction wrt whole cell if (flag_cice(i)) then adjsfculw(i) = adjsfculw_lnd(i) * frland(i) & + ulwsfc_cice(i) * tem & + adjsfculw_wat(i) * (one - frland(i) - tem) else adjsfculw(i) = adjsfculw_lnd(i) * frland(i) & + adjsfculw_ice(i) * tem & + adjsfculw_wat(i) * (one - frland(i) - tem) endif enddo else do i=1,im if (dry(i)) then ! all land adjsfculw(i) = adjsfculw_lnd(i) elseif (icy(i)) then ! ice (and water) tem = one - cice(i) if (flag_cice(i)) then if (wet(i) .and. adjsfculw_wat(i) /= huge) then adjsfculw(i) = ulwsfc_cice(i)*cice(i) + adjsfculw_wat(i)*tem else adjsfculw(i) = ulwsfc_cice(i) endif else if (wet(i) .and. adjsfculw_wat(i) /= huge) then adjsfculw(i) = adjsfculw_ice(i)*cice(i) + adjsfculw_wat(i)*tem else adjsfculw(i) = adjsfculw_ice(i) endif endif else ! all water adjsfculw(i) = adjsfculw_wat(i) endif enddo endif do i=1,im dlwsfc(i) = dlwsfc(i) + adjsfcdlw(i)*dtf ulwsfc(i) = ulwsfc(i) + adjsfculw(i)*dtf psmean(i) = psmean(i) + pgr(i)*dtf ! mean surface pressure enddo if (ldiag3d) then if (lsidea) then idtend = dtidx(index_of_temperature,index_of_process_longwave) if(idtend>=1) then dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,1)*dtf endif idtend = dtidx(index_of_temperature,index_of_process_shortwave) if(idtend>=1) then dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,2)*dtf endif idtend = dtidx(index_of_temperature,index_of_process_pbl) if(idtend>=1) then dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,3)*dtf endif idtend = dtidx(index_of_temperature,index_of_process_dcnv) if(idtend>=1) then dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,4)*dtf endif idtend = dtidx(index_of_temperature,index_of_process_scnv) if(idtend>=1) then dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,5)*dtf endif idtend = dtidx(index_of_temperature,index_of_process_mp) if(idtend>=1) then dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,6)*dtf endif else idtend = dtidx(index_of_temperature,index_of_process_longwave) if(idtend>=1) then if (use_LW_jacobian) then dtend(:,:,idtend) = dtend(:,:,idtend) + htrlwu(:,:)*dtf else dtend(:,:,idtend) = dtend(:,:,idtend) + htrlw(:,:)*dtf endif endif idtend = dtidx(index_of_temperature,index_of_process_shortwave) if(idtend>=1) then do k=1,levs do i=1,im dtend(i,k,idtend) = dtend(i,k,idtend) + htrsw(i,k)*dtf*xmu(i) enddo enddo endif endif endif endif ! end if_lssav_block do i=1, im invrsn(i) = .false. tx1(i) = zero tx2(i) = 10.0_kind_phys ctei_r(i) = 10.0_kind_phys enddo if ((((imfshalcnv == 0 .and. shal_cnv) .or. old_monin) .and. mstrat) & .or. do_shoc) then ctei_rml(:) = ctei_rm(1)*work1(:) + ctei_rm(2)*work2(:) do k=1,levs/2 do i=1,im if (prsi(i,1)-prsi(i,k+1) < 0.35_kind_phys*prsi(i,1) & .and. (.not. invrsn(i))) then tem = (tgrs(i,k+1) - tgrs(i,k)) & / (prsl(i,k) - prsl(i,k+1)) if (((tem > 0.0001_kind_phys) .and. (tx1(i) < zero)) .or. & ((tem-abs(tx1(i)) > zero) .and. (tx2(i) < zero))) then invrsn(i) = .true. if (qgrs_water_vapor(i,k) > qgrs_water_vapor(i,k+1)) then tem1 = tgrs(i,k+1) + hocp*max(qgrs_water_vapor(i,k+1),qmin) tem2 = tgrs(i,k) + hocp*max(qgrs_water_vapor(i,k),qmin) tem1 = tem1 / prslk(i,k+1) - tem2 / prslk(i,k) ! --- ... (cp/l)(deltathetae)/(deltatwater) > ctei_rm -> conditon for CTEI ctei_r(i) = (one/hocp)*tem1/(qgrs_water_vapor(i,k+1)-qgrs_water_vapor(i,k) & + qgrs_cloud_water(i,k+1)-qgrs_cloud_water(i,k)) else ctei_r(i) = 10.0_kind_phys endif if ( ctei_rml(i) > ctei_r(i) ) then kinver(i) = k else kinver(i) = levs endif endif tx2(i) = tx1(i) tx1(i) = tem endif enddo enddo endif end subroutine GFS_suite_interstitial_2_run end module GFS_suite_interstitial_2