CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_surface_generic_pre.F90
1
3
5
6 use machine, only: kind_phys
7
8 implicit none
9
10 private
11
13
14 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys, one = 1.0_kind_phys
15
16 contains
17
24 subroutine gfs_surface_generic_pre_init (nthreads, im, slmsk, isot, ivegsrc, stype,scolor, vtype, slope, &
25 vtype_save, stype_save,scolor_save, slope_save, errmsg, errflg)
26
27 implicit none
28
29 ! Interface variables
30 integer, intent(in) :: nthreads, im, isot, ivegsrc
31 real(kind_phys), dimension(:), intent(in) :: slmsk
32 integer, dimension(:), intent(inout) :: vtype, stype, scolor,slope
33 integer, dimension(:), intent(out) :: vtype_save, stype_save,scolor_save, slope_save
34
35 ! CCPP error handling
36 character(len=*), intent(out) :: errmsg
37 integer, intent(out) :: errflg
38
39 ! Local variables
40 integer, dimension(1:im) :: islmsk
41 integer :: i
42
43 ! Initialize CCPP error handling variables
44 errmsg = ''
45 errflg = 0
46
47 islmsk = nint(slmsk)
48
49 ! Save current values of vegetation, soil and slope type
50 vtype_save(:) = vtype(:)
51 stype_save(:) = stype(:)
52 scolor_save(:) = scolor(:)
53 slope_save(:) = slope(:)
54
55 call update_vegetation_soil_slope_type(nthreads, im, isot, ivegsrc, islmsk, vtype, stype,scolor, slope)
56
57 end subroutine gfs_surface_generic_pre_init
58
62 subroutine gfs_surface_generic_pre_run (nthreads, im, levs, vfrac, islmsk, isot, ivegsrc, stype, scolor,vtype, slope, &
63 prsik_1, prslk_1, tsfc, phil, con_g, sigmaf, work3, zlvl, &
64 lndp_type, n_var_lndp, sfc_wts, lndp_var_list, lndp_prt_list, &
65 z01d, zt1d, bexp1d, xlai1d, vegf1d, lndp_vgf, &
66 cplflx, flag_cice, islmsk_cice, slimskin_cpl, &
67 wind, u1, v1, cnvwind, smcwlt2, smcref2, vtype_save, stype_save,scolor_save, slope_save, &
68 errmsg, errflg)
69
71
72 implicit none
73
74 ! Interface variables
75 integer, intent(in) :: nthreads, im, levs, isot, ivegsrc
76 integer, dimension(:), intent(in) :: islmsk
77
78 real(kind=kind_phys), intent(in) :: con_g
79 real(kind=kind_phys), dimension(:), intent(in) :: vfrac, prsik_1, prslk_1
80 integer, dimension(:), intent(inout) :: vtype, stype,scolor, slope
81 integer, dimension(:), intent(out) :: vtype_save(:), stype_save(:),scolor_save(:), slope_save(:)
82
83 real(kind=kind_phys), dimension(:), intent(inout) :: tsfc
84 real(kind=kind_phys), dimension(:,:), intent(in) :: phil
85
86 real(kind=kind_phys), dimension(:), intent(inout) :: sigmaf, work3, zlvl
87
88 ! Stochastic physics / surface perturbations
89 integer, intent(in) :: lndp_type, n_var_lndp
90 character(len=3), dimension(:), intent(in), optional :: lndp_var_list
91 real(kind=kind_phys), dimension(:), intent(in), optional :: lndp_prt_list
92 real(kind=kind_phys), dimension(:,:), intent(in), optional :: sfc_wts
93 real(kind=kind_phys), dimension(:), intent(out) :: z01d
94 real(kind=kind_phys), dimension(:), intent(out) :: zt1d
95 real(kind=kind_phys), dimension(:), intent(out) :: bexp1d
96 real(kind=kind_phys), dimension(:), intent(out) :: xlai1d
97 real(kind=kind_phys), dimension(:), intent(out) :: vegf1d
98 real(kind=kind_phys), intent(out) :: lndp_vgf
99
100 logical, intent(in) :: cplflx
101 real(kind=kind_phys), dimension(:), intent(in), optional :: slimskin_cpl
102 logical, dimension(:), intent(inout) :: flag_cice
103 integer, dimension(:), intent(out) :: islmsk_cice
104
105 real(kind=kind_phys), dimension(:), intent(out) :: wind
106 real(kind=kind_phys), dimension(:), intent(in ) :: u1, v1
107 ! surface wind enhancement due to convection
108 real(kind=kind_phys), dimension(:), intent(inout ), optional :: cnvwind
109 !
110 real(kind=kind_phys), dimension(:), intent(out) :: smcwlt2, smcref2
111
112 ! CCPP error handling
113 character(len=*), intent(out) :: errmsg
114 integer, intent(out) :: errflg
115
116 ! Local variables
117 integer :: i, k
118 real(kind=kind_phys) :: onebg, cdfz
119
120 ! Set constants
121 onebg = 1.0/con_g
122
123 ! Initialize CCPP error handling variables
124 errmsg = ''
125 errflg = 0
126
127 ! Scale random patterns for surface perturbations with perturbation size
128 ! Turn vegetation fraction pattern into percentile pattern
129 lndp_vgf=-999.
130
131 if (lndp_type==1) then
132 do k =1,n_var_lndp
133 select case(lndp_var_list(k))
134 case ('rz0')
135 z01d(:) = lndp_prt_list(k)* sfc_wts(:,k)
136 case ('rzt')
137 zt1d(:) = lndp_prt_list(k)* sfc_wts(:,k)
138 case ('shc')
139 bexp1d(:) = lndp_prt_list(k) * sfc_wts(:,k)
140 case ('lai')
141 xlai1d(:) = lndp_prt_list(k)* sfc_wts(:,k)
142 case ('vgf')
143 ! note that the pertrubed vegfrac is being used in sfc_drv, but not sfc_diff
144 do i=1,im
145 call cdfnor(sfc_wts(i,k),cdfz)
146 vegf1d(i) = cdfz
147 enddo
148 lndp_vgf = lndp_prt_list(k)
149 end select
150 enddo
151 endif
152
153 ! End of stochastic physics / surface perturbation
154
155 ! Save current values of vegetation, soil and slope type
156 vtype_save(:) = vtype(:)
157 stype_save(:) = stype(:)
158 scolor_save(:) = scolor(:)
159 slope_save(:) = slope(:)
160
161 call update_vegetation_soil_slope_type(nthreads, im, isot, ivegsrc, islmsk, vtype, stype,scolor, slope)
162
163 do i=1,im
164 sigmaf(i) = max(vfrac(i), 0.01_kind_phys)
165 islmsk_cice(i) = islmsk(i)
166
167 work3(i) = prsik_1(i) / prslk_1(i)
168
169 zlvl(i) = phil(i,1) * onebg
170 smcwlt2(i) = zero
171 smcref2(i) = zero
172
173 wind(i) = max(sqrt(u1(i)*u1(i) + v1(i)*v1(i)) &
174 + max(zero, min(cnvwind(i), 30.0_kind_phys)), one)
175 !wind(i) = max(sqrt(Statein%ugrs(i,1)*Statein%ugrs(i,1) + &
176 ! Statein%vgrs(i,1)*Statein%vgrs(i,1)) &
177 ! + max(zero, min(Tbd%phy_f2d(i,Model%num_p2d), 30.0)), one)
178 cnvwind(i) = zero
179
180 enddo
181
182 if (cplflx) then
183 do i=1,im
184 islmsk_cice(i) = nint(slimskin_cpl(i))
185 flag_cice(i) = (islmsk_cice(i) == 4)
186 enddo
187 endif
188
189 end subroutine gfs_surface_generic_pre_run
190
191 subroutine update_vegetation_soil_slope_type(nthreads, im, isot, ivegsrc, islmsk, vtype, stype,scolor, slope)
192
193 implicit none
194
195 integer, intent(in) :: nthreads, im, isot, ivegsrc, islmsk(:)
196 integer, intent(inout) :: vtype(:), stype(:),scolor(:), slope(:)
197 integer :: i
198
199!$OMP parallel do num_threads(nthreads) default(none) private(i) &
200!$OMP shared(im, isot, ivegsrc, islmsk, vtype, stype,scolor, slope)
201
202! scolor is a place holder now, how to update soil color based on the mask/veg/sot src
203
204 do i=1,im
205 if (islmsk(i) == 2) then
206 if (isot == 1) then
207 stype(i) = 16
208 else
209 stype(i) = 9
210 endif
211 if (ivegsrc == 0 .or. ivegsrc == 4) then
212 vtype(i) = 24
213 elseif (ivegsrc == 1) then
214 vtype(i) = 15
215 elseif (ivegsrc == 2) then
216 vtype(i) = 13
217 elseif (ivegsrc == 3 .or. ivegsrc == 5) then
218 vtype(i) = 15
219 endif
220 slope(i) = 9
221 else
222 if (vtype(i) < 1) vtype(i) = 17
223 if (slope(i) < 1) slope(i) = 1
224 endif
225 enddo
226!$OMP end parallel do
227
230
231 end module gfs_surface_generic_pre
subroutine, public cdfnor(z, cdfz)
This subrtouine calculates the CDF of the standard normal distribution evaluated at z.
subroutine, public gfs_surface_generic_pre_run(nthreads, im, levs, vfrac, islmsk, isot, ivegsrc, stype, scolor, vtype, slope, prsik_1, prslk_1, tsfc, phil, con_g, sigmaf, work3, zlvl, lndp_type, n_var_lndp, sfc_wts, lndp_var_list, lndp_prt_list, z01d, zt1d, bexp1d, xlai1d, vegf1d, lndp_vgf, cplflx, flag_cice, islmsk_cice, slimskin_cpl, wind, u1, v1, cnvwind, smcwlt2, smcref2, vtype_save, stype_save, scolor_save, slope_save, errmsg, errflg)
subroutine, public gfs_surface_generic_pre_init(nthreads, im, slmsk, isot, ivegsrc, stype, scolor, vtype, slope, vtype_save, stype_save, scolor_save, slope_save, errmsg, errflg)
subroutine update_vegetation_soil_slope_type(nthreads, im, isot, ivegsrc, islmsk, vtype, stype, scolor, slope)