34 subroutine mp_fer_hires_init(ncol, nlev, dtp, imp_physics, &
35 imp_physics_fer_hires, restart, &
36 mpicomm, mpirank,mpiroot, &
37 threads, errmsg, errflg)
44 integer,
intent(in) :: ncol
45 integer,
intent(in) :: nlev
46 real(kind_phys),
intent(in) :: dtp
47 integer,
intent(in) :: imp_physics
48 integer,
intent(in) :: imp_physics_fer_hires
49 type(mpi_comm),
intent(in) :: mpicomm
50 integer,
intent(in) :: mpirank
51 integer,
intent(in) :: mpiroot
52 integer,
intent(in) :: threads
53 logical,
intent(in) :: restart
54 character(len=*),
intent(out) :: errmsg
55 integer,
intent(out) :: errflg
58 integer :: ims, ime, lm,i,k
65 if (is_initialized)
return
73 if (mpirank==mpiroot)
then
74 write(0,*)
' -----------------------------------------------'
75 write(0,*)
' --- !!! WARNING !!! ---'
76 write(0,*)
' --- the CCPP Ferrier-Aligo MP scheme is ---'
77 write(0,*)
' --- currently under development, use at ---'
78 write(0,*)
' --- your own risk . ---'
79 write(0,*)
' -----------------------------------------------'
83 if (imp_physics /= imp_physics_fer_hires )
then
84 write(errmsg,
'(*(a))')
"Logic error: namelist choice of microphysics is different from Ferrier-Aligo MP"
89 if (mpirank==mpiroot)
write (0,*)
'F-A: calling FERRIER_INIT_HR ...'
92 if (mpirank==mpiroot)
write (0,*)
'F-A: FERRIER_INIT_HR finished ...'
93 if (errflg /= 0 )
return
95 is_initialized = .true.
104 SUBROUTINE mp_fer_hires_run(NCOL, NLEV, DT ,SPEC_ADV &
111 ,mpirank, mpiroot, threads &
114 ,EPSQ,R_D,P608,CP,G &
124 INTEGER,
PARAMETER :: d_ss=1
130 integer,
intent(in ) :: ncol
131 integer,
intent(in ) :: nlev
132 real(kind_phys),
intent(in ) :: dt
133 integer,
intent(in ) :: threads
134 logical,
intent(in ) :: spec_adv
135 integer,
intent(in ) :: mpirank
136 integer,
intent(in ) :: mpiroot
137 real(kind_phys),
intent(in ) :: slmsk(:)
138 real(kind_phys),
intent(in ) :: prsi(:,:)
139 real(kind_phys),
intent(in ) :: p_phy(:,:)
140 real(kind_phys),
intent(in ) :: epsq,r_d,p608,cp,g
141 real(kind_phys),
intent(inout) :: t(:,:)
142 real(kind_phys),
intent(inout) :: q(:,:)
143 real(kind_phys),
intent(inout),
optional :: train(:,:)
144 real(kind_phys),
intent(out ) :: sr(:)
145 real(kind_phys),
intent(inout) :: qc(:,:)
146 real(kind_phys),
intent(inout) :: qr(:,:)
147 real(kind_phys),
intent(inout) :: qi(:,:)
148 real(kind_phys),
intent(inout) :: qg(:,:)
150 real(kind_phys),
intent(inout) :: prec(:)
151 real(kind_phys),
intent(inout) :: refl_10cm(:,:)
152 real(kind_phys),
intent(in ),
optional :: rhgrd
153 real(kind_phys),
intent(in ) :: dx(:)
154 character(len=*),
intent(out) :: errmsg
155 integer,
intent(out) :: errflg
162 integer :: lowlyr(1:ncol)
164 real(kind_phys) :: pcpcol
165 real(kind_phys) :: ql(1:nlev),tl(1:nlev)
166 real(kind_phys) :: rainnc(1:ncol),rainncv(1:ncol)
167 real(kind_phys) :: snownc(1:ncol),snowncv(1:ncol)
168 real(kind_phys) :: graupelncv(1:ncol)
169 real(kind_phys) :: train_phy(1:ncol,1:nlev)
170 real(kind_phys) :: f_ice(1:ncol,1:nlev)
171 real(kind_phys) :: f_rain(1:ncol,1:nlev)
172 real(kind_phys) :: f_rimef(1:ncol,1:nlev)
173 real(kind_phys) :: cwm(1:ncol,1:nlev)
176 integer :: ims, ime, lm
186 if (.not. is_initialized)
then
187 write(errmsg, fmt=
'((a))')
'mp_fer_hires_run called before mp_fer_hires_init'
238 cwm(i,k)=qc(i,k)+qr(i,k)+qi(i,k)
239 IF (qi(i,k) <= epsq)
THEN
242 IF (t(i,k) < t_icek) f_ice(i,k)=1.
244 f_ice(i,k)=max( 0., min(1., qi(i,k)/cwm(i,k) ) )
247 IF (qr(i,k) <= epsq)
THEN
250 f_rain(i,k)=qr(i,k)/(qr(i,k)+qc(i,k))
261 cwm(i,k) = cwm(i,k)/(1.0_kind_phys-q(i,k))
262 qr(i,k) = qr(i,k)/(1.0_kind_phys-q(i,k))
263 qi(i,k) = qi(i,k)/(1.0_kind_phys-q(i,k))
264 qc(i,k) = qc(i,k)/(1.0_kind_phys-q(i,k))
272 ,prsi=prsi,p_phy=p_phy,t_phy=t &
274 ,lowlyr=lowlyr,sr=sr,train_phy=train_phy &
275 ,f_ice_phy=f_ice,f_rain_phy=f_rain &
276 ,f_rimef_phy=f_rimef &
278 ,rainnc=rainnc,rainncv=rainncv &
280 ,ims=ims,ime=ime,lm=lm &
282 ,refl_10cm=refl_10cm,dx1=dx1)
292 qc(i,k) = qc(i,k)/(1.0_kind_phys+q(i,k))
293 qi(i,k) = qi(i,k)/(1.0_kind_phys+q(i,k))
294 qr(i,k) = qr(i,k)/(1.0_kind_phys+q(i,k))
308 qg(i,k)=qi(i,k)*f_rimef(i,k)
316 train(i,k)=train(i,k)+train_phy(i,k)
328 pcpcol=rainncv(i)*1.e-3
329 prec(i)=prec(i)+pcpcol
subroutine fer_hires(dt, rhgrd, prsi, p_phy, t_phy, q, qt, lowlyr, sr, train_phy, f_ice_phy, f_rain_phy, f_rimef_phy, qc, qr, qs, rainnc, rainncv, threads, ims, ime, lm, d_ss, refl_10cm, dx1)
This is the driver scheme of Ferrier-Aligo microphysics scheme. NOTE: The only differences between FE...