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)
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(:)
71 character(len=*),
intent(out) :: errmsg
72 integer,
intent(out) :: errflg
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
79 real(kind=kind_phys),
parameter :: timin = 173.0_kind_phys
88 sh_flux(i) = spec_sh_flux(i)
89 lh_flux(i) = spec_lh_flux(i)
90 sh_flux_chs(i) = sh_flux(i)
92 roughness_length_m = 0.01*roughness_length(i)
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))
98 sfc_stress(i) = u_star(i)*u_star(i)
100 if(wind1(i) > 0.0)
then
101 cm(i) = sfc_stress(i)/(wind1(i)*wind1(i))
105 fm(i) = sqrt((vonkarman*vonkarman)/cm(i))
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)
113 thv1 = t1(i) * exner_inverse(i) * (1.0 + fvirt*q1_non_neg)
116 w_thv1 = sh_flux(i)*exner_inverse(i) + fvirt*t1(i)*lh_flux(i)
118 obukhov_length = -u_star(i)*u_star(i)*u_star(i)*thv1/(vonkarman*grav*w_thv1)
121 tvs = t_surf(i)*(1.0 + fvirt*q1_non_neg)
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)))
128 fh(i) = rb(i)*fm(i)*fm(i)*obukhov_length/z1(i)
129 ch(i) = vonkarman*vonkarman/(fm(i)*fh(i))
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)
147 if (islmsk(i) == 1)
then
149 frland(i) = 1.0_kind_phys
150 cice(i) = 0.0_kind_phys
153 frland(i) = 0.0_kind_phys
154 if (oceanfrac(i) > 0.0_kind_phys)
then
155 if (cice(i) >= min_seaice)
then
164 if (cplice .and. cplflx)
then
165 flag_cice(i) = .true.
167 flag_cice(i) = .false.
171 cice(i) = 0.0_kind_phys
172 flag_cice(i) = .false.
176 if (cice(i) < 1.0_kind_phys)
then
180 if (cice(i) >= min_lakeice)
then
184 cice(i) = 0.0_kind_phys
188 flag_cice(i) = .false.
189 if (cice(i) < 1.0_kind_phys)
then
194 if (nint(slmsk(i)) /= 1) slmsk(i) = islmsk(i)
199 tsfc_wat(i) = t_surf(i)
206 tisfc(i) = max(timin, min(tisfc(i), tgice))
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
216 use_lake_model(i) = 0
219 use_lake_model(i) = 0