!>\file mp_fer_hires.F90 !! This file contains the Ferrier-Aligo microphysics scheme driver. ! module mp_fer_hires use machine, only : kind_phys use module_mp_fer_hires, only : ferrier_init_hr, FER_HIRES, & ferhires_finalize implicit none public :: mp_fer_hires_init, mp_fer_hires_run, mp_fer_hires_finalize private logical :: is_initialized = .False. ! * T_ICE - temperature (C) threshold at which all remaining liquid water ! is glaciated to ice ! * T_ICE_init - maximum temperature (C) at which ice nucleation occurs REAL, PUBLIC, PARAMETER :: T_ICE=-40., & T0C=273.15, & T_ICEK=T0C+T_ICE contains !> This subroutine initialize constants & lookup tables for Ferrier-Aligao MP !! scheme. !> \section arg_table_mp_fer_hires_init Argument Table !! \htmlinclude mp_fer_hires_init.html !! subroutine mp_fer_hires_init(ncol, nlev, dtp, imp_physics, & imp_physics_fer_hires, restart, & mpicomm, mpirank,mpiroot, & threads, errmsg, errflg) USE machine, ONLY : kind_phys USE MODULE_MP_FER_HIRES, ONLY : FERRIER_INIT_HR implicit none integer, intent(in) :: ncol integer, intent(in) :: nlev real(kind_phys), intent(in) :: dtp integer, intent(in) :: imp_physics integer, intent(in) :: imp_physics_fer_hires integer, intent(in) :: mpicomm integer, intent(in) :: mpirank integer, intent(in) :: mpiroot integer, intent(in) :: threads logical, intent(in) :: restart character(len=*), intent(out) :: errmsg integer, intent(out) :: errflg ! Local variables integer :: ims, ime, lm,i,k !real(kind=kind_phys) :: DT_MICRO ! Initialize the CCPP error handling variables errmsg = '' errflg = 0 if (is_initialized) return ! Set internal dimensions ims = 1 ime = ncol lm = nlev ! MZ* temporary if (mpirank==mpiroot) then write(0,*) ' -----------------------------------------------' write(0,*) ' --- !!! WARNING !!! ---' write(0,*) ' --- the CCPP Ferrier-Aligo MP scheme is ---' write(0,*) ' --- currently under development, use at ---' write(0,*) ' --- your own risk . ---' write(0,*) ' -----------------------------------------------' end if ! MZ* temporary if (imp_physics /= imp_physics_fer_hires ) then write(errmsg,'(*(a))') "Logic error: namelist choice of microphysics is different from Ferrier-Aligo MP" errflg = 1 return end if if (mpirank==mpiroot) write (0,*) 'F-A: calling FERRIER_INIT_HR ...' CALL FERRIER_INIT_HR(dtp,mpicomm,mpirank,mpiroot,threads,errmsg,errflg) if (mpirank==mpiroot) write (0,*)'F-A: FERRIER_INIT_HR finished ...' if (errflg /= 0 ) return is_initialized = .true. end subroutine mp_fer_hires_init !>\defgroup hafs_famp HWRF Ferrier-Aligo Microphysics Scheme !> This is the CCPP-compliant FER_HIRES driver module. !> \section arg_table_mp_fer_hires_run Argument Table !! \htmlinclude mp_fer_hires_run.html !! SUBROUTINE mp_fer_hires_run(NCOL, NLEV, DT ,SPEC_ADV & ,SLMSK & ,PRSI,P_PHY & ,T,Q & ,TRAIN,SR & ,QC,QR,QI,QG & ,PREC & ,mpirank, mpiroot, threads & ,refl_10cm & ,RHGRD,dx & ,EPSQ,R_D,P608,CP,G & ,errmsg,errflg) !----------------------------------------------------------------------- USE MACHINE, ONLY: kind_phys ! IMPLICIT NONE ! !----------------------------------------------------------------------- ! INTEGER,PARAMETER :: D_SS=1 ! !------------------------ !*** Argument Variables !------------------------ integer, intent(in ) :: ncol integer, intent(in ) :: nlev real(kind_phys), intent(in ) :: dt integer, intent(in ) :: threads logical, intent(in ) :: spec_adv integer, intent(in ) :: mpirank integer, intent(in ) :: mpiroot real(kind_phys), intent(in ) :: slmsk(:) real(kind_phys), intent(in ) :: prsi(:,:) real(kind_phys), intent(in ) :: p_phy(:,:) real(kind_phys), intent(in ) :: epsq,r_d,p608,cp,g real(kind_phys), intent(inout) :: t(:,:) real(kind_phys), intent(inout) :: q(:,:) real(kind_phys), intent(inout) :: train(:,:) real(kind_phys), intent(out ) :: sr(:) real(kind_phys), intent(inout) :: qc(:,:) real(kind_phys), intent(inout) :: qr(:,:) real(kind_phys), intent(inout) :: qi(:,:) real(kind_phys), intent(inout) :: qg(:,:) ! QRIMEF real(kind_phys), intent(inout) :: prec(:) real(kind_phys), intent(inout) :: refl_10cm(:,:) real(kind_phys), intent(in ) :: rhgrd real(kind_phys), intent(in ) :: dx(:) character(len=*), intent(out) :: errmsg integer, intent(out) :: errflg ! !--------------------- !*** Local Variables !--------------------- ! integer :: I,J,K,N integer :: lowlyr(1:ncol) integer :: dx1 real(kind_phys) :: PCPCOL real(kind_phys) :: ql(1:nlev),tl(1:nlev) real(kind_phys) :: rainnc(1:ncol),rainncv(1:ncol) real(kind_phys) :: snownc(1:ncol),snowncv(1:ncol) real(kind_phys) :: graupelncv(1:ncol) real(kind_phys) :: train_phy(1:ncol,1:nlev) real(kind_phys) :: f_ice(1:ncol,1:nlev) real(kind_phys) :: f_rain(1:ncol,1:nlev) real(kind_phys) :: f_rimef(1:ncol,1:nlev) real(kind_phys) :: cwm(1:ncol,1:nlev) ! Dimension integer :: ims, ime, lm !----------------------------------------------------------------------- !*********************************************************************** !----------------------------------------------------------------------- ! Initialize the CCPP error handling variables errmsg = '' errflg = 0 ! Check initialization state if (.not. is_initialized) then write(errmsg, fmt='((a))') 'mp_fer_hires_run called before mp_fer_hires_init' errflg = 1 return end if ! Set internal dimensions ims = 1 ime = ncol lm = nlev ! Use the dx of the 1st i point to set an integer value of dx to be used for ! determining where RHgrd should be set to 0.98 in the coarse domain when running HAFS. DX1=NINT(DX(1)) !----------------------------------------------------------------------- !*** NOTE: THE NMMB HAS IJK STORAGE WITH LAYER 1 AT THE TOP. !*** THE WRF PHYSICS DRIVERS HAVE IKJ STORAGE WITH LAYER 1 !*** AT THE BOTTOM. !----------------------------------------------------------------------- !....................................................................... DO I=IMS,IME ! LOWLYR(I)=1 ! !----------------------------------------------------------------------- !*** FILL RAINNC WITH ZERO (NORMALLY CONTAINS THE NONCONVECTIVE !*** ACCUMULATED RAIN BUT NOT YET USED BY NMM) !*** COULD BE OBTAINED FROM ACPREC AND CUPREC (ACPREC-CUPREC) !----------------------------------------------------------------------- !..The NC variables were designed to hold simulation total accumulations !.. whereas the NCV variables hold timestep only values, so change below !.. to zero out only the timestep amount preparing to go into each !.. micro routine while allowing NC vars to accumulate continually. !.. But, the fact is, the total accum variables are local, never saved !.. nor written so they go nowhere at the moment. ! RAINNC (I)=0. ! NOT YET USED BY NMM RAINNCv(I)=0. SNOWNCv(I)=0. graupelncv(i) = 0.0 ! !----------------------------------------------------------------------- !*** FILL THE SINGLE-COLUMN INPUT !----------------------------------------------------------------------- ! DO K=LM,1,-1 !mz* We are moving down from the top in the flipped arrays !*** CALL MICROPHYSICS !MZ* in HWRF !-- 6/11/2010: Update cwm, F_ice, F_rain and F_rimef arrays cwm(I,K)=QC(I,K)+QR(I,K)+QI(I,K) IF (QI(I,K) <= EPSQ) THEN F_ICE(I,K)=0. F_RIMEF(I,K)=1. IF (T(I,K) < T_ICEK) F_ICE(I,K)=1. ELSE F_ICE(I,K)=MAX( 0., MIN(1., QI(I,K)/cwm(I,K) ) ) F_RIMEF(I,K)=QG(I,K)!/QI(I,K) ENDIF IF (QR(I,K) <= EPSQ) THEN F_RAIN(I,K)=0. ELSE F_RAIN(I,K)=QR(I,K)/(QR(I,K)+QC(I,K)) ENDIF ENDDO ENDDO !--------------------------------------------------------------------- !aligo DO K = 1, LM DO I= IMS, IME cwm(i,k) = cwm(i,k)/(1.0_kind_phys-q(i,k)) qr(i,k) = qr(i,k)/(1.0_kind_phys-q(i,k)) qi(i,k) = qi(i,k)/(1.0_kind_phys-q(i,k)) qc(i,k) = qc(i,k)/(1.0_kind_phys-q(i,k)) ENDDO ENDDO !aligo !--------------------------------------------------------------------- CALL FER_HIRES( & DT=DT,RHgrd=RHGRD & ,PRSI=prsi,P_PHY=p_phy,T_PHY=t & ,Q=Q,QT=cwm & ,LOWLYR=LOWLYR,SR=SR,TRAIN_PHY=train_phy & ,F_ICE_PHY=F_ICE,F_RAIN_PHY=F_RAIN & ,F_RIMEF_PHY=F_RIMEF & ,QC=QC,QR=QR,QS=QI & ,RAINNC=rainnc,RAINNCV=rainncv & ,threads=threads & ,IMS=IMS,IME=IME,LM=LM & ,D_SS=d_ss & ,refl_10cm=refl_10cm,DX1=DX1) !....................................................................... !MZ* !Aligo Oct-23-2019 ! - Convert dry qc,qr,qi back to wet mixing ratio DO K = 1, LM DO I= IMS, IME qc(i,k) = qc(i,k)/(1.0_kind_phys+q(i,k)) qi(i,k) = qi(i,k)/(1.0_kind_phys+q(i,k)) qr(i,k) = qr(i,k)/(1.0_kind_phys+q(i,k)) ENDDO ENDDO !----------------------------------------------------------- DO K=1,LM DO I=IMS,IME !--------------------------------------------------------------------- !*** Calculate graupel from total ice array and rime factor !--------------------------------------------------------------------- !MZ IF (SPEC_ADV) then QG(I,K)=QI(I,K)*F_RIMEF(I,K) ENDIF ! !----------------------------------------------------------------------- !*** UPDATE TEMPERATURE, SPECIFIC HUMIDITY, CLOUD WATER, AND HEATING. !----------------------------------------------------------------------- ! TRAIN(I,K)=TRAIN(I,K)+TRAIN_PHY(I,K) ENDDO ENDDO !....................................................................... ! !----------------------------------------------------------------------- !*** UPDATE PRECIPITATION !----------------------------------------------------------------------- ! DO I=IMS,IME PCPCOL=RAINNCV(I)*1.E-3 !MZ:unit:m PREC(I)=PREC(I)+PCPCOL !MZ ACPREC(I)=ACPREC(I)+PCPCOL !MZ: not used ! ! NOTE: RAINNC IS ACCUMULATED INSIDE MICROPHYSICS BUT NMM ZEROES IT OUT ABOVE ! SINCE IT IS ONLY A LOCAL ARRAY FOR NOW ! ENDDO !----------------------------------------------------------------------- ! end subroutine mp_fer_hires_run !> \section arg_table_mp_fer_hires_finalize Argument Table !! \htmlinclude mp_fer_hires_finalize.html !! subroutine mp_fer_hires_finalize (errmsg,errflg) implicit none character(len=*), intent( out) :: errmsg integer, intent( out) :: errflg ! Initialize the CCPP error handling variables errmsg = '' errflg = 0 if (.not.is_initialized) return call ferhires_finalize() is_initialized = .false. end subroutine mp_fer_hires_finalize end module mp_fer_hires