-
Notifications
You must be signed in to change notification settings - Fork 156
/
Copy pathsmooth_tends.f90
executable file
·263 lines (214 loc) · 6.68 KB
/
smooth_tends.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
subroutine smooth_tends(u,v,t,q,oz,cw,ps,nlevs,mype)
!
! prgram history log:
! 2009-10-1 kleist - define filtering of tendencies
! 2010-9-30 rancic - start modifying code with new coding infrastructure
!
! [At this stage, the idea is to exploit existing coding structure and
! quickly generate some results. To this end, sub2grid2, etc., will be
! used in smoothing of all fields. At some later stage, the new sub2grid,
! etc., using bundle structure will be applied. - Rancic, Jan 2011]
! [Also, note that here we use surface pressure instead of 3d preasure.
! This also may be changed later on - but for now hydrostatic modeling
! paradigm should be more efficient (once the bundle is applied). - Rancic, Jan 2011]
!
use kinds, only: r_kind,i_kind
use constants, only: zero
use gridmod, only: nlat,nlon,lat2,lon2,nsig,sp_a,grd_a
use mpimod, only: nnnvsbal
implicit none
! Passed variables
integer(i_kind),intent(in):: mype
integer(i_kind),intent(in):: nlevs
real(r_kind),dimension(lat2,lon2),intent(inout):: ps
real(r_kind),dimension(lat2,lon2,nsig),intent(inout):: t,q,cw,oz,u,v
! Local Variables
real(r_kind),dimension(lat2,lon2,nsig+1):: pri_dum
real(r_kind),dimension(lat2,lon2):: t_dum
real(r_kind),dimension(nlat,nlon,nnnvsbal):: hwork
real(r_kind),dimension(sp_a%nc):: spc1
integer(i_kind) iflg,k,i,j
iflg=1
!
! First pass (smooth tendencies of balanced variables)
!
hwork=zero
pri_dum(:,:,1:nsig)=q(:,:,1:nsig)
pri_dum(:,:,nsig+1)=ps(:,:)
! Put tendencies on grid
call sub2grid2(hwork,u,v,pri_dum,t,iflg)
! Transform to spectral coefficients
do k=1,nnnvsbal
call general_g2s0(grd_a,sp_a,spc1,hwork(1,1,k))
! Truncate
call general_jcaptrans(sp_a,spc1)
! Back to the grid
call general_s2g0(grd_a,sp_a,spc1,hwork(1,1,k))
end do !End do k
! Put back on subdomains
call grid2sub2(hwork,u,v,pri_dum,t)
q(:,:,1:nsig)=pri_dum(:,:,1:nsig)
ps(:,:)=pri_dum(:,:,nsig+1)
!
! Second pass (smooth tendencies of remaining variables)
!
hwork=zero
pri_dum=zero
t_dum=zero
! Put tendencies on grid
call sub2grid2(hwork,oz,cw,pri_dum,t_dum,iflg)
! Transform to spectral coefficients
do k=1,nnnvsbal
call general_g2s0(grd_a,sp_a,spc1,hwork(1,1,k))
! Truncate
call general_jcaptrans(sp_a,spc1)
! Back to the grid
call general_s2g0(grd_a,sp_a,spc1,hwork(1,1,k))
end do !End do k
! Put back on subdomains
call grid2sub2(hwork,oz,cw,pri_dum,t_dum)
return
end subroutine smooth_tends
subroutine smooth_tends_ad(u,v,t,q,oz,cw,ps,nlevs,mype)
!
! prgram history log:
! 2009-10-1 kleist - define filtering of tendencies
! 2010-9-30 rancic - start modifying code with new coding infrastructure
!
use kinds, only: r_kind,i_kind
use constants, only: zero
use gridmod, only: nlat,nlon,lat2,lon2,nsig,sp_a,grd_a
use mpimod, only: nnnvsbal
implicit none
! Passed variables
integer(i_kind),intent(in):: mype
integer(i_kind),intent(in):: nlevs
real(r_kind),dimension(lat2,lon2),intent(inout):: ps
real(r_kind),dimension(lat2,lon2,nsig),intent(inout):: t,q,cw,oz,u,v
! Local Variables
real(r_kind),dimension(lat2,lon2,nsig+1):: pri_dum
real(r_kind),dimension(lat2,lon2):: t_dum
real(r_kind),dimension(nlat,nlon,nnnvsbal):: hwork
real(r_kind),dimension(sp_a%nc):: spc1
integer(i_kind) iflg,k,i,j
iflg=1
!
! First pass (smooth tendencies of balanced variables)
!
hwork=zero
pri_dum(:,:,1:nsig)=q(:,:,1:nsig)
pri_dum(:,:,nsig+1)=ps(:,:)
! Put tendencies on grid
call sub2grid2(hwork,u,v,pri_dum,t,iflg)
! Transform to spectral coefficients
do k=1,nnnvsbal
call general_s2g0_ad(grd_a,sp_a,spc1,hwork(1,1,k))
! Truncate
call general_jcaptrans_ad(sp_a,spc1)
! Back to the grid
call general_g2s0_ad(grd_a,sp_a,spc1,hwork(1,1,k))
end do !End do k
! Put back on subdomains
call grid2sub2(hwork,u,v,pri_dum,t)
q(:,:,1:nsig)=pri_dum(:,:,1:nsig)
ps(:,:)=pri_dum(:,:,nsig+1)
!
! Second pass (smooth tendencies of remaining variables)
!
hwork=zero
pri_dum=zero
t_dum=zero
! Put tendencies on grid
call sub2grid2(hwork,oz,cw,pri_dum,t_dum,iflg)
! Transform to spectral coefficients
do k=1,nnnvsbal
call general_s2g0_ad(grd_a,sp_a,spc1,hwork(1,1,k))
! Truncate
call general_jcaptrans_ad(sp_a,spc1)
! Back to the grid
call general_g2s0_ad(grd_a,sp_a,spc1,hwork(1,1,k))
end do !End do k
! Put back on subdomains
call grid2sub2(hwork,oz,cw,pri_dum,t_dum)
return
end subroutine smooth_tends_ad
subroutine general_jcaptrans(sp,z)
use kinds, only: r_kind
use constants, only: zero
use general_specmod, only: spec_vars
implicit none
type(spec_vars),intent(in ) :: sp
real(r_kind),dimension(sp%nc),intent(inout) :: z
real(r_kind),dimension(sp%nc_trunc):: z_trunc
integer j,l,n,ks1,ks2
! First load temp array
do l=0,sp%jcap_trunc
do n=l,sp%jcap_trunc
ks2=l*(2*sp%jcap_trunc+(-1)*(l-1))+2*n
if (l.LE.sp%jcap.AND.n.LE.sp%jcap) then
ks1=l*(2*sp%jcap+(-1)*(l-1))+2*n
z_trunc(ks2+1)=z(ks1+1)
z_trunc(ks2+2)=z(ks1+2)
else
z_trunc(ks2+1)=zero
z_trunc(ks2+2)=zero
endif
end do
end do
! Now put temp array into full truncation and pad
do l=0,sp%jcap
do n=l,sp%jcap
ks2=l*(2*sp%jcap+(-1)*(l-1))+2*n
if (l.LE.sp%jcap_trunc.AND.n.LE.sp%jcap_trunc) then
ks1=l*(2*sp%jcap_trunc+(-1)*(l-1))+2*n
z(ks2+1)=z_trunc(ks1+1)
z(ks2+2)=z_trunc(ks1+2)
else
z(ks2+1)=zero
z(ks2+2)=zero
endif
end do
end do
return
end subroutine general_jcaptrans
subroutine general_jcaptrans_ad(sp,z)
use kinds, only: r_kind
use constants, only: zero
use general_specmod, only: spec_vars
implicit none
type(spec_vars),intent(in ) :: sp
real(r_kind),dimension(sp%nc),intent(inout) :: z
real(r_kind),dimension(sp%nc_trunc):: z_trunc
integer j,l,n,ks1,ks2
z_trunc=zero
! Now put temp array into full truncation and pad
do l=0,sp%jcap
do n=l,sp%jcap
ks2=l*(2*sp%jcap+(-1)*(l-1))+2*n
if (l.LE.sp%jcap_trunc.AND.n.LE.sp%jcap_trunc) then
ks1=l*(2*sp%jcap_trunc+(-1)*(l-1))+2*n
z_trunc(ks1+1)=z_trunc(ks1+1)+z(ks2+1)
z_trunc(ks1+2)=z_trunc(ks1+2)+z(ks2+2)
else
z_trunc(ks2+1)=zero
z_trunc(ks2+2)=zero
endif
end do
end do
z=zero
! First load temp array
do l=0,sp%jcap_trunc
do n=l,sp%jcap_trunc
ks2=l*(2*sp%jcap_trunc+(-1)*(l-1))+2*n
if (l.LE.sp%jcap.AND.n.LE.sp%jcap) then
ks1=l*(2*sp%jcap+(-1)*(l-1))+2*n
z(ks1+1)=z(ks1+1)+z_trunc(ks2+1)
z(ks1+2)=z(ks1+2)+z_trunc(ks2+2)
else
z(ks1+1)=zero
z(ks2+2)=zero
endif
end do
end do
return
end subroutine general_jcaptrans_ad