source: palm/trunk/SOURCE/diffusion_w.f90 @ 1350

Last change on this file since 1350 was 1341, checked in by kanani, 11 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 22.3 KB
RevLine 
[1]1 MODULE diffusion_w_mod
2
[1036]3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
[1310]17! Copyright 1997-2014 Leibniz Universitaet Hannover
[1036]18!--------------------------------------------------------------------------------!
19!
[484]20! Current revisions:
[1]21! -----------------
[1341]22!
23!
[1321]24! Former revisions:
25! -----------------
26! $Id: diffusion_w.f90 1341 2014-03-25 19:48:09Z maronga $
27!
[1341]28! 1340 2014-03-25 19:45:13Z kanani
29! REAL constants defined as wp-kind
30!
[1323]31! 1322 2014-03-20 16:38:49Z raasch
32! REAL constants defined as wp-kind
33!
[1321]34! 1320 2014-03-20 08:40:49Z raasch
[1320]35! ONLY-attribute added to USE-statements,
36! kind-parameters added to all INTEGER and REAL declaration statements,
37! kinds are defined in new module kinds,
38! revision history before 2012 removed,
39! comment fields (!:) to be used for variable explanations added to
40! all variable declaration statements
[1]41!
[1258]42! 1257 2013-11-08 15:18:40Z raasch
43! openacc loop and loop vector clauses removed, declare create moved after
44! the FORTRAN declaration statement
45!
[1132]46! 1128 2013-04-12 06:19:32Z raasch
47! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
48! j_north
49!
[1037]50! 1036 2012-10-22 13:43:42Z raasch
51! code put under GPL (PALM 3.9)
52!
[1017]53! 1015 2012-09-27 09:23:24Z raasch
54! accelerator version (*_acc) added
55!
[1002]56! 1001 2012-09-13 14:08:46Z raasch
57! arrays comunicated by module instead of parameter list
58!
[979]59! 978 2012-08-09 08:28:32Z fricke
60! outflow damping layer removed
61! kmxm_x/_z and kmxp_x/_z change to kmxm and kmxp
62! kmym_y/_z and kmyp_y/_z change to kmym and kmyp
63!
[1]64! Revision 1.1  1997/09/12 06:24:11  raasch
65! Initial revision
66!
67!
68! Description:
69! ------------
70! Diffusion term of the w-component
71!------------------------------------------------------------------------------!
72
[1320]73    USE wall_fluxes_mod,                                                       &
74        ONLY :  wall_fluxes, wall_fluxes_acc
[56]75
[1]76    PRIVATE
[1015]77    PUBLIC diffusion_w, diffusion_w_acc
[1]78
79    INTERFACE diffusion_w
80       MODULE PROCEDURE diffusion_w
81       MODULE PROCEDURE diffusion_w_ij
82    END INTERFACE diffusion_w
83
[1015]84    INTERFACE diffusion_w_acc
85       MODULE PROCEDURE diffusion_w_acc
86    END INTERFACE diffusion_w_acc
87
[1]88 CONTAINS
89
90
91!------------------------------------------------------------------------------!
92! Call for all grid points
93!------------------------------------------------------------------------------!
[1001]94    SUBROUTINE diffusion_w
[1]95
[1320]96       USE arrays_3d,                                                          &         
97           ONLY :  ddzu, ddzw, km, tend, u, v, w
98           
99       USE control_parameters,                                                 & 
100           ONLY :  topography
101           
102       USE grid_variables,                                                     &     
103           ONLY :  ddx, ddy, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y
104           
105       USE indices,                                                            &           
106           ONLY :  nxl, nxr, nyn, nys, nzb, nzb_w_inner, nzb_w_outer, nzt
107           
108       USE kinds
[1]109
110       IMPLICIT NONE
111
[1320]112       INTEGER(iwp) ::  i     !:
113       INTEGER(iwp) ::  j     !:
114       INTEGER(iwp) ::  k     !:
115       
116       REAL(wp) ::  kmxm  !:
117       REAL(wp) ::  kmxp  !:
118       REAL(wp) ::  kmym  !:
119       REAL(wp) ::  kmyp  !:
[1001]120
[1320]121       REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) ::  wsus  !:
122       REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) ::  wsvs  !:
[1]123
124
[56]125!
126!--    First calculate horizontal momentum flux w'u' and/or w'v' at vertical
127!--    walls, if neccessary
128       IF ( topography /= 'flat' )  THEN
[1320]129          CALL wall_fluxes( wsus, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, nzb_w_inner,             &
[56]130                            nzb_w_outer, wall_w_x )
[1320]131          CALL wall_fluxes( wsvs, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, nzb_w_inner,             &
[56]132                            nzb_w_outer, wall_w_y )
133       ENDIF
134
[1]135       DO  i = nxl, nxr
136          DO  j = nys, nyn
137             DO  k = nzb_w_outer(j,i)+1, nzt-1
138!
139!--             Interpolate eddy diffusivities on staggered gridpoints
[1322]140                kmxp = 0.25_wp *                                               &
[978]141                       ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
[1322]142                kmxm = 0.25_wp *                                               &
[978]143                       ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
[1322]144                kmyp = 0.25_wp *                                               &
[978]145                       ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
[1322]146                kmym = 0.25_wp *                                               &
[978]147                       ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1]148
149                tend(k,j,i) = tend(k,j,i)                                      &
[978]150                      & + ( kmxp * ( w(k,j,i+1)   - w(k,j,i)   ) * ddx         &
151                      &   + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1)   &
152                      &   - kmxm * ( w(k,j,i)   - w(k,j,i-1) ) * ddx           &
153                      &   - kmxm * ( u(k+1,j,i) - u(k,j,i)   ) * ddzu(k+1)     &
[1]154                      &   ) * ddx                                              &
[978]155                      & + ( kmyp * ( w(k,j+1,i)   - w(k,j,i)   ) * ddy         &
156                      &   + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1)   &
157                      &   - kmym * ( w(k,j,i)   - w(k,j-1,i) ) * ddy           &
158                      &   - kmym * ( v(k+1,j,i) - v(k,j,i)   ) * ddzu(k+1)     &
[1]159                      &   ) * ddy                                              &
[1322]160                      & + 2.0_wp * (                                           &
[1]161                      &   km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1)  &
162                      & - km(k,j,i)   * ( w(k,j,i)   - w(k-1,j,i) ) * ddzw(k)  &
[1340]163                      &            ) * ddzu(k+1)
[1]164             ENDDO
165
166!
167!--          Wall functions at all vertical walls, where necessary
[1340]168             IF ( wall_w_x(j,i) /= 0.0_wp  .OR.  wall_w_y(j,i) /= 0.0_wp )  THEN
[51]169
[1]170                DO  k = nzb_w_inner(j,i)+1, nzb_w_outer(j,i)
171!
172!--                Interpolate eddy diffusivities on staggered gridpoints
[1322]173                   kmxp = 0.25_wp *                                            &
[978]174                          ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
[1322]175                   kmxm = 0.25_wp *                                            &
[978]176                          ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
[1322]177                   kmyp = 0.25_wp *                                            &
[978]178                          ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
[1322]179                   kmym = 0.25_wp *                                            &
[978]180                          ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1]181
182                   tend(k,j,i) = tend(k,j,i)                                   &
183                                 + (   fwxp(j,i) * (                           &
[978]184                            kmxp * ( w(k,j,i+1)   - w(k,j,i)   ) * ddx         &
185                          + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1)   &
[1]186                                                   )                           &
187                                     - fwxm(j,i) * (                           &
[978]188                            kmxm * ( w(k,j,i)     - w(k,j,i-1) ) * ddx         &
189                          + kmxm * ( u(k+1,j,i)   - u(k,j,i)   ) * ddzu(k+1)   &
[1]190                                                   )                           &
[56]191                                     + wall_w_x(j,i) * wsus(k,j,i)             &
[1]192                                   ) * ddx                                     &
193                                 + (   fwyp(j,i) * (                           &
[978]194                            kmyp * ( w(k,j+1,i)   - w(k,j,i)   ) * ddy         &
195                          + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1)   &
[1]196                                                   )                           &
197                                     - fwym(j,i) * (                           &
[978]198                            kmym * ( w(k,j,i)     - w(k,j-1,i) ) * ddy         &
199                          + kmym * ( v(k+1,j,i)   - v(k,j,i)   ) * ddzu(k+1)   &
[1]200                                                   )                           &
[56]201                                     + wall_w_y(j,i) * wsvs(k,j,i)             &
[1]202                                   ) * ddy                                     &
[1340]203                                 + 2.0_wp * (                                  &
[1]204                           km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
205                         - km(k,j,i)   * ( w(k,j,i)   - w(k-1,j,i) ) * ddzw(k) &
[1340]206                                            ) * ddzu(k+1)
[1]207                ENDDO
208             ENDIF
209
210          ENDDO
211       ENDDO
212
213    END SUBROUTINE diffusion_w
214
215
216!------------------------------------------------------------------------------!
[1015]217! Call for all grid points - accelerator version
218!------------------------------------------------------------------------------!
219    SUBROUTINE diffusion_w_acc
220
[1320]221       USE arrays_3d,                                                          &
222           ONLY :  ddzu, ddzw, km, tend, u, v, w
223           
224       USE control_parameters,                                                 &
225           ONLY :  topography
226           
227       USE grid_variables,                                                     &
228           ONLY : ddx, ddy, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y
229           
230       USE indices,                                                            &
231           ONLY :  i_left, i_right, j_north, j_south, nxl, nxr, nyn, nys, nzb, &
232                   nzb_w_inner, nzb_w_outer, nzt
233                   
234       USE kinds
[1015]235
236       IMPLICIT NONE
237
[1320]238       INTEGER(iwp) ::  i     !:
239       INTEGER(iwp) ::  j     !:
240       INTEGER(iwp) ::  k     !:
241       
242       REAL(wp) ::  kmxm  !:
243       REAL(wp) ::  kmxp  !:
244       REAL(wp) ::  kmym  !:
245       REAL(wp) ::  kmyp  !:
[1015]246
[1320]247       REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) ::  wsus  !:
248       REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) ::  wsvs  !:
[1015]249       !$acc declare create ( wsus, wsvs )
250
251!
252!--    First calculate horizontal momentum flux w'u' and/or w'v' at vertical
253!--    walls, if neccessary
254       IF ( topography /= 'flat' )  THEN
[1320]255          CALL wall_fluxes_acc( wsus, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp,          &
256                                nzb_w_inner, nzb_w_outer, wall_w_x )
257          CALL wall_fluxes_acc( wsvs, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp,          &
258                                nzb_w_inner, nzb_w_outer, wall_w_y )
[1015]259       ENDIF
260
261       !$acc kernels present ( u, v, w, km, tend, vsws, vswst )    &
262       !$acc         present ( ddzu, ddzw, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y )           &
263       !$acc         present ( nzb_w_inner, nzb_w_outer )
[1128]264       DO  i = i_left, i_right
265          DO  j = j_south, j_north
[1015]266             DO  k = 1, nzt
267                IF ( k > nzb_w_outer(j,i) )  THEN
268!
269!--                Interpolate eddy diffusivities on staggered gridpoints
[1322]270                   kmxp = 0.25_wp *                                            &
[1015]271                          ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
[1322]272                   kmxm = 0.25_wp *                                            &
[1015]273                          ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
[1322]274                   kmyp = 0.25_wp *                                            &
[1015]275                          ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
[1322]276                   kmym = 0.25_wp *                                            &
[1015]277                          ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
278
279                   tend(k,j,i) = tend(k,j,i)                                     &
280                         & + ( kmxp * ( w(k,j,i+1)   - w(k,j,i)   ) * ddx        &
281                         &   + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1)  &
282                         &   - kmxm * ( w(k,j,i)   - w(k,j,i-1) ) * ddx          &
283                         &   - kmxm * ( u(k+1,j,i) - u(k,j,i)   ) * ddzu(k+1)    &
284                         &   ) * ddx                                             &
285                         & + ( kmyp * ( w(k,j+1,i)   - w(k,j,i)   ) * ddy        &
286                         &   + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1)  &
287                         &   - kmym * ( w(k,j,i)   - w(k,j-1,i) ) * ddy          &
288                         &   - kmym * ( v(k+1,j,i) - v(k,j,i)   ) * ddzu(k+1)    &
289                         &   ) * ddy                                             &
[1322]290                         & + 2.0_wp * (                                          &
[1015]291                         &   km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
292                         & - km(k,j,i)   * ( w(k,j,i)   - w(k-1,j,i) ) * ddzw(k) &
[1340]293                         &            ) * ddzu(k+1)
[1015]294                ENDIF
295             ENDDO
296
297!
298!--          Wall functions at all vertical walls, where necessary
299             DO  k = 1,nzt
300
301                IF ( k > nzb_w_inner(j,i)  .AND.  k <= nzb_w_outer(j,i)  .AND. &
[1340]302                     wall_w_x(j,i) /= 0.0_wp  .AND.  wall_w_y(j,i) /= 0.0_wp )  THEN
[1015]303!
304!--                Interpolate eddy diffusivities on staggered gridpoints
[1322]305                   kmxp = 0.25_wp *                                            &
[1015]306                          ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
[1322]307                   kmxm = 0.25_wp *                                            &
[1015]308                          ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
[1322]309                   kmyp = 0.25_wp *                                            &
[1015]310                          ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
[1322]311                   kmym = 0.25_wp *                                            &
[1015]312                          ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
313
314                   tend(k,j,i) = tend(k,j,i)                                   &
315                                 + (   fwxp(j,i) * (                           &
316                            kmxp * ( w(k,j,i+1)   - w(k,j,i)   ) * ddx         &
317                          + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1)   &
318                                                   )                           &
319                                     - fwxm(j,i) * (                           &
320                            kmxm * ( w(k,j,i)     - w(k,j,i-1) ) * ddx         &
321                          + kmxm * ( u(k+1,j,i)   - u(k,j,i)   ) * ddzu(k+1)   &
322                                                   )                           &
323                                     + wall_w_x(j,i) * wsus(k,j,i)             &
324                                   ) * ddx                                     &
325                                 + (   fwyp(j,i) * (                           &
326                            kmyp * ( w(k,j+1,i)   - w(k,j,i)   ) * ddy         &
327                          + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1)   &
328                                                   )                           &
329                                     - fwym(j,i) * (                           &
330                            kmym * ( w(k,j,i)     - w(k,j-1,i) ) * ddy         &
331                          + kmym * ( v(k+1,j,i)   - v(k,j,i)   ) * ddzu(k+1)   &
332                                                   )                           &
333                                     + wall_w_y(j,i) * wsvs(k,j,i)             &
334                                   ) * ddy                                     &
[1340]335                                 + 2.0_wp * (                                  &
[1015]336                           km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
337                         - km(k,j,i)   * ( w(k,j,i)   - w(k-1,j,i) ) * ddzw(k) &
[1340]338                                            ) * ddzu(k+1)
[1015]339                ENDIF
340             ENDDO
341
342          ENDDO
343       ENDDO
344       !$acc end kernels
345
346    END SUBROUTINE diffusion_w_acc
347
348
349!------------------------------------------------------------------------------!
[1]350! Call for grid point i,j
351!------------------------------------------------------------------------------!
[1001]352    SUBROUTINE diffusion_w_ij( i, j )
[1]353
[1320]354       USE arrays_3d,                                                          &         
355           ONLY :  ddzu, ddzw, km, tend, u, v, w
356           
357       USE control_parameters,                                                 & 
358           ONLY :  topography
359           
360       USE grid_variables,                                                     &     
361           ONLY :  ddx, ddy, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y
362           
363       USE indices,                                                            &           
364           ONLY :  nxl, nxr, nyn, nys, nzb, nzb_w_inner, nzb_w_outer, nzt
365           
366       USE kinds
[1]367
368       IMPLICIT NONE
369
[1320]370       INTEGER(iwp) ::  i     !:
371       INTEGER(iwp) ::  j     !:
372       INTEGER(iwp) ::  k     !:
373       
374       REAL(wp) ::  kmxm  !:
375       REAL(wp) ::  kmxp  !:
376       REAL(wp) ::  kmym  !:
377       REAL(wp) ::  kmyp  !:
[1]378
[1320]379       REAL(wp), DIMENSION(nzb:nzt+1) ::  wsus
380       REAL(wp), DIMENSION(nzb:nzt+1) ::  wsvs
[1]381
[1001]382
[1]383       DO  k = nzb_w_outer(j,i)+1, nzt-1
384!
385!--       Interpolate eddy diffusivities on staggered gridpoints
[1322]386          kmxp = 0.25_wp * ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
387          kmxm = 0.25_wp * ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
388          kmyp = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
389          kmym = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1]390
391          tend(k,j,i) = tend(k,j,i)                                            &
[978]392                      & + ( kmxp * ( w(k,j,i+1)   - w(k,j,i)   ) * ddx         &
393                      &   + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1)   &
394                      &   - kmxm * ( w(k,j,i)   - w(k,j,i-1) ) * ddx           &
395                      &   - kmxm * ( u(k+1,j,i) - u(k,j,i)   ) * ddzu(k+1)     &
[1]396                      &   ) * ddx                                              &
[978]397                      & + ( kmyp * ( w(k,j+1,i)   - w(k,j,i)   ) * ddy         &
398                      &   + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1)   &
399                      &   - kmym * ( w(k,j,i)   - w(k,j-1,i) ) * ddy           &
400                      &   - kmym * ( v(k+1,j,i) - v(k,j,i)   ) * ddzu(k+1)     &
[1]401                      &   ) * ddy                                              &
[1322]402                      & + 2.0_wp * (                                           &
[1]403                      &   km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1)  &
404                      & - km(k,j,i)   * ( w(k,j,i)   - w(k-1,j,i) ) * ddzw(k)  &
[1340]405                      &            ) * ddzu(k+1)
[1]406       ENDDO
407
408!
409!--    Wall functions at all vertical walls, where necessary
[1340]410       IF ( wall_w_x(j,i) /= 0.0_wp  .OR.  wall_w_y(j,i) /= 0.0_wp )  THEN
[51]411
412!
413!--       Calculate the horizontal momentum fluxes w'u' and/or w'v'
[1340]414          IF ( wall_w_x(j,i) /= 0.0_wp )  THEN
[1320]415             CALL wall_fluxes( i, j, nzb_w_inner(j,i)+1, nzb_w_outer(j,i),     &
416                               wsus, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp )
[51]417          ELSE
[1340]418             wsus = 0.0_wp
[51]419          ENDIF
420
[1340]421          IF ( wall_w_y(j,i) /= 0.0_wp )  THEN
[1320]422             CALL wall_fluxes( i, j, nzb_w_inner(j,i)+1, nzb_w_outer(j,i),     &
423                               wsvs, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp )
[51]424          ELSE
[1340]425             wsvs = 0.0_wp
[51]426          ENDIF
427
[1]428          DO  k = nzb_w_inner(j,i)+1, nzb_w_outer(j,i)
429!
430!--          Interpolate eddy diffusivities on staggered gridpoints
[1322]431             kmxp = 0.25_wp * ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
432             kmxm = 0.25_wp * ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
433             kmyp = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
434             kmym = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1]435
436             tend(k,j,i) = tend(k,j,i)                                         &
437                                 + (   fwxp(j,i) * (                           &
[978]438                            kmxp * ( w(k,j,i+1)   - w(k,j,i)   ) * ddx         &
439                          + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1)   &
[1]440                                                   )                           &
441                                     - fwxm(j,i) * (                           &
[978]442                            kmxm * ( w(k,j,i)     - w(k,j,i-1) ) * ddx         &
443                          + kmxm * ( u(k+1,j,i)   - u(k,j,i)   ) * ddzu(k+1)   &
[1]444                                                   )                           &
[51]445                                     + wall_w_x(j,i) * wsus(k)                 &
[1]446                                   ) * ddx                                     &
447                                 + (   fwyp(j,i) * (                           &
[978]448                            kmyp * ( w(k,j+1,i)   - w(k,j,i)   ) * ddy         &
449                          + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1)   &
[1]450                                                   )                           &
451                                     - fwym(j,i) * (                           &
[978]452                            kmym * ( w(k,j,i)     - w(k,j-1,i) ) * ddy         &
453                          + kmym * ( v(k+1,j,i)   - v(k,j,i)   ) * ddzu(k+1)   &
[1]454                                                   )                           &
[51]455                                     + wall_w_y(j,i) * wsvs(k)                 &
[1]456                                   ) * ddy                                     &
[1340]457                                 + 2.0_wp * (                                  &
[1]458                           km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
459                         - km(k,j,i)   * ( w(k,j,i)   - w(k-1,j,i) ) * ddzw(k) &
[1340]460                                            ) * ddzu(k+1)
[1]461          ENDDO
462       ENDIF
463
464    END SUBROUTINE diffusion_w_ij
465
[1323]466 END MODULE diffusion_w_mod
Note: See TracBrowser for help on using the repository browser.