BUGSrad
 All Classes Files Functions Variables
two_rt_lw.f
Go to the documentation of this file.
1 
2 
3 ! CVS: $Id: two_rt_lw.F,v 1.7 2003/11/11 21:55:13 norm Exp $
4 ! CVS: $Name: $
5 
6 !-----------------------------------------------------------------------
7 
8  subroutine two_rt_lw
9  + ( ncol , nlm, mbs , mbir
10  +, ib , wc, asym , tau
11  +, es , bf, fu , fd
12  +, sel_rules
13  + )
14 
15  use kinds
16 
17 
18 
19  implicit none
20 !-----------------------------------------------------------------------
21 ! REFERENCES:
22 ! two_rt_lw replaces two_rt and add written by G. Stephens. two_rt_lw
23 ! computes the spectral fluxes using a two-stream approximation method.
24 ! Philip Partain, Philip Gabriel, and Laura D. Fowler/graben (09-08-99).
25 
26 ! MODIFICATIONS:
27 ! * changed declarations to adapt the code from BUGS4 to BUGS5.
28 ! Laura D. Fowler/slikrock (02-01-00).
29 
30 ! SUBROUTINES CALLED:
31 ! none.
32 
33 ! FUNCTIONS CALLED:
34 ! none.
35 
36 ! INCLUDED COMMONS:
37 ! none.
38 
39 ! ARGUMENT LIST VARIABLES:
40 ! All arrays indexed as nlm correspond to variables defined in the
41 ! middle of layers. All arrays indexed as nlm+1 correspond to variables
42 ! defined at levels at the top and bottom of layers.
43 
44 ! INPUT ARGUMENTS:
45 ! ----------------
46  logical (kind=log_kind), intent(in)::
47  & sel_rules
48 
49  integer (kind=int_kind), intent(in)::
50  & ncol !Length of sub-domain.
51  &, nlm !Number of layers.
52  &, mbs !Number of SW spectral intervals.
53  &, mbir !Number of IR spectral intervals.
54  &, ib !Index of spectral interval.
55 
56  real (kind=dbl_kind), intent(in), dimension(ncol,mbir)::
57  & es !Spectral surface emissivity (-).
58 
59  real (kind=dbl_kind), intent(in), dimension(ncol,nlm)::
60  & wc !Single scattering albedo (-).
61  &, asym !Asymmetry factor (-).
62  &, tau !Optical depth (-).
63 
64  real (kind=dbl_kind), intent(in), dimension(ncol,nlm+1)::
65  & bf !Planck function (W/m^2).
66 
67 ! OUTPUT ARGUMENTS:
68 ! -----------------
69  real (kind=dbl_kind), intent(out), dimension(ncol,nlm+1)::
70  & fd !Spectral downward flux (W/m^2).
71  &, fu !Spectral upward flux (W/m^2).
72 
73 ! LOCAL VARIABLES:
74 
75  integer(kind=int_kind)
76  & i !Horizontal index.
77  &, l !Vertical index.
78  &, ibms !Index of spectral interval.
79 
80  real (kind=dbl_kind), dimension(nlm)::
81  & rr !
82  &, tr !
83  &, sigu !
84  &, sigd !
85 
86  real(kind=dbl_kind)
87  & aa , bb , beta0 , cc
88  &, diffac , denom , fact , eggtau
89  &, ggtau
90  &, kappa , oms , prop ,r, rinf
91  &, t , taus
92  data diffac /2./
93 
94  real (kind=dbl_kind), dimension(nlm)::
95  & td , vu , exptau
96 
97  real (kind=dbl_kind), dimension(nlm+1)::
98  & re, vd
99 
100 ! SELECTION RULE VARIABLES
101 
102  logical (kind=log_kind)::
103  & fail
104 
105  real (kind=dbl_kind)::
106  & tausthresh
107  &, wcthresh
108  &, tauscat
109 
110  data tausthresh / 0.001 /
111  data wcthresh / 0.975 /
112 
113 !----------------------------------------------------------------------
114 
115  fd(:,1) = 0.
116 
117  ibms = ib - mbs
118 
119  do 1000 i = 1, ncol
120 
121  if(sel_rules) then
122 
123  fail = .false.
124  tauscat = 0.0
125  do l = nlm, 1, -1
126  if (wc(i,l).gt.wcthresh) fail = .true.
127  tauscat = tauscat + wc(i,l)*tau(i,l)
128  enddo
129  if (fail.and.tauscat.ge.tausthresh) goto 2000
130 
132 ! print *,'selection rules'
133  do l = 1, nlm
134  exptau(l) = exp(-2*tau(i,l))
135  if(tau(i,l) .lt. .8e-2) then
136  sigu(l) = (bf(i,l)+bf(i,l+1))*tau(i,l)
137  sigd(l) = sigu(l)
138  else
139  prop = (1.-exptau(l))/tau(i,l)
140  aa = 2.-prop
141  bb = -2.*exptau(l)+prop
142  cc = 0.5
143  sigu(l) = (aa*bf(i,l)+bb*bf(i,l+1))*cc
144  sigd(l) = (bb*bf(i,l)+aa*bf(i,l+1))*cc
145  endif
146  fd(i,l+1) = sigd(l) + exptau(l) * fd(i,l)
147  enddo
148 
149  fu(i,nlm+1) = bf(i,nlm+1)*es(i,ibms)
150 ! & + fd(i,nlm+1)*(1.0-es(i,ibms))
151 
152  do l = nlm , 1, -1
153  fu(i,l) = sigu(l) + exptau(l) * fu(i,l+1)
154  enddo
155 
156  cycle
157 
158  endif
159 
161 
163 2000 re(1) = 0.
164  vd(1) = 0.
165 ! print *,'full up calculation'
166 
167  do l = 1, nlm
168  fact = asym(i,l)*asym(i,l)
169  oms = ((1.-fact)*wc(i,l))/(1.-fact*wc(i,l))
170  taus = (1.-fact*wc(i,l))*tau(i,l)
171 
172  beta0 = (4.+asym(i,l))/(8.*(1.+asym(i,l)))
173  t = diffac*(1.-oms*(1.-beta0)) !-0.25
174  r = diffac*oms*beta0 !-0.25
175  kappa = sqrt(t**2-r**2)
176  rinf = r/(kappa+t)
177  ggtau = kappa*taus
178  eggtau = exp(-ggtau)
179  denom = (1.-rinf**2*eggtau**2)
180  tr(l) = (1.-rinf**2)*eggtau/denom
181  rr(l) = rinf*(1.-eggtau**2)/denom
182 
183  if(taus .lt. .8e-2) then
184  sigu(l) = 0.5*diffac*(bf(i,l)+bf(i,l+1))*taus
185  sigd(l) = sigu(l)
186  else
187  aa = (t+r)*(1.-rr(l))-(1.+rr(l)-tr(l))/taus
188  bb = -(t+r)*tr(l)+(1.+rr(l)-tr(l))/taus
189  cc = diffac*(1.-oms)/kappa**2
190  sigu(l) = cc*(aa*bf(i,l)+bb*bf(i,l+1))
191  sigd(l) = cc*(bb*bf(i,l)+aa*bf(i,l+1))
192  endif
193  enddo
194 
195 !---- 1. do adding, going from top down:
196 
197  do l = 1, nlm
198  prop = 1. / (1. - re(l)*rr(l))
199  re(l+1) = rr(l) + tr(l)**2*re(l)*prop
200  vd(l+1) = sigd(l) + (tr(l)*vd(l)
201  + + tr(l)*re(l)*sigu(l))*prop
202  vu(l) = (rr(l)*vd(l) + sigu(l))*prop
203  td(l) = prop
204  enddo
205 
206 !---- 2. calculate fluxes going up through the layers:
207 
208  fu(i,nlm+1) = es(i,ibms)*bf(i,nlm+1)
209 
210  do l = nlm+1, 2, -1
211  fd(i,l) = re(l)*fu(i,l) + vd(l)
212  fu(i,l-1) = tr(l-1)*fu(i,l)*td(l-1) + vu(l-1)
213  enddo
214 
216 
217  1000 continue
218 
219  return
220  end subroutine two_rt_lw
221 
222 !------------------------------------------------------------------------
Generated on Tue Aug 13 2013 10:38:25 for BUGSrad by   doxygen 1.8.3.1