CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
m_micro_pre.F90
1
5
6 implicit none
7
8 contains
9
10!! \section arg_table_m_micro_pre_run Argument Table
11!! \htmlinclude m_micro_pre_run.html
12!!
13 subroutine m_micro_pre_run (im, levs, do_shoc, skip_macro, fprcp, mg3_as_mg2, gq0_ice, gq0_water, gq0_rain, &
14 gq0_snow, gq0_graupel, gq0_rain_nc, gq0_snow_nc, gq0_graupel_nc, cld_shoc, cnvc, cnvw, tcr, tcrf, gt0, &
15 qrn, qsnw, qgl, ncpr, ncps, ncgl, cld_frc_MG, clw_water, clw_ice, clcn, errmsg, errflg )
16
17 use machine, only : kind_phys
18 implicit none
19
20 integer, intent(in) :: im, levs, fprcp
21 logical, intent(in) :: do_shoc, mg3_as_mg2
22 logical, intent(inout) :: skip_macro
23 real(kind=kind_phys), intent(in) :: tcr, tcrf
24
25 real(kind=kind_phys), intent(in) :: &
26 gq0_ice(:,:), gq0_water(:,:), gq0_rain(:,:), gq0_snow(:,:), &
27 gq0_graupel(:,:), gq0_rain_nc(:,:), gq0_snow_nc(:,:), &
28 gq0_graupel_nc(:,:), cnvc(:,:), cnvw(:,:), gt0(:,:)
29 real(kind=kind_phys), intent(in), optional :: cld_shoc(:,:)
30 real(kind=kind_phys), intent(inout), optional :: &
31 qrn(:,:), qsnw(:,:), qgl(:,:), ncpr(:,:), ncps(:,:), ncgl(:,:), &
32 cld_frc_mg(:,:)
33
34 real(kind=kind_phys), intent(out) :: clw_ice(:,:), clw_water(:,:)
35
36 real(kind=kind_phys), intent(in), optional :: clcn(:,:)
37
38 character(len=*), intent(out) :: errmsg
39 integer, intent(out) :: errflg
40
41 integer :: i, k
42 real(kind=kind_phys) :: tem
43
44 ! Initialize CCPP error handling variables
45 errmsg = ''
46 errflg = 0
47
48 ! Acheng used clw here for other code to run smoothly and minimum change
49 ! to make the code work. However, the nc and clw should be treated
50 ! in other procceses too. August 28/2015; Hope that can be done next
51 ! year. I believe this will make the physical interaction more reasonable
52 ! Anning 12/5/2015 changed ntcw hold liquid only
53 skip_macro = do_shoc
54 if (do_shoc) then
55 if (fprcp == 0) then
56 do k=1,levs
57 do i=1,im
58 clw_ice(i,k) = gq0_ice(i,k)
59 clw_water(i,k) = gq0_water(i,k)
60 cld_frc_mg(i,k) = cld_shoc(i,k)
61 enddo
62 enddo
63 else if ((abs(fprcp) == 1) .or. mg3_as_mg2) then
64 do k=1,levs
65 do i=1,im
66 clw_ice(i,k) = gq0_ice(i,k)
67 clw_water(i,k) = gq0_water(i,k)
68 qrn(i,k) = gq0_rain(i,k)
69 qsnw(i,k) = gq0_snow(i,k)
70 ncpr(i,k) = gq0_rain_nc(i,k)
71 ncps(i,k) = gq0_snow_nc(i,k)
72 cld_frc_mg(i,k) = cld_shoc(i,k)
73 enddo
74 enddo
75 else
76 do k=1,levs
77 do i=1,im
78 clw_ice(i,k) = gq0_ice(i,k)
79 clw_water(i,k) = gq0_water(i,k)
80 qrn(i,k) = gq0_rain(i,k)
81 qsnw(i,k) = gq0_snow(i,k)
82 qgl(i,k) = gq0_graupel(i,k)
83 ncpr(i,k) = gq0_rain_nc(i,k)
84 ncps(i,k) = gq0_snow_nc(i,k)
85 ncgl(i,k) = gq0_graupel_nc(i,k)
86 cld_frc_mg(i,k) = cld_shoc(i,k)
87 enddo
88 enddo
89 end if
90 else
91 if (fprcp == 0 ) then
92 do k=1,levs
93 do i=1,im
94 clw_ice(i,k) = gq0_ice(i,k)
95 clw_water(i,k) = gq0_water(i,k)
96 enddo
97 enddo
98 elseif (abs(fprcp) == 1 .or. mg3_as_mg2) then
99 do k=1,levs
100 do i=1,im
101 clw_ice(i,k) = gq0_ice(i,k)
102 clw_water(i,k) = gq0_water(i,k)
103 qrn(i,k) = gq0_rain(i,k)
104 qsnw(i,k) = gq0_snow(i,k)
105 ncpr(i,k) = gq0_rain_nc(i,k)
106 ncps(i,k) = gq0_snow_nc(i,k)
107 enddo
108 enddo
109 else
110 do k=1,levs
111 do i=1,im
112 clw_ice(i,k) = gq0_ice(i,k)
113 clw_water(i,k) = gq0_water(i,k)
114 qrn(i,k) = gq0_rain(i,k)
115 qsnw(i,k) = gq0_snow(i,k)
116 qgl(i,k) = gq0_graupel(i,k)
117 ncpr(i,k) = gq0_rain_nc(i,k)
118 ncps(i,k) = gq0_snow_nc(i,k)
119 ncgl(i,k) = gq0_graupel_nc(i,k)
120 enddo
121 enddo
122 endif
123 end if
124
125 ! add convective cloud fraction
126 do k = 1,levs
127 do i = 1,im
128 cld_frc_mg(i,k) = min(1.0, cld_frc_mg(i,k) + clcn(i,k))
129 enddo
130 enddo
131
132 end subroutine m_micro_pre_run
133
134 end module m_micro_pre