CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
mp_fer_hires.F90
1
3
4!
6
7 use machine, only : kind_phys
8
10 ferhires_finalize
11
12 implicit none
13
14 public :: mp_fer_hires_init, mp_fer_hires_run, mp_fer_hires_finalize
15
16 private
17
18 logical :: is_initialized = .false.
19
20 ! * T_ICE - temperature (C) threshold at which all remaining liquid water
21 ! is glaciated to ice
22 ! * T_ICE_init - maximum temperature (C) at which ice nucleation occurs
23 REAL, PUBLIC, PARAMETER :: t_ice=-40., &
24 t0c=273.15, &
25 t_icek=t0c+t_ice
26
27 contains
28
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)
38
39 USE mpi_f08
40 USE machine, ONLY : kind_phys
42 implicit none
43
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
56
57 ! Local variables
58 integer :: ims, ime, lm,i,k
59 !real(kind=kind_phys) :: DT_MICRO
60
61 ! Initialize the CCPP error handling variables
62 errmsg = ''
63 errflg = 0
64
65 if (is_initialized) return
66
67 ! Set internal dimensions
68 ims = 1
69 ime = ncol
70 lm = nlev
71
72 ! MZ* temporary
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,*) ' -----------------------------------------------'
80 end if
81 ! MZ* temporary
82
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"
85 errflg = 1
86 return
87 end if
88
89 if (mpirank==mpiroot) write (0,*) 'F-A: calling FERRIER_INIT_HR ...'
90 CALL ferrier_init_hr(dtp,mpicomm,mpirank,mpiroot,threads,errmsg,errflg)
91
92 if (mpirank==mpiroot) write (0,*)'F-A: FERRIER_INIT_HR finished ...'
93 if (errflg /= 0 ) return
94
95 is_initialized = .true.
96
97 end subroutine mp_fer_hires_init
98
104 SUBROUTINE mp_fer_hires_run(NCOL, NLEV, DT ,SPEC_ADV &
105 ,SLMSK &
106 ,PRSI,P_PHY &
107 ,T,Q &
108 ,TRAIN,SR &
109 ,QC,QR,QI,QG &
110 ,PREC &
111 ,mpirank, mpiroot, threads &
112 ,refl_10cm &
113 ,RHGRD,dx &
114 ,EPSQ,R_D,P608,CP,G &
115 ,errmsg,errflg)
116
117!-----------------------------------------------------------------------
118 USE machine, ONLY: kind_phys
119!
120 IMPLICIT NONE
121!
122!-----------------------------------------------------------------------
123!
124 INTEGER,PARAMETER :: d_ss=1
125!
126!------------------------
127!*** Argument Variables
128!------------------------
129
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(:,:) ! QRIMEF
149
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
156!
157!---------------------
158!*** Local Variables
159!---------------------
160!
161 integer :: i,j,k,n
162 integer :: lowlyr(1:ncol)
163 integer :: dx1
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)
174
175! Dimension
176 integer :: ims, ime, lm
177
178!-----------------------------------------------------------------------
179!***********************************************************************
180!-----------------------------------------------------------------------
181 ! Initialize the CCPP error handling variables
182 errmsg = ''
183 errflg = 0
184
185 ! Check initialization state
186 if (.not. is_initialized) then
187 write(errmsg, fmt='((a))') 'mp_fer_hires_run called before mp_fer_hires_init'
188 errflg = 1
189 return
190 end if
191
192! Set internal dimensions
193 ims = 1
194 ime = ncol
195 lm = nlev
196
197! Use the dx of the 1st i point to set an integer value of dx to be used for
198! determining where RHgrd should be set to 0.98 in the coarse domain when running HAFS.
199 dx1=nint(dx(1))
200
201!-----------------------------------------------------------------------
202!*** NOTE: THE NMMB HAS IJK STORAGE WITH LAYER 1 AT THE TOP.
203!*** THE WRF PHYSICS DRIVERS HAVE IKJ STORAGE WITH LAYER 1
204!*** AT THE BOTTOM.
205!-----------------------------------------------------------------------
206!.......................................................................
207 DO i=ims,ime
208!
209 lowlyr(i)=1
210!
211!-----------------------------------------------------------------------
212!*** FILL RAINNC WITH ZERO (NORMALLY CONTAINS THE NONCONVECTIVE
213!*** ACCUMULATED RAIN BUT NOT YET USED BY NMM)
214!*** COULD BE OBTAINED FROM ACPREC AND CUPREC (ACPREC-CUPREC)
215!-----------------------------------------------------------------------
216!..The NC variables were designed to hold simulation total accumulations
217!.. whereas the NCV variables hold timestep only values, so change below
218!.. to zero out only the timestep amount preparing to go into each
219!.. micro routine while allowing NC vars to accumulate continually.
220!.. But, the fact is, the total accum variables are local, never saved
221!.. nor written so they go nowhere at the moment.
222!
223 rainnc(i)=0. ! NOT YET USED BY NMM
224 rainncv(i)=0.
225 snowncv(i)=0.
226 graupelncv(i) = 0.0
227!
228!-----------------------------------------------------------------------
229!*** FILL THE SINGLE-COLUMN INPUT
230!-----------------------------------------------------------------------
231!
232 DO k=lm,1,-1 !mz* We are moving down from the top in the flipped arrays
233
234!*** CALL MICROPHYSICS
235
236!MZ* in HWRF
237!-- 6/11/2010: Update cwm, F_ice, F_rain and F_rimef arrays
238 cwm(i,k)=qc(i,k)+qr(i,k)+qi(i,k)
239 IF (qi(i,k) <= epsq) THEN
240 f_ice(i,k)=0.
241 f_rimef(i,k)=1.
242 IF (t(i,k) < t_icek) f_ice(i,k)=1.
243 ELSE
244 f_ice(i,k)=max( 0., min(1., qi(i,k)/cwm(i,k) ) )
245 f_rimef(i,k)=qg(i,k)!/QI(I,K)
246 ENDIF
247 IF (qr(i,k) <= epsq) THEN
248 f_rain(i,k)=0.
249 ELSE
250 f_rain(i,k)=qr(i,k)/(qr(i,k)+qc(i,k))
251 ENDIF
252
253 ENDDO
254
255 ENDDO
256
257!---------------------------------------------------------------------
258!aligo
259 DO k = 1, lm
260 DO i= ims, ime
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))
265 ENDDO
266 ENDDO
267!aligo
268!---------------------------------------------------------------------
269
270 CALL fer_hires( &
271 dt=dt,rhgrd=rhgrd &
272 ,prsi=prsi,p_phy=p_phy,t_phy=t &
273 ,q=q,qt=cwm &
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 &
277 ,qc=qc,qr=qr,qs=qi &
278 ,rainnc=rainnc,rainncv=rainncv &
279 ,threads=threads &
280 ,ims=ims,ime=ime,lm=lm &
281 ,d_ss=d_ss &
282 ,refl_10cm=refl_10cm,dx1=dx1)
283
284
285!.......................................................................
286
287!MZ*
288!Aligo Oct-23-2019
289! - Convert dry qc,qr,qi back to wet mixing ratio
290 DO k = 1, lm
291 DO i= ims, ime
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))
295 ENDDO
296 ENDDO
297
298!-----------------------------------------------------------
299 DO k=1,lm
300 DO i=ims,ime
301
302!---------------------------------------------------------------------
303!*** Calculate graupel from total ice array and rime factor
304!---------------------------------------------------------------------
305
306!MZ
307 IF (spec_adv) then
308 qg(i,k)=qi(i,k)*f_rimef(i,k)
309 ENDIF
310
311!
312!-----------------------------------------------------------------------
313!*** UPDATE TEMPERATURE, SPECIFIC HUMIDITY, CLOUD WATER, AND HEATING.
314!-----------------------------------------------------------------------
315!
316 train(i,k)=train(i,k)+train_phy(i,k)
317 ENDDO
318 ENDDO
319
320!.......................................................................
321
322!
323!-----------------------------------------------------------------------
324!*** UPDATE PRECIPITATION
325!-----------------------------------------------------------------------
326!
327 DO i=ims,ime
328 pcpcol=rainncv(i)*1.e-3 !MZ:unit:m
329 prec(i)=prec(i)+pcpcol
330!MZ ACPREC(I)=ACPREC(I)+PCPCOL !MZ: not used
331!
332! NOTE: RAINNC IS ACCUMULATED INSIDE MICROPHYSICS BUT NMM ZEROES IT OUT ABOVE
333! SINCE IT IS ONLY A LOCAL ARRAY FOR NOW
334!
335 ENDDO
336!-----------------------------------------------------------------------
337!
338 end subroutine mp_fer_hires_run
339
340
344 subroutine mp_fer_hires_finalize (errmsg,errflg)
345 implicit none
346
347 character(len=*), intent( out) :: errmsg
348 integer, intent( out) :: errflg
349
350 ! Initialize the CCPP error handling variables
351 errmsg = ''
352 errflg = 0
353
354 if (.not.is_initialized) return
355
356 call ferhires_finalize()
357
358 is_initialized = .false.
359 end subroutine mp_fer_hires_finalize
360
361end module mp_fer_hires
subroutine, public ferrier_init_hr(gsmdt, mpi_comm_comp, mpirank, mpiroot, threads, errmsg, errflg)
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...