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 !---------------------------------------------------------------------------
29 ! [5.0] For Global application, set Wast-Eest Periodic boundary
30 !---------------------------------------------------------------------------
31
32 if (global) then
33 call da_set_boundary_3d(u)
34 call da_set_boundary_3d(v)
35 end if
36
37 !---------------------------------------------------------------------------
38 ! Computation to check for edge of domain:
39 !---------------------------------------------------------------------------
40
41 is = its-1
42 ie = ite+1
43 js = jts-1
44 je = jte+1
45 if (jts == jds) js = jds+1
46 if (jte == jde) je = jde-1
47
48 if (.not. global) then
49 if (its == ids) is = ids+1
50 if (ite == ide) ie = ide-1
51
52 do k = kts, kte
53
54 !---------------------------------------------------------------------
55 ! [4.0] Corner points (assume average of surrounding points - poor?):
56 !---------------------------------------------------------------------
57
58 ! [4.1] Bottom-left point:
59
60 if (its == ids .AND. jts == jds) then
61 u(its+1,jts,k) = u(its+1,jts,k) + 0.5 * u(its,jts,k)
62 u(its,jts+1,k) = u(its,jts+1,k) + 0.5 * u(its,jts,k)
63 v(its+1,jts,k) = v(its+1,jts,k) + 0.5 * v(its,jts,k)
64 v(its,jts+1,k) = v(its,jts+1,k) + 0.5 * v(its,jts,k)
65 end if
66
67 ! [4.2] Top-left point:
68
69 if (ite == ide .AND. jts == jds) then
70 u(ite-1,jts,k) = u(ite-1,jts,k) + 0.5 * u(ite,jts,k)
71 u(ite,jts+1,k) = u(ite,jts+1,k) + 0.5 * u(ite,jts,k)
72 v(ite-1,jts,k) = v(ite-1,jts,k) + 0.5 * v(ite,jts,k)
73 v(ite,jts+1,k) = v(ite,jts+1,k) + 0.5 * v(ite,jts,k)
74 end if
75
76 ! [4.3] Bottom-right point:
77
78 if (its == ids .AND. jte == jde) then
79 u(its+1,jte,k) = u(its+1,jte,k) + 0.5 * u(its,jte,k)
80 u(its,jte-1,k) = u(its,jte-1,k) + 0.5 * u(its,jte,k)
81 v(its+1,jte,k) = v(its+1,jte,k) + 0.5 * v(its,jte,k)
82 v(its,jte-1,k) = v(its,jte-1,k) + 0.5 * v(its,jte,k)
83 end if
84
85 ! [4.4] Top-right point:
86
87 if (ite == ide .AND. jte == jde) then
88 u(ite-1,jte,k) = u(ite-1,jte,k) + 0.5 * u(ite,jte,k)
89 u(ite,jte-1,k) = u(ite,jte-1,k) + 0.5 * u(ite,jte,k)
90 v(ite-1,jte,k) = v(ite-1,jte,k) + 0.5 * v(ite,jte,k)
91 v(ite,jte-1,k) = v(ite,jte-1,k) + 0.5 * v(ite,jte,k)
92 end if
93 end do
94 end if
95
96 ! [3.0] Compute u, v at domain boundaries:
97
98 do k = kts, kte
99 if (.not. global) then
100 ! [3.4] Northern boundaries:
101
102 if (jte == jde) then
103 j = jte
104
105 do i = ie, is, -1
106 psi(i+1,j,k) = psi(i+1,j,k) + coefx(i,j) * v(i,j,k)
107 psi(i-1,j,k) = psi(i-1,j,k) - coefx(i,j) * v(i,j,k)
108 chi(i,j ,k) = chi(i,j ,k) + coefy(i,j) * v(i,j,k)
109 chi(i,j-2,k) = chi(i,j-2,k) - coefy(i,j) * v(i,j,k)
110
111 psi(i,j ,k) = psi(i,j ,k) - coefy(i,j) * u(i,j,k)
112 psi(i,j-2,k) = psi(i,j-2,k) + coefy(i,j) * u(i,j,k)
113 chi(i+1,j,k) = chi(i+1,j,k) + coefx(i,j) * u(i,j,k)
114 chi(i-1,j,k) = chi(i-1,j,k) - coefx(i,j) * u(i,j,k)
115 end do
116 end if
117
118 ! [3.3] Southern boundaries:
119
120 if (jts == jds) then
121 j = jts
122
123 do i = ie, is, -1
124
125
126 psi(i+1,j,k) = psi(i+1,j,k) + coefx(i,j) * v(i,j,k)
127 psi(i-1,j,k) = psi(i-1,j,k) - coefx(i,j) * v(i,j,k)
128 chi(i,j+2,k) = chi(i,j+2,k) + coefy(i,j) * v(i,j,k)
129 chi(i,j ,k) = chi(i,j ,k) - coefy(i,j) * v(i,j,k)
130
131 psi(i,j+2,k) = psi(i,j+2,k) - coefy(i,j) * u(i,j,k)
132 psi(i,j ,k) = psi(i,j ,k) + coefy(i,j) * u(i,j,k)
133 chi(i+1,j,k) = chi(i+1,j,k) + coefx(i,j) * u(i,j,k)
134 chi(i-1,j,k) = chi(i-1,j,k) - coefx(i,j) * u(i,j,k)
135
136 end do
137 end if
138
139 ! [3.2] Eastern boundaries:
140
141 if (ite == ide) then
142 i = ite
143 do j = je, js, -1
144 psi(i ,j,k) = psi(i ,j,k) + coefx(i,j) * v(i,j,k)
145 psi(i-2,j,k) = psi(i-2,j,k) - coefx(i,j) * v(i,j,k)
146 chi(i,j+1,k) = chi(i,j+1,k) + coefy(i,j) * v(i,j,k)
147 chi(i,j-1,k) = chi(i,j-1,k) - coefy(i,j) * v(i,j,k)
148
149 psi(i,j+1,k) = psi(i,j+1,k) - coefy(i,j) * u(i,j,k)
150 psi(i,j-1,k) = psi(i,j-1,k) + coefy(i,j) * u(i,j,k)
151 chi(i ,j,k) = chi(i ,j,k) + coefx(i,j) * u(i,j,k)
152 chi(i-2,j,k) = chi(i-2,j,k) - coefx(i,j) * u(i,j,k)
153 end do
154 end if
155
156 ! [3.1] Western boundaries:
157 if (its == ids) then
158 i = its
159
160 do j = je, js, -1
161 psi(i+2,j,k) = psi(i+2,j,k) + coefx(i,j) * v(i,j,k)
162 psi(i ,j,k) = psi(i ,j,k) - coefx(i,j) * v(i,j,k)
163 chi(i,j+1,k) = chi(i,j+1,k) + coefy(i,j) * v(i,j,k)
164 chi(i,j-1,k) = chi(i,j-1,k) - coefy(i,j) * v(i,j,k)
165
166 psi(i,j+1,k) = psi(i,j+1,k) - coefy(i,j) * u(i,j,k)
167 psi(i,j-1,k) = psi(i,j-1,k) + coefy(i,j) * u(i,j,k)
168 chi(i+2,j,k) = chi(i+2,j,k) + coefx(i,j) * u(i,j,k)
169 chi(i, j,k) = chi(i, j,k) - coefx(i,j) * u(i,j,k)
170 end do
171 end if
172 end if
173
174 !------------------------------------------------------------------------
175 ! [2.0] Compute u, v at interior points (2nd order central finite diffs):
176 !------------------------------------------------------------------------
177
178 do j = je, js, -1
179 do i = ie, is, -1
180 psi(i+1,j,k) = psi(i+1,j,k) + coefx(i,j) * v(i,j,k)
181 psi(i-1,j,k) = psi(i-1,j,k) - coefx(i,j) * v(i,j,k)
182 chi(i,j+1,k) = chi(i,j+1,k) + coefy(i,j) * v(i,j,k)
183 chi(i,j-1,k) = chi(i,j-1,k) - coefy(i,j) * v(i,j,k)
184
185 psi(i,j+1,k) = psi(i,j+1,k) - coefy(i,j) * u(i,j,k)
186 psi(i,j-1,k) = psi(i,j-1,k) + coefy(i,j) * u(i,j,k)
187 chi(i+1,j,k) = chi(i+1,j,k) + coefx(i,j) * u(i,j,k)
188 chi(i-1,j,k) = chi(i-1,j,k) - coefx(i,j) * u(i,j,k)
189 end do
190 end do
191 end do ! loop over levels
192
193 !---------------------------------------------------------------------------
194 ! [5.0] For Global application, set Wast-Eest Periodic boundary
195 !---------------------------------------------------------------------------
196
197 if (global) then
198 call da_set_boundary_3d(psi)
199 call da_set_boundary_3d(chi)
200 end if
201
202 end subroutine da_psichi_to_uv_adj
203
204