source: palm/trunk/SOURCE/diffusion_s.f90 @ 23

Last change on this file since 23 was 20, checked in by raasch, 18 years ago

new parameter use_top_fluxes, Bugfix: ddzw dimensioned 1:nzt+1

  • Property svn:keywords set to Id
File size: 10.4 KB
Line 
1 MODULE diffusion_s_mod
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6! Bugfix: ddzw dimensioned 1:nzt"+1"
7! Calculation extended for gridpoint nzt, fluxes can be given at top,
8! +s_flux_t in parameter list, s_flux renamed s_flux_b,
9!
10! Former revisions:
11! -----------------
12! $Id: diffusion_s.f90 20 2007-02-26 00:12:32Z raasch $
13! RCS Log replace by Id keyword, revision history cleaned up
14!
15! Revision 1.8  2006/02/23 10:34:17  raasch
16! nzb_2d replaced by nzb_s_outer in horizontal diffusion and by nzb_s_inner
17! or nzb_diff_s_inner, respectively, in vertical diffusion, prescribed surface
18! fluxes at vertically oriented topography
19!
20! Revision 1.1  2000/04/13 14:54:02  schroeter
21! Initial revision
22!
23!
24! Description:
25! ------------
26! Diffusion term of scalar quantities (temperature and water content)
27!------------------------------------------------------------------------------!
28
29    PRIVATE
30    PUBLIC diffusion_s
31
32    INTERFACE diffusion_s
33       MODULE PROCEDURE diffusion_s
34       MODULE PROCEDURE diffusion_s_ij
35    END INTERFACE diffusion_s
36
37 CONTAINS
38
39
40!------------------------------------------------------------------------------!
41! Call for all grid points
42!------------------------------------------------------------------------------!
43    SUBROUTINE diffusion_s( ddzu, ddzw, kh, s, s_flux_b, s_flux_t, tend )
44
45       USE control_parameters
46       USE grid_variables
47       USE indices
48
49       IMPLICIT NONE
50
51       INTEGER ::  i, j, k
52       REAL    ::  vertical_gridspace
53       REAL    ::  ddzu(1:nzt+1), ddzw(1:nzt+1)
54       REAL    ::  tend(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1)
55       REAL, DIMENSION(:,:),   POINTER ::  s_flux_b, s_flux_t
56       REAL, DIMENSION(:,:,:), POINTER ::  kh, s
57
58       DO  i = nxl, nxr
59          DO  j = nys,nyn
60!
61!--          Compute horizontal diffusion
62             DO  k = nzb_s_outer(j,i)+1, nzt
63
64                tend(k,j,i) = tend(k,j,i)                                     &
65                                          + 0.5 * (                           &
66                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
67                      - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
68                                                  ) * ddx2                    &
69                                          + 0.5 * (                           &
70                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
71                      - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
72                                                  ) * ddy2
73             ENDDO
74
75!
76!--          Apply prescribed horizontal wall heatflux where necessary
77             IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0 ) ) &
78             THEN
79                DO  k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
80
81                   tend(k,j,i) = tend(k,j,i)                                  &
82                                          + 0.5 * ( fwxp(j,i) *               &
83                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
84                        - ( 1.0 - fwxp(j,i) ) * wall_heatflux(1)              &
85                                                   -fwxm(j,i) *               &
86                        ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
87                        + ( 1.0 - fwxm(j,i) ) * wall_heatflux(3)              &
88                                                  ) * ddx2                    &
89                                          + 0.5 * ( fwyp(j,i) *               &
90                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
91                        - ( 1.0 - fwyp(j,i) ) * wall_heatflux(2)              &
92                                                   -fwym(j,i) *               &
93                        ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
94                        + ( 1.0 - fwym(j,i) ) * wall_heatflux(4)              &
95                                                  ) * ddy2
96                ENDDO
97             ENDIF
98
99!
100!--          Compute vertical diffusion. In case that surface fluxes have been
101!--          prescribed or computed at bottom and/or top, index k starts/ends at
102!--          nzb+2 or nzt-1, respectively.
103             DO  k = nzb_diff_s_inner(j,i), nzt_diff
104
105                tend(k,j,i) = tend(k,j,i)                                     &
106                                       + 0.5 * (                              &
107            ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
108          - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k)   &
109                                               ) * ddzw(k)
110             ENDDO
111
112!
113!--          Vertical diffusion at the first computational gridpoint along
114!--          z-direction
115             IF ( use_surface_fluxes )  THEN
116
117                k = nzb_s_inner(j,i)+1
118
119                tend(k,j,i) = tend(k,j,i)                                     &
120                                       + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )    &
121                                               * ( s(k+1,j,i)-s(k,j,i) )      &
122                                               * ddzu(k+1)                    &
123                                           + s_flux_b(j,i)                    &
124                                         ) * ddzw(k)
125
126             ENDIF
127
128!
129!--          Vertical diffusion at the last computational gridpoint along
130!--          z-direction
131             IF ( use_top_fluxes )  THEN
132
133                k = nzt
134
135                tend(k,j,i) = tend(k,j,i)                                     &
136                                       + ( - s_flux_t(j,i)                    &
137                                           - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )  &
138                                                 * ( s(k,j,i)-s(k-1,j,i) )    &
139                                                 * ddzu(k)                    &
140                                         ) * ddzw(k)
141
142             ENDIF
143
144          ENDDO
145       ENDDO
146
147    END SUBROUTINE diffusion_s
148
149
150!------------------------------------------------------------------------------!
151! Call for grid point i,j
152!------------------------------------------------------------------------------!
153    SUBROUTINE diffusion_s_ij( i, j, ddzu, ddzw, kh, s, s_flux_b, s_flux_t, &
154                               tend )
155
156       USE control_parameters
157       USE grid_variables
158       USE indices
159
160       IMPLICIT NONE
161
162       INTEGER ::  i, j, k
163       REAL    ::  vertical_gridspace
164       REAL    ::  ddzu(1:nzt+1), ddzw(1:nzt+1)
165       REAL    ::  tend(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1)
166       REAL, DIMENSION(:,:),   POINTER ::  s_flux_b, s_flux_t
167       REAL, DIMENSION(:,:,:), POINTER ::  kh, s
168
169!
170!--    Compute horizontal diffusion
171       DO  k = nzb_s_outer(j,i)+1, nzt
172
173          tend(k,j,i) = tend(k,j,i)                                           &
174                                          + 0.5 * (                           &
175                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
176                      - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
177                                                  ) * ddx2                    &
178                                          + 0.5 * (                           &
179                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
180                      - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
181                                                  ) * ddy2
182       ENDDO
183
184!
185!--    Apply prescribed horizontal wall heatflux where necessary
186       IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0 ) ) &
187       THEN
188          DO  k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
189
190             tend(k,j,i) = tend(k,j,i)                                        &
191                                          + 0.5 * ( fwxp(j,i) *               &
192                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
193                        - ( 1.0 - fwxp(j,i) ) * wall_heatflux(1)              &
194                                                   -fwxm(j,i) *               &
195                        ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
196                        + ( 1.0 - fwxm(j,i) ) * wall_heatflux(3)              &
197                                                  ) * ddx2                    &
198                                          + 0.5 * ( fwyp(j,i) *               &
199                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
200                        - ( 1.0 - fwyp(j,i) ) * wall_heatflux(2)              &
201                                                   -fwym(j,i) *               &
202                        ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
203                        + ( 1.0 - fwym(j,i) ) * wall_heatflux(4)              &
204                                                  ) * ddy2
205          ENDDO
206       ENDIF
207
208!
209!--    Compute vertical diffusion. In case that surface fluxes have been
210!--    prescribed or computed at bottom and/or top, index k starts/ends at
211!--    nzb+2 or nzt-1, respectively.
212       DO  k = nzb_diff_s_inner(j,i), nzt_diff
213
214          tend(k,j,i) = tend(k,j,i)                                           &
215                                       + 0.5 * (                              &
216            ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
217          - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k)   &
218                                               ) * ddzw(k)
219       ENDDO
220
221!
222!--    Vertical diffusion at the first computational gridpoint along z-direction
223       IF ( use_surface_fluxes )  THEN
224
225          k = nzb_s_inner(j,i)+1
226
227          tend(k,j,i) = tend(k,j,i) + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )  &
228                                            * ( s(k+1,j,i)-s(k,j,i) )    &
229                                            * ddzu(k+1)                  &
230                                        + s_flux_b(j,i)                  &
231                                      ) * ddzw(k)
232
233       ENDIF
234
235!
236!--    Vertical diffusion at the last computational gridpoint along z-direction
237       IF ( use_top_fluxes )  THEN
238
239          k = nzt
240
241          tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(j,i)                  &
242                                      - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )  &
243                                            * ( s(k,j,i)-s(k-1,j,i) )    &
244                                            * ddzu(k)                    &
245                                      ) * ddzw(k)
246
247       ENDIF
248
249    END SUBROUTINE diffusion_s_ij
250
251 END MODULE diffusion_s_mod
Note: See TracBrowser for help on using the repository browser.