source: palm/trunk/SOURCE/advec_w_up.f90 @ 4834

Last change on this file since 4834 was 4828, checked in by Giersch, 4 years ago

Copyright updated to year 2021, interface pmc_sort removed to accelarate the nesting code

  • Property svn:keywords set to Id
File size: 6.2 KB
Line 
1!> @file advec_w_up.f90
2!--------------------------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
6! Public License as published by the Free Software Foundation, either version 3 of the License, or
7! (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11! Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with PALM. If not, see
14! <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2021 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------------------------!
18!
19! Current revisions:
20! -----------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: advec_w_up.f90 4828 2021-01-05 11:21:41Z raasch $
26! file re-formatted to follow the PALM coding standard
27!
28! 4360 2020-01-07 11:25:50Z suehring
29! Corrected "Former revisions" section
30!
31! 3655 2019-01-07 16:51:22Z knoop
32! variables documented
33!
34! Revision 1.1  1997/08/29 08:56:33  raasch
35! Initial revision
36!
37!
38! Description:
39! ------------
40!> Advection term for the w velocity-component using upstream scheme.
41!> NOTE: vertical advection at k=1 still has wrong grid spacing for w>0
42!>       The same problem occurs for all topography boundaries!
43!--------------------------------------------------------------------------------------------------!
44 MODULE advec_w_up_mod
45 
46
47    PRIVATE
48    PUBLIC advec_w_up
49
50    INTERFACE advec_w_up
51       MODULE PROCEDURE advec_w_up
52       MODULE PROCEDURE advec_w_up_ij
53    END INTERFACE advec_w_up
54
55 CONTAINS
56
57
58!--------------------------------------------------------------------------------------------------!
59! Description:
60! ------------
61!> Call for all grid points
62!--------------------------------------------------------------------------------------------------!
63 SUBROUTINE advec_w_up
64
65    USE arrays_3d,                                                                                 &
66        ONLY:  ddzw, tend, u, v, w
67
68    USE control_parameters,                                                                        &
69        ONLY:  u_gtrans, v_gtrans
70
71    USE grid_variables,                                                                            &
72        ONLY:  ddx, ddy
73
74    USE indices,                                                                                   &
75        ONLY:  nxl, nxr, nyn, nys, nzb, nzt
76
77    USE kinds
78
79
80    IMPLICIT NONE
81
82    INTEGER(iwp) ::  i !< grid index along x-direction
83    INTEGER(iwp) ::  j !< grid index along y-direction
84    INTEGER(iwp) ::  k !< grid index along z-direction
85
86    REAL(wp) ::  ukomp !< advection velocity along x-direction
87    REAL(wp) ::  vkomp !< advection velocity along y-direction
88
89
90    DO  i = nxl, nxr
91       DO  j = nys, nyn
92          DO  k = nzb+1, nzt-1
93!
94!--          x-direction
95             ukomp = 0.25_wp * ( u(k,j,i)   + u(k,j,i+1) + u(k+1,j,i) + u(k+1,j,i+1) ) - u_gtrans
96             IF ( ukomp > 0.0_wp )  THEN
97                tend(k,j,i) = tend(k,j,i) - ukomp * ( w(k,j,i) - w(k,j,i-1) ) * ddx
98             ELSE
99                tend(k,j,i) = tend(k,j,i) - ukomp * ( w(k,j,i+1) - w(k,j,i) ) * ddx
100             ENDIF
101!
102!--          y-direction
103             vkomp = 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k+1,j,i) + v(k+1,j+1,i) ) - v_gtrans
104             IF ( vkomp > 0.0_wp )  THEN
105                tend(k,j,i) = tend(k,j,i) - vkomp * ( w(k,j,i) - w(k,j-1,i) ) * ddy
106             ELSE
107                tend(k,j,i) = tend(k,j,i) - vkomp * ( w(k,j+1,i) - w(k,j,i) ) * ddy
108             ENDIF
109!
110!--          z-direction
111             IF ( w(k,j,i) > 0.0_wp )  THEN
112                tend(k,j,i) = tend(k,j,i) - w(k,j,i) * ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k)
113             ELSE
114                tend(k,j,i) = tend(k,j,i) - w(k,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1)
115             ENDIF
116
117          ENDDO
118       ENDDO
119    ENDDO
120
121 END SUBROUTINE advec_w_up
122
123
124!--------------------------------------------------------------------------------------------------!
125! Description:
126! ------------
127!> Call for grid point i,j
128!--------------------------------------------------------------------------------------------------!
129 SUBROUTINE advec_w_up_ij( i, j )
130
131    USE arrays_3d,                                                                                 &
132        ONLY:  ddzw, tend, u, v, w
133
134    USE control_parameters,                                                                        &
135        ONLY:  u_gtrans, v_gtrans
136
137    USE grid_variables,                                                                            &
138        ONLY:  ddx, ddy
139
140    USE indices,                                                                                   &
141        ONLY:  nzb, nzt
142
143    USE kinds
144
145
146    IMPLICIT NONE
147
148    INTEGER(iwp) ::  i !< grid index along x-direction
149    INTEGER(iwp) ::  j !< grid index along y-direction
150    INTEGER(iwp) ::  k !< grid index along z-direction
151
152    REAL(wp) ::  ukomp !< advection velocity along x-direction
153    REAL(wp) ::  vkomp !< advection velocity along y-direction
154
155
156    DO  k = nzb+1, nzt-1
157!
158!--    x-direction
159       ukomp = 0.25_wp * ( u(k,j,i) + u(k,j,i+1) + u(k+1,j,i) + u(k+1,j,i+1) ) - u_gtrans
160       IF ( ukomp > 0.0_wp )  THEN
161          tend(k,j,i) = tend(k,j,i) - ukomp * ( w(k,j,i) - w(k,j,i-1) ) * ddx
162       ELSE
163          tend(k,j,i) = tend(k,j,i) - ukomp * ( w(k,j,i+1) - w(k,j,i) ) * ddx
164       ENDIF
165!
166!--    y-direction
167       vkomp = 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k+1,j,i) + v(k+1,j+1,i) ) - v_gtrans
168       IF ( vkomp > 0.0_wp )  THEN
169          tend(k,j,i) = tend(k,j,i) - vkomp * ( w(k,j,i) - w(k,j-1,i) ) * ddy
170       ELSE
171          tend(k,j,i) = tend(k,j,i) - vkomp * ( w(k,j+1,i) - w(k,j,i) ) * ddy
172       ENDIF
173!
174!--    z-direction
175       IF ( w(k,j,i) > 0.0_wp )  THEN
176          tend(k,j,i) = tend(k,j,i) - w(k,j,i) * ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k)
177       ELSE
178          tend(k,j,i) = tend(k,j,i) - w(k,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1)
179       ENDIF
180
181    ENDDO
182
183 END SUBROUTINE advec_w_up_ij
184
185 END MODULE advec_w_up_mod
Note: See TracBrowser for help on using the repository browser.