CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
rrtmg_sw_post.F90
1
4 contains
5
13 subroutine rrtmg_sw_post_run (im, levr, levs, ltp, nday, lm, kd, lsswr, &
14 swhtr, sfcalb1, sfcalb2, sfcalb3, sfcalb4, htswc, htsw0, &
15 nirbmdi, nirdfdi, visbmdi, visdfdi, nirbmui, nirdfui, visbmui,&
16 visdfui, sfcdsw, sfcnsw, htrsw, swhc, scmpsw, sfcfsw, topfsw, &
17 errmsg, errflg)
18
19 use machine, only: kind_phys
22
23 implicit none
24
25 integer, intent(in) :: im, levr, levs, &
26 ltp, nday, lm, kd
27 logical, intent(in) :: lsswr, swhtr
28 real(kind=kind_phys), dimension(:), intent(in) :: sfcalb1, sfcalb2, &
29 sfcalb3, sfcalb4
30 real(kind=kind_phys), dimension(:,:), intent(in) :: htswc, htsw0
31
32 real(kind=kind_phys), dimension(:), intent(inout) :: nirbmdi, nirdfdi, &
33 visbmdi, visdfdi, &
34 nirbmui, nirdfui, &
35 visbmui, visdfui, &
36 sfcdsw, sfcnsw
37 real(kind=kind_phys), dimension(:,:), intent(inout) :: htrsw, swhc
38
39 type(cmpfsw_type), dimension(:), intent(inout) :: scmpsw
40 type(sfcfsw_type), dimension(:), intent(inout) :: sfcfsw
41 type(topfsw_type), dimension(:), intent(inout) :: topfsw
42
43 character(len=*), intent(out) :: errmsg
44 integer, intent(out) :: errflg
45 ! Local variables
46 integer :: i, k1, k
47
48 ! Initialize CCPP error handling variables
49 errmsg = ''
50 errflg = 0
51
52 if (lsswr) then
53 if (nday > 0) then
54 do k = 1, lm
55 k1 = k + kd
56 htrsw(1:im,k) = htswc(1:im,k1)
57 enddo
58 ! We are assuming that radiative tendencies are from bottom to top
59 ! --- repopulate the points above levr i.e. LM
60 if (lm < levs) then
61 do k = lm+1, levs
62 htrsw(1:im,k) = htrsw(1:im,lm)
63 enddo
64 endif
65
66 if (swhtr) then
67 do k = 1, lm
68 k1 = k + kd
69 swhc(1:im,k) = htsw0(1:im,k1)
70 enddo
71 ! --- repopulate the points above levr i.e. LM
72 if (lm < levs) then
73 do k = lm+1, levs
74 swhc(1:im,k) = swhc(1:im,lm)
75 enddo
76 endif
77 endif
78
79! --- surface down and up spectral component fluxes
80! Save two spectral bands' surface downward and upward fluxes for
81! output.
82
83 do i=1,im
84 nirbmdi(i) = scmpsw(i)%nirbm
85 nirdfdi(i) = scmpsw(i)%nirdf
86 visbmdi(i) = scmpsw(i)%visbm
87 visdfdi(i) = scmpsw(i)%visdf
88
89 nirbmui(i) = scmpsw(i)%nirbm * sfcalb1(i)
90 nirdfui(i) = scmpsw(i)%nirdf * sfcalb2(i)
91 visbmui(i) = scmpsw(i)%visbm * sfcalb3(i)
92 visdfui(i) = scmpsw(i)%visdf * sfcalb4(i)
93 enddo
94
95 else ! if_nday_block
96
97 htrsw(:,:) = 0.0
98
99 sfcfsw = sfcfsw_type( 0.0, 0.0, 0.0, 0.0 )
100 topfsw = topfsw_type( 0.0, 0.0, 0.0 )
101 scmpsw = cmpfsw_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 )
102
103 do i=1,im
104 nirbmdi(i) = 0.0
105 nirdfdi(i) = 0.0
106 visbmdi(i) = 0.0
107 visdfdi(i) = 0.0
108
109 nirbmui(i) = 0.0
110 nirdfui(i) = 0.0
111 visbmui(i) = 0.0
112 visdfui(i) = 0.0
113 enddo
114
115 if (swhtr) then
116 swhc(:,:) = 0
117 endif
118
119 endif ! end_if_nday
120
121! --- radiation fluxes for other physics processes
122 do i=1,im
123 sfcnsw(i) = sfcfsw(i)%dnfxc - sfcfsw(i)%upfxc
124 sfcdsw(i) = sfcfsw(i)%dnfxc
125 enddo
126
127 endif ! end_if_lsswr
128
129 end subroutine rrtmg_sw_post_run
131 end module rrtmg_sw_post
subroutine rrtmg_sw_post_run(im, levr, levs, ltp, nday, lm, kd, lsswr, swhtr, sfcalb1, sfcalb2, sfcalb3, sfcalb4, htswc, htsw0, nirbmdi, nirdfdi, visbmdi, visdfdi, nirbmui, nirdfui, visbmui, visdfui, sfcdsw, sfcnsw, htrsw, swhc, scmpsw, sfcfsw, topfsw, errmsg, errflg)
This module is for specifying the band structures and program parameters used by the RRTMG-SW scheme.
Definition radsw_param.f:62
derived type for special components of surface SW fluxes
derived type for SW fluxes at surface
Definition radsw_param.f:89
derived type for SW fluxes at TOA
Definition radsw_param.f:79