da_effht_tl.inc
References to this file elsewhere.
1 subroutine da_effht_tl(ho,hv,sigo,sigv,mu,zcld,hdn,hup,hdninf,hupinf, &
2 TGL_ho,TGL_hv,TGL_sigo,TGL_sigv,TGL_mu, &
3 TGL_zcld,TGL_hdn,TGL_hup,TGL_hdninf,TGL_hupinf)
4
5 !--------------------------------------------------------------------
6 ! Purpose: TBD
7 ! Input : TGL_ho, TGL_hv, TGL_sigo, TGL_sigv, TGL_mu, TGL_zcld
8 ! Output : TGL_hdn, hdn, TGL_hup, hup,
9 ! TGL_hdninf, hdninf, TGL_hupinf, hupinf
10 !--------------------------------------------------------------------
11
12 implicit none
13
14 real, intent(in) :: ho,hv,sigo,sigv,mu,zcld
15 real, intent(in) :: TGL_ho,TGL_hv,TGL_sigo,TGL_sigv,TGL_zcld, TGL_mu
16 real, intent(out) :: hdn,hup,hdninf,hupinf
17 real, intent(out) :: TGL_hdn,TGL_hup,TGL_hdninf,TGL_hupinf
18
19 real :: gint,zgint,hint,zhint
20 real :: ginf,zginf,hinf,zhinf
21 real :: TGL_gint,TGL_zgint,TGL_hint,TGL_zhint
22 real :: TGL_ginf,TGL_zginf,TGL_hinf,TGL_zhinf
23 real :: TGL_mu2,TGL_halfmu,TGL_sixthmu2,TGL_etnthmu2
24 real :: TGL_quartmu,TGL_halfmu2
25
26 real :: hoinv,hvinv,chio,chiv,ezho,ezhv,alpha,alph2,alph3
27 real :: beta,beta2,beta3,mu2,mualph,cplus,cmin,dplus,dmin
28 real :: chiov,chivv,chioo,chioov,chiovv,chiooo,chivvv
29 real :: h11,h21,h12,newh11
30 real :: sigoo,sigov,sigvv,sigooo,sigoov,sigovv,sigvvv
31 real :: ezhoo,ezhov,ezhvv,ezhooo,ezhoov,ezhovv,ezhvvv
32 real :: s,sprim,t,tprim,u,uprim,term1,term2,term3
33 real :: halfmu,halfmu2,sixthmu2,etnthmu2,quartmu
34
35 real :: TGL_hoinv,TGL_hvinv,TGL_chio,TGL_chiv,TGL_ezho
36 real :: TGL_ezhv,TGL_alpha,TGL_alph2,TGL_alph3
37 real :: TGL_beta,TGL_beta2,TGL_beta3,TGL_mualph
38 real :: TGL_cplus,TGL_cmin,TGL_dplus,TGL_dmin
39 real :: TGL_chiov,TGL_chivv,TGL_chioo,TGL_chioov
40 real :: TGL_chiovv,TGL_chiooo,TGL_chivvv
41 real :: TGL_h11,TGL_h21,TGL_h12,TGL_newh11
42 real :: TGL_sigoo,TGL_sigov,TGL_sigvv,TGL_sigooo
43 real :: TGL_sigoov,TGL_sigovv,TGL_sigvvv
44 real :: TGL_ezhoo,TGL_ezhov,TGL_ezhvv,TGL_ezhooo
45 real :: TGL_ezhoov,TGL_ezhovv,TGL_ezhvvv
46 real :: TGL_s,TGL_sprim,TGL_t,TGL_tprim
47 real :: TGL_u,TGL_uprim,TGL_term1,TGL_term2,TGL_term3
48
49 if (trace_use) call da_trace_entry("da_effht_tl")
50
51 hoinv = 1.0d0/ho
52 TGL_hoinv = -1.0d0*hoinv*hoinv*TGL_ho
53
54 hvinv = 1.0d0/hv
55 TGL_hvinv = -1.0d0*hvinv*hvinv*TGL_hv
56
57 chio = zcld*hoinv
58 TGL_chio = TGL_zcld*hoinv + zcld*TGL_hoinv
59
60 chiv = zcld*hvinv
61 TGL_chiv = TGL_zcld*hvinv + zcld*TGL_hvinv
62
63 ezho = sigo*exp(-chio)
64 TGL_ezho = TGL_sigo*exp(-chio)-TGL_chio*ezho
65
66 ezhv = sigv*exp(-chiv)
67 TGL_ezhv = TGL_sigv*exp(-chiv)-TGL_chiv*ezhv
68
69 alpha = sigo + sigv
70 TGL_alpha = TGL_sigo + TGL_sigv
71
72 alph2 = alpha*alpha
73 TGL_alph2 = 2.0*alpha*TGL_alpha
74
75 alph3 = alpha*alph2
76 TGL_alph3 = TGL_alpha*alph2+alpha*TGL_alph2
77
78 beta = ezho + ezhv
79 TGL_beta = TGL_ezho + TGL_ezhv
80
81 beta2 = beta*beta
82 TGL_beta2 = 2.0*beta*TGL_beta
83
84 beta3 = beta*beta2
85 TGL_beta3 = TGL_beta*beta2+beta*TGL_beta2
86
87 mu2 = mu*mu
88 TGL_mu2 = 2.0*mu*TGL_mu
89 halfmu = 0.5d0* mu
90 TGL_halfmu = 0.5d0*TGL_mu
91 sixthmu2 = mu2/6.0d0
92 TGL_sixthmu2 = TGL_mu2/6.0d0
93 etnthmu2 = mu2/18.0d0
94 TGL_etnthmu2 = TGL_mu2/18.0d0
95 quartmu = 0.25d0* mu
96 TGL_quartmu = 0.25d0*TGL_mu
97 halfmu2 = 0.5d0* mu2
98 TGL_halfmu2 = 0.5d0*TGL_mu2
99
100 mualph = mu*alpha
101 TGL_mualph = TGL_mu*alpha + mu*TGL_alpha
102
103 cplus = 1.0d0 + mualph
104 TGL_cplus = TGL_mualph
105
106 cmin = 1.0d0 - mualph
107 TGL_cmin = - TGL_mualph
108
109 dplus = halfmu2*alph2
110 TGL_dplus = TGL_halfmu2*alph2 + halfmu2*TGL_alph2
111
112 dmin = dplus
113 TGL_dmin = TGL_dplus
114
115 TGL_dplus = TGL_cplus + TGL_dplus
116 dplus = cplus + dplus
117
118 TGL_dmin = TGL_cmin + TGL_dmin
119 dmin = cmin + dmin
120
121
122 h11 = hoinv + hvinv
123 TGL_h11 = TGL_hoinv + TGL_hvinv
124
125 h21 = 1.0d0/(h11 + hvinv)
126 TGL_h21 = -1.0d0*h21*h21*(TGL_h11+TGL_hvinv)
127
128 h12 = 1.0d0/(h11 + hoinv)
129 TGL_h12 = -1.0d0*h12*h12*(TGL_h11 + TGL_hoinv)
130
131 newh11 = 1.0d0/h11
132 TGL_newh11 = -1.0d0*newh11*newh11*TGL_h11
133
134 chiov = 1.0d0 + chio + chiv
135 TGL_chiov = TGL_chio + TGL_chiv
136
137 chioo = 1.0d0 + chio + chio
138 TGL_chioo = TGL_chio + TGL_chio
139
140 chivv = 1.0d0 + chiv + chiv
141 TGL_chivv = TGL_chiv + TGL_chiv
142
143 chioov = chioo + chiv
144 TGL_chioov = TGL_chioo + TGL_chiv
145
146 chiovv = chio + chivv
147 TGL_chiovv = TGL_chio + TGL_chivv
148
149 chiooo = chioo + chio
150 TGL_chiooo = TGL_chioo + TGL_chio
151
152 chivvv = chivv + chiv
153 TGL_chivvv = TGL_chivv + TGL_chiv
154
155 TGL_chio = TGL_chio
156 chio = 1.0d0 + chio
157
158 TGL_chiv = TGL_chiv
159 chiv = 1.0d0 + chiv
160
161 sigov = sigo*sigv
162 TGL_sigov = TGL_sigo*sigv + sigo*TGL_sigv
163
164 sigoo = sigo*sigo
165 TGL_sigoo = 2.0*sigo*TGL_sigo
166
167 sigvv = sigv*sigv
168 TGL_sigvv = 2.0*sigv*TGL_sigv
169
170 sigooo = sigoo*sigo
171 TGL_sigooo = TGL_sigoo*sigo + sigoo*TGL_sigo
172
173 sigoov = sigoo*sigv
174 TGL_sigoov = TGL_sigoo*sigv + sigoo*TGL_sigv
175
176 sigovv = sigo*sigvv
177 TGL_sigovv = TGL_sigo*sigvv + sigo*TGL_sigvv
178
179 sigvvv = sigvv*sigv
180 TGL_sigvvv = TGL_sigvv*sigv + sigvv*TGL_sigv
181
182 ezhoo = ezho*ezho
183 TGL_ezhoo = 2.0*ezho*TGL_ezho
184
185 ezhov = ezho*ezhv
186 TGL_ezhov = TGL_ezho*ezhv + ezho*TGL_ezhv
187
188 ezhvv = ezhv*ezhv
189 TGL_ezhvv = 2.0*ezhv*TGL_ezhv
190
191 ezhovv = ezho*ezhvv
192 TGL_ezhovv = TGL_ezho*ezhvv + ezho*TGL_ezhvv
193
194 ezhoov = ezhoo*ezhv
195 TGL_ezhoov = TGL_ezhoo*ezhv + ezhoo*TGL_ezhv
196
197 ezhooo = ezhoo*ezho
198 TGL_ezhooo = TGL_ezhoo*ezho + ezhoo*TGL_ezho
199
200 ezhvvv = ezhvv*ezhv
201 TGL_ezhvvv = TGL_ezhvv*ezhv + ezhvv*TGL_ezhv
202
203 s = sigo*ho + sigv*hv
204 TGL_s = TGL_sigo*ho + sigo*TGL_ho + TGL_sigv*hv + sigv*TGL_hv
205
206 sprim = ezho*ho*chio + ezhv*hv*chiv
207 TGL_sprim = TGL_ezho*ho*chio + ezho*TGL_ho*chio + ezho*ho*TGL_chio + &
208 TGL_ezhv*hv*chiv + ezhv*TGL_hv*chiv + ezhv*hv*TGL_chiv
209
210 t = sigoo*ho + 4.0d0*sigov*newh11 + sigvv*hv
211 TGL_t = TGL_sigoo*ho + sigoo*TGL_ho + &
212 4.0d0*(TGL_sigov*newh11 + sigov*TGL_newh11) + &
213 TGL_sigvv*hv + sigvv*TGL_hv
214
215 tprim = ezhoo*ho*chioo + 4.0d0*ezhov*newh11*chiov + ezhvv*hv*chivv
216 TGL_tprim = TGL_ezhoo*ho*chioo +ezhoo*TGL_ho*chioo + ezhoo*ho*TGL_chioo + &
217 4.0d0*(TGL_ezhov*newh11*chiov + &
218 ezhov*TGL_newh11*chiov + &
219 ezhov*newh11*TGL_chiov ) + &
220 TGL_ezhvv*hv*chivv + ezhvv*TGL_hv*chivv + ezhvv*hv*TGL_chivv
221
222 u = sigooo*ho + 9.0d0*(sigovv*h21+sigoov*h12) + sigvvv*hv
223 TGL_u = TGL_sigooo*ho + sigooo*TGL_ho + &
224 9.0d0*(TGL_sigovv*h21 + sigovv*TGL_h21 + &
225 TGL_sigoov*h12 + sigoov*TGL_h12 ) + &
226 TGL_sigvvv*hv + sigvvv*TGL_hv
227
228 uprim = ezhvvv*hv*chivvv + &
229 9.0d0*(ezhovv*h21*chiovv + ezhoov*h12*chioov) + &
230 ezhooo*ho*chiooo
231 TGL_uprim = TGL_ezhvvv*hv*chivvv +ezhvvv*TGL_hv*chivvv+ ezhvvv*hv*TGL_chivvv+ &
232 9.0d0*(TGL_ezhovv*h21*chiovv + &
233 ezhovv*TGL_h21*chiovv + &
234 ezhovv*h21*TGL_chiovv + &
235 TGL_ezhoov*h12*chioov + &
236 ezhoov*TGL_h12*chioov + &
237 ezhoov*h12*TGL_chioov ) + &
238 TGL_ezhooo*ho*chiooo + ezhooo*TGL_ho*chiooo + ezhooo*ho*TGL_chiooo
239
240 term1 = s - sprim
241 TGL_term1 = TGL_s - TGL_sprim
242
243 term2 = quartmu*(t - tprim)
244 TGL_term2 = TGL_quartmu*(t - tprim) + quartmu*(TGL_t - TGL_tprim)
245
246 term3 = etnthmu2*( u - uprim)
247 TGL_term3 = TGL_etnthmu2*(u - uprim) + etnthmu2*(TGL_u - TGL_uprim)
248
249 zgint = dmin*term1 + cmin*term2 + term3
250 TGL_zgint = TGL_dmin*term1 + dmin*TGL_term1 + &
251 TGL_cmin*term2 + cmin*TGL_term2 + TGL_term3
252
253 zhint = -dplus*term1 + cplus*term2 - term3
254 TGL_zhint = -TGL_dplus*term1 - dplus*TGL_term1 + &
255 TGL_cplus*term2 + cplus*TGL_term2 - TGL_term3
256
257 term2 = quartmu * t
258 TGL_term2 = TGL_quartmu*t + quartmu*TGL_t
259
260 term3 = etnthmu2*u
261 TGL_term3 = TGL_etnthmu2*u + etnthmu2*TGL_u
262
263 zginf = dmin*s + cmin*term2 + term3
264 TGL_zginf = TGL_dmin*s + dmin*TGL_s + &
265 TGL_cmin*term2 + cmin*TGL_term2 + TGL_term3
266
267 zhinf = -dplus*s + cplus*term2 - term3
268 TGL_zhinf = -TGL_dplus*s - dplus*TGL_s + &
269 TGL_cplus*term2 + cplus*TGL_term2 - TGL_term3
270
271 term1 = alpha - beta
272 TGL_term1 = TGL_alpha - TGL_beta
273
274 term2 = halfmu*( alph2 - beta2)
275 TGL_term2 = TGL_halfmu*(alph2 - beta2) + halfmu*(TGL_alph2 - TGL_beta2)
276
277 term3 = sixthmu2*( alph3 - beta3)
278 TGL_term3 = TGL_sixthmu2*(alph3 - beta3) + sixthmu2*(TGL_alph3 - TGL_beta3)
279
280 gint = dmin*term1 + cmin*term2 + term3
281 TGL_gint = TGL_dmin*term1 + dmin*TGL_term1 + &
282 TGL_cmin*term2 + cmin*TGL_term2 + TGL_term3
283
284 hint = -dplus*term1 + cplus*term2 - term3
285 TGL_hint = -TGL_dplus*term1 - dplus*TGL_term1 + &
286 TGL_cplus*term2 + cplus*TGL_term2 - TGL_term3
287
288 term2 = halfmu*alph2
289 TGL_term2 = TGL_halfmu*alph2 + halfmu*TGL_alph2
290
291 term3 = sixthmu2*alph3
292 TGL_term3 = TGL_sixthmu2*alph3 + sixthmu2*TGL_alph3
293
294 ginf = dmin*alpha + cmin*term2 + term3
295 TGL_ginf = TGL_dmin*alpha + dmin*TGL_alpha + &
296 TGL_cmin*term2 + cmin*TGL_term2 + TGL_term3
297
298 hinf = -dplus*alpha + cplus*term2 - term3
299 TGL_hinf = -TGL_dplus*alpha - dplus*TGL_alpha + &
300 TGL_cplus*term2 + cplus*TGL_term2 - TGL_term3
301
302 hdn = zgint/gint
303 TGL_hdn = TGL_zgint/gint - hdn * TGL_gint/gint
304
305 hup = zhint/hint
306 TGL_hup = TGL_zhint/hint - hup*TGL_hint/hint
307
308 hdninf = zginf/ginf
309 TGL_hdninf = TGL_zginf/ginf - hdninf*TGL_ginf/ginf
310
311 hupinf = zhinf/hinf
312 TGL_hupinf = TGL_zhinf/hinf - hupinf*TGL_hinf/hinf
313
314 if (trace_use) call da_trace_exit("da_effht_tl")
315
316 end subroutine da_effht_tl
317
318