da_find_layer_adj.inc
References to this file elsewhere.
1 subroutine da_find_layer_adj(layer,tv,pre,pre_ma,tv_ma,ks,ke, &
2 ADJ_tv,ADJ_pre_ma,ADJ_tv_ma)
3
4 !-----------------------------------------------------------------------
5 ! Purpose: adjoint routine for da_find_layer
6 !-----------------------------------------------------------------------
7
8 implicit none
9
10 integer, intent(in) :: ks, ke
11 integer, intent(out) :: layer
12 real, dimension(ks-1:ke+1), intent(in) :: pre_ma, tv_ma
13 real, dimension(ks-1:ke+1), intent(inout) :: ADJ_pre_ma, ADJ_tv_ma
14 real, intent(in) :: pre, tv
15 real, intent(inout) :: ADJ_tv
16
17 integer :: k
18 real :: alpha, coef1, coef2
19 real :: ADJ_alpha
20
21 if (pre >= pre_ma(ks)) then
22
23 layer = ks
24 coef1=log(pre/pre_ma(ks+1))/(pre_ma(ks)* &
25 (log(pre_ma(ks)/pre_ma(ks+1)))**2)
26 coef2=log(pre_ma(ks)/pre)/(pre_ma(ks+1)* &
27 (log(pre_ma(ks)/pre_ma(ks+1)))**2)
28 alpha = log(pre_ma(ks)/pre)/log(pre_ma(ks)/pre_ma(ks+1))
29
30 ADJ_pre_ma(ks-1)= 0.
31 ADJ_tv_ma(ks) = ADJ_tv_ma(ks) + (1.-alpha)*ADJ_tv
32 ADJ_alpha = (tv_ma(ks+1)-tv_ma(ks))*ADJ_tv
33 ADJ_tv_ma(ks+1) = ADJ_tv_ma(ks+1) + alpha*ADJ_tv
34
35 ADJ_pre_ma(ks) = ADJ_pre_ma(ks) + coef1 * ADJ_alpha
36 ADJ_pre_ma(ks+1) = ADJ_pre_ma(ks+1) + coef2 * ADJ_alpha
37
38 else if (pre <= pre_ma(ke)) then
39
40 layer = ke+1
41 coef1=log(pre/pre_ma(ke))/(pre_ma(ke-1)* &
42 (log(pre_ma(ke-1)/pre_ma(ke)))**2)
43 coef2=log(pre_ma(ke-1)/pre)/(pre_ma(ke)* &
44 (log(pre_ma(ke-1)/pre_ma(ke)))**2)
45 alpha = log(pre_ma(ke-1)/pre)/log(pre_ma(ke-1)/pre_ma(ke))
46
47 ADJ_pre_ma(ke+1) = 0.
48 ADJ_tv_ma(ke-1) = ADJ_tv_ma(ke-1) + (1.-alpha)*ADJ_tv
49 ADJ_alpha = (tv_ma(ke)-tv_ma(ke-1))*ADJ_tv
50 ADJ_tv_ma(ke) = ADJ_tv_ma(ke) + alpha*ADJ_tv
51
52 ADJ_pre_ma(ke-1) = ADJ_pre_ma(ke-1) + coef1 * ADJ_alpha
53 ADJ_pre_ma(ke) = ADJ_pre_ma(ke) + coef2 * ADJ_alpha
54
55 else
56
57 do k=ks,ke-1
58 if (pre>=pre_ma(k+1) .and. pre<pre_ma(k)) then
59 layer = k+1
60 coef1=log(pre/pre_ma(k+1))/(pre_ma(k)* &
61 (log(pre_ma(k)/pre_ma(k+1)))**2)
62 coef2=log(pre_ma(k)/pre)/(pre_ma(k+1)* &
63 (log(pre_ma(k)/pre_ma(k+1)))**2)
64 alpha = log(pre_ma(k)/pre)/log(pre_ma(k)/pre_ma(k+1))
65
66 ADJ_tv_ma(k) = ADJ_tv_ma(k) + (1.-alpha)*ADJ_tv
67 ADJ_alpha = (tv_ma(k+1)-tv_ma(k))*ADJ_tv
68 ADJ_tv_ma(k+1) = ADJ_tv_ma(k+1) + alpha * ADJ_tv
69
70 ADJ_pre_ma(k) = ADJ_pre_ma(k) + coef1 * ADJ_alpha
71 ADJ_pre_ma(k+1) = ADJ_pre_ma(k+1) + coef2 * ADJ_alpha
72 exit
73 end if
74 end do
75
76 end if
77 ADJ_tv = 0.
78
79 end subroutine da_find_layer_adj
80
81