da_psichi_to_uv_adj.inc
References to this file elsewhere.
1 subroutine da_psichi_to_uv_adj(u, v, coefx, coefy, psi, chi)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Adjoint code of da_psichi_to_uv
5 !
6 ! Method: u = coefx * (-dpsi/dy + dchi/dx)
7 ! v = coefy * ( dpsi/dx + dchi/dy)
8 !
9 ! Assumptions: Unstaggered grid.
10 ! Lateral boundary conditions - dpsi/dn, dchi/dn = 0 (FCT)
11 ! grid size may or may not be equal
12 !
13 !---------------------------------------------------------------------------
14
15 implicit none
16
17 real, intent(inout) :: u(ims:ime,jms:jme,kms:kme) ! u wind comp.
18 real, intent(inout) :: v(ims:ime,jms:jme,kms:kme) ! v wind comp.
19 real, intent(inout) :: psi(ims:ime,jms:jme,kms:kme) ! Stream function
20 real, intent(inout) :: chi(ims:ime,jms:jme,kms:kme) ! Velocity potential
21 real, intent(in) :: coefx(ims:ime,jms:jme) ! Multiplicative coeff.
22 real, intent(in) :: coefy(ims:ime,jms:jme) ! Multiplicative coeff.
23
24 integer :: i, j, k ! Loop counters.
25 integer :: is, ie ! 1st dim. end points.
26 integer :: js, je ! 2nd dim. end points.
27
28 if (trace_use) call da_trace_entry("da_psichi_to_uv_adj")
29
30 !---------------------------------------------------------------------------
31 ! [5.0] For Global application, set Wast-Eest Periodic boundary
32 !---------------------------------------------------------------------------
33
34 if (global) then
35 call da_set_boundary_3d(u)
36 call da_set_boundary_3d(v)
37 end if
38
39 !---------------------------------------------------------------------------
40 ! Computation to check for edge of domain:
41 !---------------------------------------------------------------------------
42
43 is = its-1
44 ie = ite+1
45 js = jts-1
46 je = jte+1
47 if (jts == jds) js = jds+1
48 if (jte == jde) je = jde-1
49
50 if (.not. global) then
51 if (its == ids) is = ids+1
52 if (ite == ide) ie = ide-1
53
54 do k = kts, kte
55
56 !---------------------------------------------------------------------
57 ! [4.0] Corner points (assume average of surrounding points - poor?):
58 !---------------------------------------------------------------------
59
60 ! [4.1] Bottom-left point:
61
62 if (its == ids .AND. jts == jds) then
63 u(its+1,jts,k) = u(its+1,jts,k) + 0.5 * u(its,jts,k)
64 u(its,jts+1,k) = u(its,jts+1,k) + 0.5 * u(its,jts,k)
65 v(its+1,jts,k) = v(its+1,jts,k) + 0.5 * v(its,jts,k)
66 v(its,jts+1,k) = v(its,jts+1,k) + 0.5 * v(its,jts,k)
67 end if
68
69 ! [4.2] Top-left point:
70
71 if (ite == ide .AND. jts == jds) then
72 u(ite-1,jts,k) = u(ite-1,jts,k) + 0.5 * u(ite,jts,k)
73 u(ite,jts+1,k) = u(ite,jts+1,k) + 0.5 * u(ite,jts,k)
74 v(ite-1,jts,k) = v(ite-1,jts,k) + 0.5 * v(ite,jts,k)
75 v(ite,jts+1,k) = v(ite,jts+1,k) + 0.5 * v(ite,jts,k)
76 end if
77
78 ! [4.3] Bottom-right point:
79
80 if (its == ids .AND. jte == jde) then
81 u(its+1,jte,k) = u(its+1,jte,k) + 0.5 * u(its,jte,k)
82 u(its,jte-1,k) = u(its,jte-1,k) + 0.5 * u(its,jte,k)
83 v(its+1,jte,k) = v(its+1,jte,k) + 0.5 * v(its,jte,k)
84 v(its,jte-1,k) = v(its,jte-1,k) + 0.5 * v(its,jte,k)
85 end if
86
87 ! [4.4] Top-right point:
88
89 if (ite == ide .AND. jte == jde) then
90 u(ite-1,jte,k) = u(ite-1,jte,k) + 0.5 * u(ite,jte,k)
91 u(ite,jte-1,k) = u(ite,jte-1,k) + 0.5 * u(ite,jte,k)
92 v(ite-1,jte,k) = v(ite-1,jte,k) + 0.5 * v(ite,jte,k)
93 v(ite,jte-1,k) = v(ite,jte-1,k) + 0.5 * v(ite,jte,k)
94 end if
95 end do
96 end if
97
98 ! [3.0] Compute u, v at domain boundaries:
99
100 do k = kts, kte
101 if (.not. global) then
102 ! [3.4] Northern boundaries:
103
104 if (jte == jde) then
105 j = jte
106
107 do i = ie, is, -1
108 psi(i+1,j,k) = psi(i+1,j,k) + coefx(i,j) * v(i,j,k)
109 psi(i-1,j,k) = psi(i-1,j,k) - coefx(i,j) * v(i,j,k)
110 chi(i,j ,k) = chi(i,j ,k) + coefy(i,j) * v(i,j,k)
111 chi(i,j-2,k) = chi(i,j-2,k) - coefy(i,j) * v(i,j,k)
112
113 psi(i,j ,k) = psi(i,j ,k) - coefy(i,j) * u(i,j,k)
114 psi(i,j-2,k) = psi(i,j-2,k) + coefy(i,j) * u(i,j,k)
115 chi(i+1,j,k) = chi(i+1,j,k) + coefx(i,j) * u(i,j,k)
116 chi(i-1,j,k) = chi(i-1,j,k) - coefx(i,j) * u(i,j,k)
117 end do
118 end if
119
120 ! [3.3] Southern boundaries:
121
122 if (jts == jds) then
123 j = jts
124
125 do i = ie, is, -1
126
127
128 psi(i+1,j,k) = psi(i+1,j,k) + coefx(i,j) * v(i,j,k)
129 psi(i-1,j,k) = psi(i-1,j,k) - coefx(i,j) * v(i,j,k)
130 chi(i,j+2,k) = chi(i,j+2,k) + coefy(i,j) * v(i,j,k)
131 chi(i,j ,k) = chi(i,j ,k) - coefy(i,j) * v(i,j,k)
132
133 psi(i,j+2,k) = psi(i,j+2,k) - coefy(i,j) * u(i,j,k)
134 psi(i,j ,k) = psi(i,j ,k) + coefy(i,j) * u(i,j,k)
135 chi(i+1,j,k) = chi(i+1,j,k) + coefx(i,j) * u(i,j,k)
136 chi(i-1,j,k) = chi(i-1,j,k) - coefx(i,j) * u(i,j,k)
137
138 end do
139 end if
140
141 ! [3.2] Eastern boundaries:
142
143 if (ite == ide) then
144 i = ite
145 do j = je, js, -1
146 psi(i ,j,k) = psi(i ,j,k) + coefx(i,j) * v(i,j,k)
147 psi(i-2,j,k) = psi(i-2,j,k) - coefx(i,j) * v(i,j,k)
148 chi(i,j+1,k) = chi(i,j+1,k) + coefy(i,j) * v(i,j,k)
149 chi(i,j-1,k) = chi(i,j-1,k) - coefy(i,j) * v(i,j,k)
150
151 psi(i,j+1,k) = psi(i,j+1,k) - coefy(i,j) * u(i,j,k)
152 psi(i,j-1,k) = psi(i,j-1,k) + coefy(i,j) * u(i,j,k)
153 chi(i ,j,k) = chi(i ,j,k) + coefx(i,j) * u(i,j,k)
154 chi(i-2,j,k) = chi(i-2,j,k) - coefx(i,j) * u(i,j,k)
155 end do
156 end if
157
158 ! [3.1] Western boundaries:
159 if (its == ids) then
160 i = its
161
162 do j = je, js, -1
163 psi(i+2,j,k) = psi(i+2,j,k) + coefx(i,j) * v(i,j,k)
164 psi(i ,j,k) = psi(i ,j,k) - coefx(i,j) * v(i,j,k)
165 chi(i,j+1,k) = chi(i,j+1,k) + coefy(i,j) * v(i,j,k)
166 chi(i,j-1,k) = chi(i,j-1,k) - coefy(i,j) * v(i,j,k)
167
168 psi(i,j+1,k) = psi(i,j+1,k) - coefy(i,j) * u(i,j,k)
169 psi(i,j-1,k) = psi(i,j-1,k) + coefy(i,j) * u(i,j,k)
170 chi(i+2,j,k) = chi(i+2,j,k) + coefx(i,j) * u(i,j,k)
171 chi(i, j,k) = chi(i, j,k) - coefx(i,j) * u(i,j,k)
172 end do
173 end if
174 end if
175
176 !------------------------------------------------------------------------
177 ! [2.0] Compute u, v at interior points (2nd order central finite diffs):
178 !------------------------------------------------------------------------
179
180 do j = je, js, -1
181 do i = ie, is, -1
182 psi(i+1,j,k) = psi(i+1,j,k) + coefx(i,j) * v(i,j,k)
183 psi(i-1,j,k) = psi(i-1,j,k) - coefx(i,j) * v(i,j,k)
184 chi(i,j+1,k) = chi(i,j+1,k) + coefy(i,j) * v(i,j,k)
185 chi(i,j-1,k) = chi(i,j-1,k) - coefy(i,j) * v(i,j,k)
186
187 psi(i,j+1,k) = psi(i,j+1,k) - coefy(i,j) * u(i,j,k)
188 psi(i,j-1,k) = psi(i,j-1,k) + coefy(i,j) * u(i,j,k)
189 chi(i+1,j,k) = chi(i+1,j,k) + coefx(i,j) * u(i,j,k)
190 chi(i-1,j,k) = chi(i-1,j,k) - coefx(i,j) * u(i,j,k)
191 end do
192 end do
193 end do ! loop over levels
194
195 !---------------------------------------------------------------------------
196 ! [5.0] For Global application, set Wast-Eest Periodic boundary
197 !---------------------------------------------------------------------------
198
199 if (global) then
200 call da_set_boundary_3d(psi)
201 call da_set_boundary_3d(chi)
202 end if
203
204 if (trace_use) call da_trace_exit("da_psichi_to_uv_adj")
205
206 end subroutine da_psichi_to_uv_adj
207
208