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

Last change on this file since 496 was 484, checked in by raasch, 15 years ago

typo in file headers removed

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