CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
m_micro_post.F90
1
4
6
7 implicit none
8
9 contains
10
11!! \section arg_table_m_micro_post_run Argument Table
12!! \htmlinclude m_micro_post_run.html
13!!
14 subroutine m_micro_post_run( &
15 im, levs, fprcp, mg3_as_mg2, ncpr, ncps, ncgl, qrn, qsnw, qgl, &
16 gq0_ice, gq0_rain, gq0_snow, gq0_graupel, gq0_rain_nc, gq0_snow_nc, &
17 gq0_graupel_nc, ice, snow, graupel, dtp, errmsg, errflg)
18
19 use machine, only : kind_phys
20 implicit none
21
22 integer, intent(in) :: im, levs, fprcp
23 logical, intent(in) :: mg3_as_mg2
24
25 real(kind=kind_phys), intent(in ),optional :: ncpr(:,:)
26 real(kind=kind_phys), intent(in ),optional :: ncps(:,:)
27 real(kind=kind_phys), intent(in ),optional :: ncgl(:,:)
28 real(kind=kind_phys), intent(inout),optional :: qrn(:,:)
29 real(kind=kind_phys), intent(inout),optional :: qsnw(:,:)
30 real(kind=kind_phys), intent(inout),optional :: qgl(:,:)
31 real(kind=kind_phys), intent(in ) :: gq0_ice(:,:)
32 real(kind=kind_phys), intent(out ) :: gq0_rain(:,:)
33 real(kind=kind_phys), intent(out ) :: gq0_snow(:,:)
34 real(kind=kind_phys), intent(out ) :: gq0_graupel(:,:)
35 real(kind=kind_phys), intent(out ) :: gq0_rain_nc(:,:)
36 real(kind=kind_phys), intent(out ) :: gq0_snow_nc(:,:)
37 real(kind=kind_phys), intent(out ) :: gq0_graupel_nc(:,:)
38 real(kind=kind_phys), intent( out) :: ice(:)
39 real(kind=kind_phys), intent( out) :: snow(:)
40 real(kind=kind_phys), intent( out) :: graupel(:)
41 real(kind=kind_phys), intent(in ) :: dtp
42
43 character(len=*), intent(out) :: errmsg
44 integer, intent(out) :: errflg
45
46 ! Local variables
47 real(kind=kind_phys), parameter :: qsmall = 1.0d-20
48 real(kind=kind_phys), parameter :: con_p001 = 0.001d0
49 real(kind=kind_phys), parameter :: con_day = 86400.0d0
50 integer :: i, k
51 real(kind=kind_phys) :: tem
52
53 ! Initialize CCPP error handling variables
54 errmsg = ''
55 errflg = 0
56! do k=1,levs
57! write(1000+me,*)' maxwatnca=',maxval(Stateout%gq0(1:im,k,ntlnc)),' k=',k,' kdt=',kdt
58! enddo
59! write(1000+me,*)' at latitude = ',lat
60! tx1 = 1000.0
61! call moist_bud(im,ix,ix,levs,me,kdt,con_g,tx1,del,rain1
62! &, txa, clw(1,1,2), clw(1,1,1)
63! &, gq0(1,1,1),gq0(1,1,ntcw),gq0(1,1,ntcw+1),' m_micro ')
64
65! if (lprnt) write(0,*) ' rain1=',rain1(ipr)*86400.0, &
66! &' rainc=',diag%rainc(ipr)*86400.0
67! &,' cn_prc=',cn_prc(ipr),' cn_snr=',cn_snr(ipr)
68! if(lprnt) write(0,*) ' aftgt0=',Stateout%gt0(ipr,:),' kdt=',kdt
69! if (lprnt) write(0,*) ' aftlsgq0=',stateout%gq0(ipr,:,1),' kdt=',kdt
70! if (lprnt) write(0,*)' clw1aft=',stateout%gq0(ipr,:,ntiw),' kdt=',kdt
71! if (ntgl > 0 .and. lprnt) &
72! write(0,*)' cgw1aft=',stateout%gq0(ipr,:,ntgl),' kdt=',kdt
73! if (lprnt) write(0,*)' cloudsm=',tbd%phy_f3d(ipr,:,1)*100,' kdt=',kdt
74! if (lprnt) write(0,*)' clw2aft=',stateout%gq0(ipr,:,ntcw),' kdt=',kdt
75! if (lprnt) write(0,*)' qrna=',qrn(ipr,:),' kdt=',kdt
76! if (lprnt) write(0,*)' qsnwa=',qsnw(ipr,:),' kdt=',kdt
77! if (lprnt) write(0,*)' qglba',qgl(ipr,:),' kdt=',kdt
78
79 tem = dtp * con_p001 / con_day
80 if (abs(fprcp) == 1 .or. mg3_as_mg2) then
81 do k=1,levs
82 do i=1,im
83 if (abs(qrn(i,k)) < qsmall) qrn(i,k) = 0.0
84 if (abs(qsnw(i,k)) < qsmall) qsnw(i,k) = 0.0
85 gq0_rain(i,k) = qrn(i,k)
86 gq0_snow(i,k) = qsnw(i,k)
87 gq0_rain_nc(i,k) = ncpr(i,k)
88 gq0_snow_nc(i,k) = ncps(i,k)
89 enddo
90 enddo
91 do i=1,im
92 ice(i) = tem * gq0_ice(i,1)
93 snow(i) = tem * qsnw(i,1)
94 enddo
95 elseif (fprcp > 1) then
96 do k=1,levs
97 do i=1,im
98 if (abs(qrn(i,k)) < qsmall) qrn(i,k) = 0.0
99 if (abs(qsnw(i,k)) < qsmall) qsnw(i,k) = 0.0
100 if (abs(qgl(i,k)) < qsmall) qgl(i,k) = 0.0
101 gq0_rain(i,k) = qrn(i,k)
102 gq0_snow(i,k) = qsnw(i,k)
103 gq0_graupel(i,k) = qgl(i,k)
104 gq0_rain_nc(i,k) = ncpr(i,k)
105 gq0_snow_nc(i,k) = ncps(i,k)
106 gq0_graupel_nc(i,k) = ncgl(i,k)
107 enddo
108 enddo
109 do i=1,im
110 ice(i) = tem * gq0_ice(i,1)
111 snow(i) = tem * qsnw(i,1)
112 graupel(i) = tem * qgl(i,1)
113 enddo
114
115 endif
116
117! if (lprnt) write(0,*)' cloudsm=',tbd%phy_f3d(ipr,:,1)*100,' kdt=',kdt
118! if (lprnt) write(0,*)' clw2aft=',stateout%gq0(ipr,:,ntcw),' kdt=',kdt
119! if (lprnt) write(0,*)' qrna=',qrn(ipr,:),' kdt=',kdt
120! if (lprnt) write(0,*)' qsnwa=',qsnw(ipr,:),' kdt=',kdt
121! if (lprnt) write(0,*)' qglba',qgl(ipr,:),' kdt=',kdt
122!
123
124
125 end subroutine m_micro_post_run
126
127 end module m_micro_post