CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_surface_composites_inter.F90
1
3
5
6 use machine, only: kind_phys
7
8 implicit none
9
10 private
11
12 public gfs_surface_composites_inter_run
13
14contains
15
19 subroutine gfs_surface_composites_inter_run (im, dry, icy, wet, semis_wat, semis_lnd, semis_ice, &
20 adjsfcdlw, gabsbdlw_lnd, gabsbdlw_ice, gabsbdlw_wat,&
21 adjsfcusw, adjsfcdsw, adjsfcnsw, use_lake_model, errmsg, errflg)
22
23 implicit none
24
25 ! Interface variables
26 integer, intent(in ) :: im
27 logical, dimension(:), intent(in ) :: dry, icy
28 logical, dimension(:), intent(inout) :: wet
29 real(kind=kind_phys), dimension(:), intent(in ) :: semis_wat, semis_lnd, semis_ice, &
30 adjsfcdlw, adjsfcdsw, adjsfcnsw
31 real(kind=kind_phys), dimension(:), intent(inout) :: gabsbdlw_lnd, gabsbdlw_ice, gabsbdlw_wat
32 real(kind=kind_phys), dimension(:), intent(out) :: adjsfcusw
33 integer, dimension(:), intent(in) :: use_lake_model
34
35 ! CCPP error handling
36 character(len=*), intent(out) :: errmsg
37 integer, intent(out) :: errflg
38!
39 ! Local variables
40 integer :: i
41
42 ! Initialize CCPP error handling variables
43 errmsg = ''
44 errflg = 0
45
46 ! --- convert lw fluxes for land/ocean/sea-ice models - requires dcyc2t3 to set adjsfcdlw
47 ! note: for sw: adjsfcdsw and adjsfcnsw are zenith angle adjusted downward/net fluxes.
48 ! for lw: adjsfcdlw is (sfc temp adjusted) downward fluxe with no emiss effect.
49 ! adjsfculw is (sfc temp adjusted) upward fluxe including emiss effect.
50 ! one needs to be aware that that the absorbed downward lw flux (used by land/ocean
51 ! models as downward flux) is not the same as adjsfcdlw but a value reduced by
52 ! the factor of emissivity. however, the net effects are the same when seeing
53 ! it either above the surface interface or below.
54 !
55 ! - flux above the interface used by atmosphere model:
56 ! down: adjsfcdlw; up: adjsfculw = sfcemis*sigma*T**4 + (1-sfcemis)*adjsfcdlw
57 ! net = up - down = sfcemis * (sigma*T**4 - adjsfcdlw)
58 ! - flux below the interface used by lnd/oc/ice models:
59 ! down: sfcemis*adjsfcdlw; up: sfcemis*sigma*T**4
60 ! net = up - down = sfcemis * (sigma*T**4 - adjsfcdlw)
61 ! surface upwelling shortwave flux at current time is in adjsfcusw
62
63 ! --- ... define the downward lw flux absorbed by ground
64 do i=1,im
65 if (dry(i)) gabsbdlw_lnd(i) = semis_lnd(i) * adjsfcdlw(i)
66 if (icy(i)) gabsbdlw_ice(i) = semis_ice(i) * adjsfcdlw(i)
67 if (wet(i)) gabsbdlw_wat(i) = semis_wat(i) * adjsfcdlw(i)
68 adjsfcusw(i) = adjsfcdsw(i) - adjsfcnsw(i)
69 enddo
70
71 end subroutine gfs_surface_composites_inter_run
72