da_balance_cycloterm_lin.inc
References to this file elsewhere.
1 subroutine da_balance_cycloterm_lin( rho, ub, vb, u, v, coefx, coefy, term_x, term_y)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Calculates linearised cyclostrophic term in balance equation.
5 !
6 ! Method: Term is rho (u.grad) u on a single level.
7 !---------------------------------------------------------------------------
8
9 implicit none
10
11 real, intent(in) :: rho(ims:ime,jms:jme) ! Density.
12 real, intent(in) :: ub(ims:ime,jms:jme) ! Background u wind
13 real, intent(in) :: vb(ims:ime,jms:jme) ! Background u wind
14 real, intent(in) :: u(ims:ime,jms:jme) ! u wind increment
15 real, intent(in) :: v(ims:ime,jms:jme) ! v wind increment
16 real, intent(in) :: coefx(ims:ime,jms:jme) ! Multiplicative const.
17 real, intent(in) :: coefy(ims:ime,jms:jme)
18 real, intent(inout) :: term_x(ims:ime,jms:jme) ! x component of term.
19 real, intent(inout) :: term_y(ims:ime,jms:jme) ! y component of term.
20
21 integer :: i, j ! Loop counters.
22 integer :: is, ie ! 1st dim. end points.
23 integer :: js, je ! 2nd dim. end points.
24
25 real :: data(ims:ime,jms:jme) ! Temporary storage.
26
27 real :: varb, var
28
29 if (trace_use) call da_trace_entry("da_balance_cycloterm_lin")
30
31 !---------------------------------------------------------------------------
32 ! [1.0] Initialise:
33 !---------------------------------------------------------------------------
34
35 ! Computation to check for edge of domain:
36 is = its; ie = ite; js = jts; je = jte
37 if (.not. global .and. its == ids) is = ids+1
38 if (.not. global .and. ite == ide) ie = ide-1
39 if (jts == jds) js = jds+1; if (jte == jde) je = jde-1
40
41 !---------------------------------------------------------------------------
42 ! [2.0] Calculate term_x = rho M ( u'du/dx + v'du/dy + udu'/dx + vdu'/dy):
43 !---------------------------------------------------------------------------
44
45 ! [2.1] Interior points:
46
47 do j = js, je
48 do i = is, ie
49 data(i,j) = coefx(i,j)*u(i,j) * ( ub(i+1,j) - ub(i-1,j)) + &
50 coefy(i,j)*v(i,j) * ( ub(i,j+1) - ub(i,j-1)) + &
51 coefx(i,j)*ub(i,j) * ( u(i+1,j) - u(i-1,j)) + &
52 coefy(i,j)*vb(i,j) * ( u(i,j+1) - u(i,j-1))
53 end do
54 end do
55
56 if (.NOT. global) then ! For global only interior points needed
57
58 ! [2.2] Bottom boundaries:
59
60 if (its == ids) then
61 i = its
62
63 do j = js, je
64 var = -3.0*u (i,j)+4.0*u(i+1,j)-u(i+2,j)
65 varb = -3.0*ub(i,j)+4.0*ub(i+1,j)-ub(i+2,j)
66
67 data(i,j) = coefx(i,j)* u(i,j) * varb + &
68 coefy(i,j)* v(i,j) * ( ub(i,j+1) - ub(i,j-1)) + &
69 coefx(i,j)*ub(i,j) * var + &
70 coefy(i,j)*vb(i,j) * ( u(i,j+1) - u(i,j-1))
71 end do
72 end if
73
74 ! [2.3] Top boundaries:
75
76 if (ite == ide) then
77 i = ite
78
79 do j = js, je
80 var = 3.0*u (i,j)-4.0*u (i-1,j)+u (i-2,j)
81 varb = 3.0*ub(i,j)-4.0*ub(i-1,j)+ub(i-2,j)
82
83 data(i,j) = coefx(i,j)*u(i,j) * varb + &
84 coefy(i,j)*v(i,j) * ( ub(i,j+1) - ub(i,j-1)) + &
85 coefx(i,j)*ub(i,j) * var + &
86 coefy(i,j)*vb(i,j) * ( u(i,j+1) - u(i,j-1))
87 end do
88 end if
89
90 ! [2.4] Left boundaries:
91
92 if (jts == jds) then
93 j = jts
94
95 do i = is, ie
96 var = -3.0*u (i,j)+4.0*u (i,j+1)-u (i,j+2)
97 varb = -3.0*ub(i,j)+4.0*ub(i,j+1)-ub(i,j+2)
98
99 data(i,j) = coefx(i,j)*u(i,j) * ( ub(i+1,j) - ub(i-1,j)) + &
100 coefy(i,j)*v(i,j) * varb + &
101 coefx(i,j)*ub(i,j) * ( u(i+1,j) - u(i-1,j)) + &
102 coefy(i,j)*vb(i,j) * var
103 end do
104 end if
105
106 ! [2.5] Right boundaries:
107
108 if (jte == jde) then
109 j = jte
110
111 do i = is, ie
112 var = 3.0*u (i,j)-4.0*u (i,j-1)+u (i,j-2)
113 varb = 3.0*ub(i,j)-4.0*ub(i,j-1)+ub(i,j-2)
114
115 data(i,j) = coefx(i,j)*u(i,j) * ( ub(i+1,j) - ub(i-1,j)) + &
116 coefy(i,j)*v(i,j) * varb + &
117 coefx(i,j)*ub(i,j) * ( u(i+1,j) - u(i-1,j)) + &
118 coefy(i,j)*vb(i,j) * var
119 end do
120 end if
121
122 ! [2.6] Corner points:
123
124 if (its == ids .AND. jts == jds) then
125 data(its,jts) = 0.5 * ( data(its,jts+1) + data(its+1,jts))
126 end if
127
128 if (ite == ide .AND. jts == jds) then
129 data(ite,jts) = 0.5 * ( data(ite-1,jts) + data(ite,jts+1))
130 end if
131
132 if (its == ids .AND. jte == jde) then
133 data(its,jde) = 0.5 * ( data(its,jde-1) + data(its+1,jde))
134 end if
135
136 if (ite == ide .AND. jte == jde) then
137 data(ite,jte) = 0.5 * ( data(ite-1,jte) + data(ite,jte-1))
138 end if
139 end if ! not global
140
141 ! [2.7] Multiply by rho and add to term_x:
142
143 term_x(its:ite,jts:jte) = rho(its:ite,jts:jte) * data(its:ite,jts:jte) + &
144 term_x(its:ite,jts:jte)
145
146 !---------------------------------------------------------------------------
147 ! [3.0] Calculate term_y = rho M ( u'dv/dx + v'dv/dy + udv'/dx + vdv'/dy):
148 !---------------------------------------------------------------------------
149
150 ! [3.1] Interior points:
151
152 do j = js, je
153 do i = is, ie
154 data(i,j) = coefx(i,j)*u(i,j) * ( vb(i+1,j) - vb(i-1,j)) + &
155 coefy(i,j)*v(i,j) * ( vb(i,j+1) - vb(i,j-1)) + &
156 coefx(i,j)*ub(i,j) * ( v(i+1,j) - v(i-1,j)) + &
157 coefy(i,j)*vb(i,j) * ( v(i,j+1) - v(i,j-1))
158 end do
159 end do
160
161 if (.NOT. global) then ! For global only interior points needed
162 ! [3.2] Bottom boundaries:
163
164 if (its == ids) then
165 i = its
166
167 do j = js, je
168 var = -3.0*v (i,j)+4.0*v (i+1,j)-v (i+2,j)
169 varb = -3.0*vb(i,j)+4.0*vb(i+1,j)-vb(i+2,j)
170
171 data(i,j) = coefx(i,j)*u(i,j) * varb + &
172 coefy(i,j)*v(i,j) * ( vb(i,j+1) - vb(i,j-1)) + &
173 coefx(i,j)*ub(i,j) * var + &
174 coefy(i,j)*vb(i,j) * ( v(i,j+1) - v(i,j-1))
175 end do
176 end if
177
178 ! [3.3] Top boundaries:
179
180 if (ite == ide) then
181 i = ite
182
183 do j = js, je
184 var = 3.0*v (i,j)-4.0*v (i-1,j)+v (i-2,j)
185 varb = 3.0*vb(i,j)-4.0*vb(i-1,j)+vb(i-2,j)
186
187 data(i,j) = coefx(i,j)*u(i,j) * varb + &
188 coefy(i,j)*v(i,j) * ( vb(i,j+1) - vb(i,j-1)) + &
189 coefx(i,j)*ub(i,j) * var + &
190 coefy(i,j)*vb(i,j) * ( v(i,j+1) - v(i,j-1))
191 end do
192 end if
193
194 ! [3.4] Left boundaries:
195
196 if (jts == jds) then
197 j = jts
198
199 do i = is, ie
200 varb = -3.0*vb(i,j)+4.0*vb(i,j+1)-vb(i,j+2)
201 var = -3.0*v (i,j)+4.0*v (i,j+1)-v (i,j+2)
202
203 data(i,j) = coefx(i,j)*u(i,j) * ( vb(i+1,j) - vb(i-1,j)) + &
204 coefy(i,j)*v(i,j) * varb + &
205 coefx(i,j)*ub(i,j) * ( v(i+1,j) - v(i-1,j)) + &
206 coefy(i,j)*vb(i,j) * var
207 end do
208 end if
209
210 ! [3.5] Right boundaries:
211
212 if (jte == jde) then
213 j = jte
214
215 do i = is, ie
216 varb = 3.0*vb(i,j)-4.0*vb(i,j-1)+vb(i,j-2)
217 var = 3.0*v (i,j)-4.0*v (i,j-1)+v (i,j-2)
218
219 data(i,j) = coefx(i,j)*u(i,j) * ( vb(i+1,j) - vb(i-1,j)) + &
220 coefy(i,j)*v(i,j) * varb + &
221 coefx(i,j)*ub(i,j) * ( v(i+1,j) - v(i-1,j)) + &
222 coefy(i,j)*vb(i,j) * var
223 end do
224 end if
225
226 ! [3.6] Corner points:
227
228 if (its == ids .AND. jts == jds) then
229 data(its,jts) = 0.5 * ( data(its,jts+1) + data(its+1,jts))
230 end if
231
232 if (ite == ide .AND. jts == jds) then
233 data(ite,jts) = 0.5 * ( data(ite-1,jts) + data(ite,jts+1))
234 end if
235
236 if (its == ids .AND. jte == jde) then
237 data(its,jde) = 0.5 * ( data(its,jde-1) + data(its+1,jde))
238 end if
239
240 if (ite == ide .AND. jte == jde) then
241 data(ite,jte) = 0.5 * ( data(ite-1,jte) + data(ite,jte-1))
242 end if
243 end if ! not global
244
245 ! [3.7] Multiply by rho and add to term_y
246
247 term_y(its:ite,jts:jte) = rho(its:ite,jts:jte) * data(its:ite,jts:jte) + term_y(its:ite,jts:jte)
248
249 if (trace_use) call da_trace_exit("da_balance_cycloterm_lin")
250
251 end subroutine da_balance_cycloterm_lin
252
253