1 | SUBROUTINE SFLUXI(PLEV,TLEV,DTAUI,TAUCUMI,UBARI,RSFI,WNOI,DWNI, |
---|
2 | * COSBI,WBARI,GWEIGHT,NFLUXTOPI,NFLUXTOPI_nu, |
---|
3 | * FMNETI,fluxupi,fluxdni,fluxupi_nu, |
---|
4 | * FZEROI,TAUGSURF) |
---|
5 | |
---|
6 | use radinc_h |
---|
7 | use radcommon_h, only: planckir, tlimit,sigma |
---|
8 | |
---|
9 | implicit none |
---|
10 | |
---|
11 | #include "comcstfi.h" |
---|
12 | |
---|
13 | integer NLEVRAD, L, NW, NG, NTS, NTT |
---|
14 | |
---|
15 | real*8 TLEV(L_LEVELS), PLEV(L_LEVELS) |
---|
16 | real*8 TAUCUMI(L_LEVELS,L_NSPECTI,L_NGAUSS) |
---|
17 | real*8 FMNETI(L_NLAYRAD) |
---|
18 | real*8 WNOI(L_NSPECTI), DWNI(L_NSPECTI) |
---|
19 | real*8 DTAUI(L_NLAYRAD,L_NSPECTI,L_NGAUSS) |
---|
20 | real*8 FMUPI(L_NLEVRAD), FMDI(L_NLEVRAD) |
---|
21 | real*8 COSBI(L_NLAYRAD,L_NSPECTI,L_NGAUSS) |
---|
22 | real*8 WBARI(L_NLAYRAD,L_NSPECTI,L_NGAUSS) |
---|
23 | real*8 GWEIGHT(L_NGAUSS), NFLUXTOPI |
---|
24 | real*8 NFLUXTOPI_nu(L_NSPECTI) |
---|
25 | real*8 fluxupi_nu(L_NLAYRAD,L_NSPECTI) |
---|
26 | real*8 FTOPUP |
---|
27 | |
---|
28 | real*8 UBARI, RSFI, TSURF, BSURF, TTOP, BTOP, TAUTOP |
---|
29 | real*8 PLANCK, PLTOP |
---|
30 | real*8 fluxupi(L_NLAYRAD), fluxdni(L_NLAYRAD) |
---|
31 | real*8 FZEROI(L_NSPECTI) |
---|
32 | real*8 taugsurf(L_NSPECTI,L_NGAUSS-1), fzero |
---|
33 | |
---|
34 | real*8 fup_tmp(L_NSPECTI),fdn_tmp(L_NSPECTI) |
---|
35 | real*8 PLANCKSUM,PLANCKREF |
---|
36 | |
---|
37 | |
---|
38 | C======================================================================C |
---|
39 | |
---|
40 | NLEVRAD = L_NLEVRAD |
---|
41 | |
---|
42 | |
---|
43 | C ZERO THE NET FLUXES |
---|
44 | |
---|
45 | NFLUXTOPI = 0.0D0 |
---|
46 | |
---|
47 | DO NW=1,L_NSPECTI |
---|
48 | NFLUXTOPI_nu(NW) = 0.0D0 |
---|
49 | DO L=1,L_NLAYRAD |
---|
50 | FLUXUPI_nu(L,NW) = 0.0D0 |
---|
51 | |
---|
52 | fup_tmp(nw)=0.0D0 |
---|
53 | fdn_tmp(nw)=0.0D0 |
---|
54 | |
---|
55 | END DO |
---|
56 | END DO |
---|
57 | |
---|
58 | DO L=1,L_NLAYRAD |
---|
59 | FMNETI(L) = 0.0D0 |
---|
60 | FLUXUPI(L) = 0.0D0 |
---|
61 | FLUXDNI(L) = 0.0D0 |
---|
62 | END DO |
---|
63 | |
---|
64 | C WE NOW ENTER A MAJOR LOOP OVER SPECTRAL INTERVALS IN THE INFRARED |
---|
65 | C TO CALCULATE THE NET FLUX IN EACH SPECTRAL INTERVAL |
---|
66 | |
---|
67 | TTOP = TLEV(2) ! JL12 why not (1) ??? |
---|
68 | TSURF = TLEV(L_LEVELS) |
---|
69 | |
---|
70 | NTS = int(TSURF*NTfac)-NTstar+1 |
---|
71 | NTT = int(TTOP *NTfac)-NTstar+1 |
---|
72 | |
---|
73 | !JL12 corrects the surface planck function so that its integral is equal to sigma Tsurf^4 |
---|
74 | !JL12 this ensure that no flux is lost due to: |
---|
75 | !JL12 -truncation of the planck function at high/low wavenumber |
---|
76 | !JL12 -numerical error during first spectral integration |
---|
77 | !JL12 -discrepancy between Tsurf and NTS/NTfac |
---|
78 | PLANCKSUM=0.d0 |
---|
79 | PLANCKREF=TSURF*TSURF |
---|
80 | PLANCKREF=sigma*PLANCKREF*PLANCKREF |
---|
81 | DO NW=1,L_NSPECTI |
---|
82 | PLANCKSUM=PLANCKSUM+PLANCKIR(NW,NTS)*DWNI(NW) |
---|
83 | ENDDO |
---|
84 | PLANCKSUM=PLANCKREF/(PLANCKSUM*Pi) |
---|
85 | !JL12 |
---|
86 | |
---|
87 | DO 501 NW=1,L_NSPECTI |
---|
88 | |
---|
89 | C SURFACE EMISSIONS - INDEPENDENT OF GAUSS POINTS |
---|
90 | BSURF = (1.-RSFI)*PLANCKIR(NW,NTS)*PLANCKSUM !JL12 plancksum see above |
---|
91 | PLTOP = PLANCKIR(NW,NTT) |
---|
92 | |
---|
93 | C If FZEROI(NW) = 1, then the k-coefficients are zero - skip to the |
---|
94 | C special Gauss point at the end. |
---|
95 | |
---|
96 | FZERO = FZEROI(NW) |
---|
97 | IF(FZERO.ge.0.99) goto 40 |
---|
98 | |
---|
99 | DO NG=1,L_NGAUSS-1 |
---|
100 | |
---|
101 | if(TAUGSURF(NW,NG).lt. TLIMIT) then |
---|
102 | fzero = fzero + (1.0D0-FZEROI(NW))*GWEIGHT(NG) |
---|
103 | goto 30 |
---|
104 | end if |
---|
105 | |
---|
106 | C SET UP THE UPPER AND LOWER BOUNDARY CONDITIONS ON THE IR |
---|
107 | C CALCULATE THE DOWNWELLING RADIATION AT THE TOP OF THE MODEL |
---|
108 | C OR THE TOP LAYER WILL COOL TO SPACE UNPHYSICALLY |
---|
109 | |
---|
110 | ! TAUTOP = DTAUI(1,NW,NG)*PLEV(2)/(PLEV(4)-PLEV(2)) |
---|
111 | TAUTOP = TAUCUMI(2,NW,NG) |
---|
112 | BTOP = (1.0D0-EXP(-TAUTOP/UBARI))*PLTOP |
---|
113 | |
---|
114 | C WE CAN NOW SOLVE FOR THE COEFFICIENTS OF THE TWO STREAM |
---|
115 | C CALL A SUBROUTINE THAT SOLVES FOR THE FLUX TERMS |
---|
116 | C WITHIN EACH INTERVAL AT THE MIDPOINT WAVENUMBER |
---|
117 | |
---|
118 | CALL GFLUXI(NLEVRAD,TLEV,NW,DWNI(NW),DTAUI(1,NW,NG), |
---|
119 | * TAUCUMI(1,NW,NG), |
---|
120 | * WBARI(1,NW,NG),COSBI(1,NW,NG),UBARI,RSFI,BTOP, |
---|
121 | * BSURF,FTOPUP,FMUPI,FMDI) |
---|
122 | |
---|
123 | |
---|
124 | |
---|
125 | C NOW CALCULATE THE CUMULATIVE IR NET FLUX |
---|
126 | |
---|
127 | NFLUXTOPI = NFLUXTOPI+FTOPUP*DWNI(NW)*GWEIGHT(NG)* |
---|
128 | * (1.0D0-FZEROI(NW)) |
---|
129 | |
---|
130 | c and same thing by spectral band... (RDW) |
---|
131 | NFLUXTOPI_nu(NW) = NFLUXTOPI_nu(NW) |
---|
132 | * +FTOPUP*DWNI(NW)*GWEIGHT(NG)*(1.0D0-FZEROI(NW)) |
---|
133 | |
---|
134 | |
---|
135 | DO L=1,L_NLEVRAD-1 |
---|
136 | |
---|
137 | C CORRECT FOR THE WAVENUMBER INTERVALS |
---|
138 | |
---|
139 | FMNETI(L) = FMNETI(L)+(FMUPI(L)-FMDI(L))*DWNI(NW)* |
---|
140 | * GWEIGHT(NG)*(1.0D0-FZEROI(NW)) |
---|
141 | FLUXUPI(L) = FLUXUPI(L) + FMUPI(L)*DWNI(NW)*GWEIGHT(NG)* |
---|
142 | * (1.0D0-FZEROI(NW)) |
---|
143 | FLUXDNI(L) = FLUXDNI(L) + FMDI(L)*DWNI(NW)*GWEIGHT(NG)* |
---|
144 | * (1.0D0-FZEROI(NW)) |
---|
145 | |
---|
146 | c and same thing by spectral band... (RW) |
---|
147 | FLUXUPI_nu(L,NW) = FLUXUPI_nu(L,NW) + |
---|
148 | * FMUPI(L)*DWNI(NW)*GWEIGHT(NG)*(1.0D0-FZEROI(NW)) |
---|
149 | |
---|
150 | END DO |
---|
151 | |
---|
152 | 30 CONTINUE |
---|
153 | |
---|
154 | END DO !End NGAUSS LOOP |
---|
155 | |
---|
156 | 40 CONTINUE |
---|
157 | |
---|
158 | C SPECIAL 17th Gauss point |
---|
159 | |
---|
160 | NG = L_NGAUSS |
---|
161 | |
---|
162 | ! TAUTOP = DTAUI(1,NW,NG)*PLEV(2)/(PLEV(4)-PLEV(2)) |
---|
163 | TAUTOP = TAUCUMI(2,NW,NG) |
---|
164 | BTOP = (1.0D0-EXP(-TAUTOP/UBARI))*PLTOP |
---|
165 | |
---|
166 | C WE CAN NOW SOLVE FOR THE COEFFICIENTS OF THE TWO STREAM |
---|
167 | C CALL A SUBROUTINE THAT SOLVES FOR THE FLUX TERMS |
---|
168 | C WITHIN EACH INTERVAL AT THE MIDPOINT WAVENUMBER |
---|
169 | |
---|
170 | |
---|
171 | CALL GFLUXI(NLEVRAD,TLEV,NW,DWNI(NW),DTAUI(1,NW,NG), |
---|
172 | * TAUCUMI(1,NW,NG), |
---|
173 | * WBARI(1,NW,NG),COSBI(1,NW,NG),UBARI,RSFI,BTOP, |
---|
174 | * BSURF,FTOPUP,FMUPI,FMDI) |
---|
175 | |
---|
176 | C NOW CALCULATE THE CUMULATIVE IR NET FLUX |
---|
177 | |
---|
178 | NFLUXTOPI = NFLUXTOPI+FTOPUP*DWNI(NW)*FZERO |
---|
179 | |
---|
180 | c and same thing by spectral band... (RW) |
---|
181 | NFLUXTOPI_nu(NW) = NFLUXTOPI_nu(NW) |
---|
182 | * +FTOPUP*DWNI(NW)*FZERO |
---|
183 | |
---|
184 | DO L=1,L_NLEVRAD-1 |
---|
185 | |
---|
186 | C CORRECT FOR THE WAVENUMBER INTERVALS |
---|
187 | |
---|
188 | FMNETI(L) = FMNETI(L)+(FMUPI(L)-FMDI(L))*DWNI(NW)*FZERO |
---|
189 | FLUXUPI(L) = FLUXUPI(L) + FMUPI(L)*DWNI(NW)*FZERO |
---|
190 | FLUXDNI(L) = FLUXDNI(L) + FMDI(L)*DWNI(NW)*FZERO |
---|
191 | |
---|
192 | c and same thing by spectral band... (RW) |
---|
193 | FLUXUPI_nu(L,NW) = FLUXUPI_nu(L,NW) + FMUPI(L)*DWNI(NW)*FZERO |
---|
194 | |
---|
195 | END DO |
---|
196 | |
---|
197 | 501 CONTINUE !End Spectral Interval LOOP |
---|
198 | |
---|
199 | C *** END OF MAJOR SPECTRAL INTERVAL LOOP IN THE INFRARED**** |
---|
200 | |
---|
201 | RETURN |
---|
202 | END |
---|