CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
scm_sfc_flux_spec.F90
1
3
5
6 implicit none
7
8 private
9
10!----------------
11! Public entities
12!----------------
13 public scm_sfc_flux_spec_init, scm_sfc_flux_spec_run
14
15 CONTAINS
16!*******************************************************************************************
17
21 subroutine scm_sfc_flux_spec_init(lheatstrg, errmsg, errflg)
22
23 logical, intent(in) :: lheatstrg
24
25 character(len=*), intent(out) :: errmsg
26 integer, intent(out) :: errflg
27
28 if (lheatstrg) then
29 errmsg = 'Using specified surface fluxes is not compatible with canopy heat storage (lheatstrg) being true. Stopping.'
30 errflg = 1
31 return
32 end if
33 end subroutine scm_sfc_flux_spec_init
34
51 subroutine scm_sfc_flux_spec_run (im, u1, v1, z1, t1, q1, p1, roughness_length, spec_sh_flux, spec_lh_flux, &
52 exner_inverse, T_surf, cp, grav, hvap, rd, fvirt, vonKarman, tgice, islmsk, dry, frland, cice, icy, tisfc,&
53 oceanfrac, min_seaice, cplflx, cplice, flag_cice, wet, min_lakeice, tsfcl, tsfc_wat, slmsk, lakefrac, lkm,&
54 lakedepth, use_lake_model, sh_flux, lh_flux, sh_flux_chs, u_star, sfc_stress, cm, ch, &
55 fm, fh, rb, u10m, v10m, wind1, qss, t2m, q2m, errmsg, errflg)
56
57 use machine, only: kind_phys
58
59 integer, intent(in) :: im, lkm
60 integer, intent(inout) :: islmsk(:), use_lake_model(:)
61 logical, intent(in) :: cplflx, cplice
62 logical, intent(inout) :: dry(:), icy(:), flag_cice(:), wet(:)
63 real(kind=kind_phys), intent(in) :: cp, grav, hvap, rd, fvirt, vonkarman, min_seaice, tgice, min_lakeice
64 real(kind=kind_phys), intent(in) :: u1(:), v1(:), z1(:), t1(:), q1(:), p1(:), roughness_length(:), &
65 spec_sh_flux(:), spec_lh_flux(:), exner_inverse(:), t_surf(:), oceanfrac(:), lakefrac(:), lakedepth(:)
66 real(kind=kind_phys), intent(inout) :: cice(:), tisfc(:), tsfcl(:), tsfc_wat(:), slmsk(:)
67 real(kind=kind_phys), intent(out) :: sh_flux(:), lh_flux(:), u_star(:), sfc_stress(:), &
68 cm(:), ch(:), fm(:), fh(:), rb(:), u10m(:), v10m(:), wind1(:), qss(:), t2m(:), q2m(:), &
69 sh_flux_chs(:), frland(:)
70
71 character(len=*), intent(out) :: errmsg
72 integer, intent(out) :: errflg
73
74 integer :: i
75
76 real(kind=kind_phys) :: rho, q1_non_neg, w_thv1, rho_cp_inverse, rho_hvap_inverse, obukhov_length, thv1, tvs, &
77 dtv, adtv, wind10m, u_fraction, roughness_length_m
78
79 real(kind=kind_phys), parameter :: timin = 173.0_kind_phys ! minimum temperature allowed for snow/ice
80
81 ! Initialize CCPP error handling variables
82 errmsg = ''
83 errflg = 0
84
85! !--- set control properties (including namelist read)
86 !calculate u_star from wind profiles (need roughness length, and wind and height at lowest model level)
87 do i=1, im
88 sh_flux(i) = spec_sh_flux(i)
89 lh_flux(i) = spec_lh_flux(i)
90 sh_flux_chs(i) = sh_flux(i)
91
92 roughness_length_m = 0.01*roughness_length(i)
93
94 wind1(i) = sqrt(u1(i)*u1(i) + v1(i)*v1(i))
95 u_star(i) = vonkarman*wind1(i)/(log(z1(i)/roughness_length_m))
96
97 !calculate variables related to u_star
98 sfc_stress(i) = u_star(i)*u_star(i)
99
100 if(wind1(i) > 0.0) then
101 cm(i) = sfc_stress(i)/(wind1(i)*wind1(i))
102 else
103 cm(i) = 0.0
104 end if
105 fm(i) = sqrt((vonkarman*vonkarman)/cm(i))
106
107 !calculate the Obukhov length from the specified surface fluxes
108 q1_non_neg = max( q1(i), 1.0e-8 )
109 rho = p1(i) / (rd*t1(i)*(1.0 + fvirt*q1_non_neg))
110 rho_cp_inverse = 1.0/(rho*cp)
111 rho_hvap_inverse = 1.0/(rho*hvap)
112
113 thv1 = t1(i) * exner_inverse(i) * (1.0 + fvirt*q1_non_neg)
114 ! sh_flux = rho_cp_inverse*sh_flux_wm2(i)
115 ! lh_flux = rho_hvap_inverse*lh_flux_wm2(i)
116 w_thv1 = sh_flux(i)*exner_inverse(i) + fvirt*t1(i)*lh_flux(i)
117
118 obukhov_length = -u_star(i)*u_star(i)*u_star(i)*thv1/(vonkarman*grav*w_thv1)
119
120 !calculate the bulk Richardson number and the M-O function for scalars
121 tvs = t_surf(i)*(1.0 + fvirt*q1_non_neg)
122
123 dtv = thv1 - tvs
124 adtv = max(abs(dtv),0.001)
125 dtv = sign(1._kind_phys,dtv) * adtv
126 rb(i) = max(-5000.0, (grav+grav) * dtv * z1(i) / ((thv1 + tvs) * wind1(i) * wind1(i)))
127
128 fh(i) = rb(i)*fm(i)*fm(i)*obukhov_length/z1(i)
129 ch(i) = vonkarman*vonkarman/(fm(i)*fh(i))
130
131 !calculate sfc_diag variables (bypassing t2m and q2m for now)
132 !should calculate qss, but it is not needed in moninedmf (q2m depends on qss - will implement later)
133
134 wind10m = u_star(i)/vonkarman*log(10.0/roughness_length_m)
135 u_fraction = sqrt(1.0 - v1(i)**2/wind1(i)**2)
136 u10m(i) = u_fraction*wind10m
137 v10m(i) = sqrt(wind10m**2 - u10m(i)**2)
138
139 qss(i) = 0.0
140 t2m(i) = 0.0
141 q2m(i) = 0.0
142 end do
143
144 !GJF: The following code is from GFS_surface_composites.F90; only statements that are used in physics schemes outside of surface schemes are kept
145 !GJF: Adding this code means that this scheme should be called before dcyc2t3
146 do i = 1, im
147 if (islmsk(i) == 1) then
148 dry(i) = .true.
149 frland(i) = 1.0_kind_phys
150 cice(i) = 0.0_kind_phys
151 icy(i) = .false.
152 else
153 frland(i) = 0.0_kind_phys
154 if (oceanfrac(i) > 0.0_kind_phys) then
155 if (cice(i) >= min_seaice) then
156 icy(i) = .true.
157 ! This cplice namelist option was added to deal with the
158 ! situation of the FV3ATM-HYCOM coupling without an active sea
159 ! ice (e.g., CICE6) component. By default, the cplice is true
160 ! when cplflx is .true. (e.g., for the S2S application).
161 ! Whereas, for the HAFS FV3ATM-HYCOM coupling, cplice is set as
162 ! .false.. In the future HAFS FV3ATM-MOM6 coupling, the cplflx
163 ! could be .true., while cplice being .false..
164 if (cplice .and. cplflx) then
165 flag_cice(i) = .true.
166 else
167 flag_cice(i) = .false.
168 endif
169 islmsk(i) = 2
170 else
171 cice(i) = 0.0_kind_phys
172 flag_cice(i) = .false.
173 islmsk(i) = 0
174 icy(i) = .false.
175 endif
176 if (cice(i) < 1.0_kind_phys) then
177 wet(i) = .true. ! some open ocean
178 endif
179 else
180 if (cice(i) >= min_lakeice) then
181 icy(i) = .true.
182 islmsk(i) = 2
183 else
184 cice(i) = 0.0_kind_phys
185 islmsk(i) = 0
186 icy(i) = .false.
187 endif
188 flag_cice(i) = .false.
189 if (cice(i) < 1.0_kind_phys) then
190 wet(i) = .true. ! some open lake
191 endif
192 endif
193 endif
194 if (nint(slmsk(i)) /= 1) slmsk(i) = islmsk(i)
195 enddo
196
197 do i = 1, im
198 if (wet(i)) then
199 tsfc_wat(i) = t_surf(i)
200 end if
201 if (dry(i)) then
202 tsfcl(i) = t_surf(i)
203 end if
204 if (icy(i)) then
205 tisfc(i) = t_surf(i)
206 tisfc(i) = max(timin, min(tisfc(i), tgice))
207 end if
208 end do
209
210! to prepare to separate lake from ocean under water category
211 do i = 1, im
212 if ((wet(i) .or. icy(i)) .and. lakefrac(i) > 0.0_kind_phys) then
213 if (lkm == 1 .and. lakefrac(i) >= 0.15 .and. lakedepth(i) > 1.0_kind_phys) then
214 use_lake_model(i) = 1
215 else
216 use_lake_model(i) = 0
217 endif
218 else
219 use_lake_model(i) = 0
220 endif
221 enddo
222!
223
224 end subroutine scm_sfc_flux_spec_run
225
226end module scm_sfc_flux_spec